[WIP] HDBSCAN Clustering

cpp/CMakeLists.txt

@@ -424,6 +424,7 @@ if(BUILD_CUML_CPP_LIBRARY)
     src/kmeans/kmeans.cu
     src/knn/knn.cu
     src/knn/knn_sparse.cu
+    src/hierarchy/linkage.cu
     src/metrics/accuracy_score.cu
     src/metrics/adjusted_rand_index.cu
     src/metrics/completeness_score.cu

cpp/bench/CMakeLists.txt

@@ -25,6 +25,7 @@ if(BUILD_CUML_BENCH)
     sg/arima_loglikelihood.cu
     sg/dbscan.cu
     sg/kmeans.cu
+    sg/linkage.cu
     sg/main.cpp
     sg/rf_classifier.cu
     # FIXME: RF Regressor is having an issue where the tests now seem to take

cpp/bench/prims/fused_l2_nn.cu

@@ -1,5 +1,5 @@
 /*
- * Copyright (c) 2019-2020, NVIDIA CORPORATION.
+ * Copyright (c) 2019-2021, NVIDIA CORPORATION.
  *
  * Licensed under the Apache License, Version 2.0 (the "License");
  * you may not use this file except in compliance with the License.
@@ -15,8 +15,8 @@
  */
 
 #include <raft/cudart_utils.h>
-#include <distance/fused_l2_nn.cuh>
 #include <limits>
+#include <raft/distance/fused_l2_nn.cuh>
 #include <raft/linalg/norm.cuh>
 #include <raft/random/rng.cuh>
 #include "../common/ml_benchmark.hpp"
@@ -52,7 +52,7 @@ struct FusedL2NN : public Fixture {
     raft::linalg::rowNorm(yn, y, params.k, params.n, raft::linalg::L2Norm, true,
                           stream);
     auto blks = raft::ceildiv(params.m, 256);
-    MLCommon::Distance::initKernel<T, cub::KeyValuePair<int, T>, int>
+    raft::distance::initKernel<T, cub::KeyValuePair<int, T>, int>
       <<<blks, 256, 0, stream>>>(out, params.m, std::numeric_limits<T>::max(),
                                  op);
   }
@@ -69,9 +69,9 @@ struct FusedL2NN : public Fixture {
   void runBenchmark(::benchmark::State& state) override {
     loopOnState(state, [this]() {
       // it is enough to only benchmark the L2-squared metric
-      MLCommon::Distance::fusedL2NN<T, cub::KeyValuePair<int, T>, int>(
+      raft::distance::fusedL2NN<T, cub::KeyValuePair<int, T>, int>(
         out, x, y, xn, yn, params.m, params.n, params.k, (void*)workspace, op,
-        false, false, stream);
+        pairRedOp, false, false, stream);
     });
   }
 
@@ -80,7 +80,8 @@ struct FusedL2NN : public Fixture {
   T *x, *y, *xn, *yn;
   cub::KeyValuePair<int, T>* out;
   int* workspace;
-  MLCommon::Distance::MinAndDistanceReduceOp<int, T> op;
+  raft::distance::KVPMinReduce<int, T> pairRedOp;
+  raft::distance::MinAndDistanceReduceOp<int, T> op;
 };  // struct FusedL2NN
 
 static std::vector<FLNParams> getInputs() {

cpp/bench/sg/linkage.cu

@@ -0,0 +1,103 @@
+/*
+ * Copyright (c) 2019-2021, NVIDIA CORPORATION.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include <raft/linalg/distance_type.h>
+#include <raft/sparse/hierarchy/common.h>
+#include <cuml/cluster/linkage.hpp>
+#include <cuml/common/logger.hpp>
+#include <cuml/cuml.hpp>
+#include <utility>
+#include "benchmark.cuh"
+
+namespace ML {
+namespace Bench {
+namespace linkage {
+
+struct Params {
+  DatasetParams data;
+  BlobsParams blobs;
+};
+
+template <typename D>
+class Linkage : public BlobsFixture<D> {
+ public:
+  Linkage(const std::string& name, const Params& p)
+    : BlobsFixture<D>(name, p.data, p.blobs) {}
+
+ protected:
+  void runBenchmark(::benchmark::State& state) override {
+    using MLCommon::Bench::CudaEventTimer;
+    if (!this->params.rowMajor) {
+      state.SkipWithError("Single-Linkage only supports row-major inputs");
+    }
+
+    this->loopOnState(state, [this]() {
+      out_arrs.labels = labels;
+      out_arrs.children = out_children;
+
+      printf("RUNNING!\n");
+
+      Logger::get().setLevel(CUML_LEVEL_WARN);
+      ML::single_linkage_pairwise(
+        *this->handle, this->data.X, this->params.nrows, this->params.ncols,
+        raft::distance::DistanceType::L2Expanded, &out_arrs, 50, 50);
+    });
+  }
+
+  void allocateTempBuffers(const ::benchmark::State& state) override {
+    this->alloc(labels, this->params.nrows);
+    this->alloc(out_children, (this->params.nrows - 1) * 2);
+  }
+
+  void deallocateTempBuffers(const ::benchmark::State& state) override {
+    this->dealloc(labels, this->params.nrows);
+    this->dealloc(out_children, (this->params.nrows - 1) * 2);
+  }
+
+ private:
+  int* labels;
+  int* out_children;
+  raft::hierarchy::linkage_output<int, D> out_arrs;
+};
+
+std::vector<Params> getInputs() {
+  std::vector<Params> out;
+  Params p;
+  p.data.rowMajor = true;
+  p.blobs.cluster_std = 5.0;
+  p.blobs.shuffle = false;
+  p.blobs.center_box_min = -10.0;
+  p.blobs.center_box_max = 10.0;
+  p.blobs.seed = 12345ULL;
+  std::vector<std::pair<int, int>> rowcols = {
+    {35000, 128}, {16384, 128}, {12288, 128}, {8192, 128}, {4096, 128},
+  };
+  for (auto& rc : rowcols) {
+    p.data.nrows = rc.first;
+    p.data.ncols = rc.second;
+    for (auto nclass : std::vector<int>({1})) {
+      p.data.nclasses = nclass;
+      out.push_back(p);
+    }
+  }
+  return out;
+}
+
+ML_BENCH_REGISTER(Params, Linkage<float>, "blobs", getInputs());
+
+}  // namespace linkage
+}  // end namespace Bench
+}  // end namespace ML

cpp/cmake/Dependencies.cmake

@@ -39,7 +39,7 @@ else(DEFINED ENV{RAFT_PATH})
 
   ExternalProject_Add(raft
     GIT_REPOSITORY    https://github.com/rapidsai/raft.git
-    GIT_TAG           4a79adcb0c0e87964dcdc9b9122f242b5235b702
+    GIT_TAG           6455e05b3889db2b495cf3189b33c2b07bfbebf2
     PREFIX            ${RAFT_DIR}
     CONFIGURE_COMMAND ""
     BUILD_COMMAND     ""

cpp/include/cuml/cluster/hdbscan.hpp

@@ -0,0 +1,65 @@
+/*
+ * Copyright (c) 2018-2020, NVIDIA CORPORATION.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#pragma once
+
+#include <raft/linalg/distance_type.h>
+
+#include <cuml/common/logger.hpp>
+#include <cuml/cuml.hpp>
+
+namespace ML {
+
+template <typename value_idx, typename value_t>
+struct hdbscan_output {
+  value_idx m;
+  value_idx *labels;       // size: m
+  value_t *probabilities;  // size: m
+
+  value_idx *mst_src;  // size: m-1
+  value_idx *mst_dst;  // size: m-1
+  value_t *mst_data;   // size: m-1
+
+  value_idx *linakage_parents;  // size: m
+  value_idx *linkage_children;  // size: m
+  value_idx *linkage_deltas;    // size: m
+  value_idx *linkage_sizes;     // size: m
+};
+
+struct hdbscan_output_float : public hdbscan_output<int, float> {};
+
+/**
+ * @defgroup HdbscanCpp C++ implementation of Dbscan algo
+ * @brief Fits an HDBSCAN model on an input feature matrix and outputs the labels,
+ *        dendrogram, and minimum spanning tree.
+
+ * @param[in] handle
+ * @param[in] X
+ * @param[in] m
+ * @param[in] n
+ * @param[in] metric
+ * @param[in] k
+ * @param[in] min_pts
+ * @param[in] alpha
+ * @param[out] out
+ */
+void hdbscan(const raft::handle_t &handle, const float *X, int m, int n,
+             raft::distance::DistanceType metric, int k, int min_pts,
+             float alpha, hdbscan_output<int, float> *out);
+
+/** @} */
+
+}  // namespace ML

cpp/include/cuml/cluster/linkage.hpp

@@ -0,0 +1,60 @@
+/*
+ * Copyright (c) 2018-2021, NVIDIA CORPORATION.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#pragma once
+
+#include <raft/linalg/distance_type.h>
+#include <raft/sparse/hierarchy/common.h>
+
+#include <cuml/cuml.hpp>
+
+namespace ML {
+
+/**
+ * @defgroup HdbscanCpp C++ implementation of Dbscan algo
+ * @brief Fits an HDBSCAN model on an input feature matrix and outputs the labels,
+ *        dendrogram, and minimum spanning tree.
+ *  TODO: Use a separate type to represent number of edges so we can scale up
+ *  number of edges without having to use 64-bit ints for vertices.
+
+ * @param[in] handle
+ * @param[in] X
+ * @param[in] m
+ * @param[in] n
+ * @param[in] metric
+ * @param[out] out
+ */
+void single_linkage_pairwise(const raft::handle_t &handle, const float *X,
+                             size_t m, size_t n,
+                             raft::distance::DistanceType metric,
+                             raft::hierarchy::linkage_output<int, float> *out,
+                             int c = 15, int n_clusters = 5);
+
+void single_linkage_neighbors(const raft::handle_t &handle, const float *X,
+                              size_t m, size_t n,
+                              raft::distance::DistanceType metric,
+                              raft::hierarchy::linkage_output<int, float> *out,
+                              int c = 15, int n_clusters = 5);
+
+void single_linkage_pairwise(
+  const raft::handle_t &handle, const float *X, size_t m, size_t n,
+  raft::distance::DistanceType metric,
+  raft::hierarchy::linkage_output<int64_t, float> *out, int c = 15,
+  int n_clusters = 5);
+
+/** @} */
+
+};  // namespace ML

cpp/src/dbscan/adjgraph/algo.cuh

@@ -23,7 +23,7 @@
 #include "../common.cuh"
 #include "pack.h"
 
-#include <sparse/convert/csr.cuh>
+#include <raft/sparse/convert/csr.cuh>
 
 using namespace thrust;
 

cpp/src/dbscan/runner.cuh

@@ -21,7 +21,7 @@
 #include <cuml/common/device_buffer.hpp>
 #include <label/classlabels.cuh>
 #include <raft/cuda_utils.cuh>
-#include <sparse/csr.cuh>
+#include <raft/sparse/csr.cuh>
 #include "adjgraph/runner.cuh"
 #include "corepoints/compute.cuh"
 #include "corepoints/exchange.cuh"

cpp/src/hdbscan/hdbscan.cu

@@ -0,0 +1,36 @@
+/*
+ * Copyright (c) 2018-2020, NVIDIA CORPORATION.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include <cuml/cuml_api.h>
+#include <common/cumlHandle.hpp>
+#include <cuml/cluster/hdbscan.hpp>
+
+#include <hdbscan/runner.h>
+
+namespace ML {
+
+template <typename value_idx = int64_t, typename value_t = float>
+void hdbscan(const raft::handle_t &handle, value_t *X, size_t m, size_t n,
+             raft::distance::DistanceType metric, int k, int min_pts,
+             float alpha, hdbscan_output<value_idx, value_t> *out) {
+  HDBSCAN::_fit<value_idx, value_t>(handle, X, m, n, metric, k, min_pts, alpha);
+}
+
+void hdbscan(const raft::handle_t &handle, const float *X, size_t m, size_t n,
+             raft::distance::DistanceType metric, int k, int min_pts,
+             float alpha, hdbscan_output<int, float> *out);
+
+};  // end namespace ML