add hdbscan engine (#247)

* add hdbscan engine to `db_clust()` (#238)

* test hdbscan engine

* document hdbscan engine

* add NEWS bullet for hdbscan engine (#238)

Author: EmilHvitfeldt

NAMESPACE

@@ -14,6 +14,7 @@ S3method(extract_cluster_assignment,cluster_fit)
 S3method(extract_cluster_assignment,cluster_spec)
 S3method(extract_cluster_assignment,dbscan)
 S3method(extract_cluster_assignment,hclust)
+S3method(extract_cluster_assignment,hdbscan)
 S3method(extract_cluster_assignment,kmeans)
 S3method(extract_cluster_assignment,kmodes)
 S3method(extract_cluster_assignment,kproto)
@@ -27,6 +28,7 @@ S3method(extract_fit_summary,cluster_fit)
 S3method(extract_fit_summary,cluster_spec)
 S3method(extract_fit_summary,dbscan)
 S3method(extract_fit_summary,hclust)
+S3method(extract_fit_summary,hdbscan)
 S3method(extract_fit_summary,kmeans)
 S3method(extract_fit_summary,kmodes)
 S3method(extract_fit_summary,kproto)
@@ -93,6 +95,7 @@ S3method(update,k_means)
 S3method(update,mean_shift)
 export("%>%")
 export(.db_clust_fit_dbscan)
+export(.db_clust_fit_hdbscan)
 export(.gm_clust_fit_mclust)
 export(.hier_clust_fit_stats)
 export(.k_means_fit_ClusterR)

NEWS.md

@@ -28,6 +28,8 @@
 
 * `mean_shift()` gains a new engine with `meanShiftR`. (#244)
 
+* `db_clust()` gains a new engine with `hdbscan`, using `dbscan::hdbscan()` to fit HDBSCAN models. (#238)
+
 * The `.config` column produced by `tune_cluster()` has changed from the
   `Preprocessor{num}_Model{num}` pattern to `pre{num}_mod{num}_post{num}` to
   align with updates in the tune package. (#220)

R/db_clust.R

@@ -10,6 +10,7 @@
 #' are listed below.
 #'
 #' - \link[=details_db_clust_dbscan]{dbscan}
+#' - \link[=details_db_clust_hdbscan]{hdbscan}
 #'
 #' @param mode A single character string for the type of model. The only
 #'   possible value for this model is `"partition"`.
@@ -263,3 +264,39 @@ dbscan_helper <- function(object, ...) {
   )
   training_data$cluster[order(training_data$overall_order)]
 }
+
+#' Simple Wrapper around hdbscan function
+#'
+#' This wrapper passes the data to `dbscan::hdbscan()` and stashes the training
+#' data on the result so it can be reused for prediction and extraction.
+#'
+#' @param x matrix or data frame.
+#' @param min_points Minimum cluster size used as the `minPts` argument of
+#'   `dbscan::hdbscan()`.
+#' @param min_cluster_size Engine-specific override for `minPts`. When supplied,
+#'   it is used in place of `min_points`.
+#'
+#' @return hdbscan object
+#' @keywords internal
+#' @export
+.db_clust_fit_hdbscan <- function(
+  x,
+  min_points = NULL,
+  min_cluster_size = NULL,
+  ...
+) {
+  min_pts <- min_cluster_size %||% min_points
+
+  if (is.null(min_pts)) {
+    cli::cli_abort(
+      "Please specify `min_points` to be able to fit specification.",
+      call = call("fit")
+    )
+  }
+
+  res <- dbscan::hdbscan(x, minPts = min_pts, ...)
+  attr(res, "min_points") <- min_pts
+  attr(res, "training_data") <- x
+
+  res
+}

R/db_clust_data.R

@@ -78,6 +78,71 @@ make_db_clust <- function() {
       )
     )
   )
+
+  # ----------------------------------------------------------------------------
+
+  modelenv::set_model_engine("db_clust", "partition", "hdbscan")
+  modelenv::set_dependency(
+    model = "db_clust",
+    mode = "partition",
+    eng = "hdbscan",
+    pkg = "dbscan"
+  )
+  modelenv::set_dependency(
+    model = "db_clust",
+    mode = "partition",
+    eng = "hdbscan",
+    pkg = "tidyclust"
+  )
+
+  modelenv::set_fit(
+    model = "db_clust",
+    eng = "hdbscan",
+    mode = "partition",
+    value = list(
+      interface = "matrix",
+      protect = c("x", "min_points"),
+      func = c(pkg = "tidyclust", fun = ".db_clust_fit_hdbscan"),
+      defaults = list()
+    )
+  )
+
+  modelenv::set_encoding(
+    model = "db_clust",
+    eng = "hdbscan",
+    mode = "partition",
+    options = list(
+      predictor_indicators = "traditional",
+      compute_intercept = TRUE,
+      remove_intercept = TRUE,
+      allow_sparse_x = FALSE
+    )
+  )
+
+  modelenv::set_model_arg(
+    model = "db_clust",
+    eng = "hdbscan",
+    exposed = "min_points",
+    original = "min_points",
+    func = list(pkg = "dials", fun = "min_points"),
+    has_submodel = TRUE
+  )
+
+  modelenv::set_pred(
+    model = "db_clust",
+    eng = "hdbscan",
+    mode = "partition",
+    type = "cluster",
+    value = list(
+      pre = NULL,
+      post = NULL,
+      func = c(fun = ".db_clust_predict_hdbscan"),
+      args = list(
+        object = rlang::expr(object$fit),
+        new_data = rlang::expr(new_data)
+      )
+    )
+  )
 }
 
 # nocov end

R/db_clust_hdbscan.R

@@ -0,0 +1,11 @@
+#' Hierarchical Density-Based Spatial Clustering (HDBSCAN) via dbscan
+#'
+#' [db_clust()] creates an HDBSCAN model.
+#'
+#' @includeRmd man/rmd/db_clust_hdbscan.md details
+#'
+#' @name details_db_clust_hdbscan
+#' @keywords internal
+NULL
+
+# See inst/README-DOCS.md for a description of how these files are processed

R/extract_cluster_assignment.R

@@ -172,6 +172,13 @@ extract_cluster_assignment.dbscan <- function(object, ...) {
   cluster_assignment_tibble_w_outliers(clusters, n_clusters, ...)
 }
 
+#' @export
+extract_cluster_assignment.hdbscan <- function(object, ...) {
+  clusters <- object$cluster
+  n_clusters <- length(unique(clusters[clusters != 0])) + 1
+  cluster_assignment_tibble_w_outliers(clusters, n_clusters, ...)
+}
+
 #' @export
 extract_cluster_assignment.ms <- function(object, ...) {
   n_clusters <- nrow(object$cluster.center)

R/extract_fit_summary.R

@@ -275,6 +275,70 @@ extract_fit_summary.dbscan <- function(object, ...) {
   summary
 }
 
+#' @export
+extract_fit_summary.hdbscan <- function(object, ...) {
+  clusts <- extract_cluster_assignment(object, ...)$.cluster
+  n_clust <- dplyr::n_distinct(clusts)
+  training_data <- attr(object, "training_data")
+
+  overall_centroid <- colMeans(training_data)
+
+  by_clust <- training_data %>%
+    tibble::as_tibble() %>%
+    dplyr::mutate(
+      .cluster = clusts
+    ) %>%
+    dplyr::group_by(.cluster) %>%
+    tidyr::nest()
+
+  centroids <- by_clust$data %>%
+    map(dplyr::summarize_all, mean) %>%
+    dplyr::bind_rows()
+
+  outlier_idx <- which(unique(clusts) == "Outlier")
+
+  sse_within_total_total <- map2_dbl(
+    by_clust$data,
+    seq_len(n_clust),
+    ~ sum(
+      philentropy::dist_many_many(
+        as.matrix(centroids[.y, ]),
+        as.matrix(.x),
+        method = "euclidean"
+      )
+    )
+  )
+
+  clust_names <- unique(clusts)[c(outlier_idx, setdiff(1:n_clust, outlier_idx))]
+  summary <- list(
+    cluster_names = clust_names,
+    centroids = centroids,
+    n_members = unname(as.integer(table(clusts))),
+    sse_within_total_total = sse_within_total_total,
+    sse_total = sum(
+      philentropy::dist_many_many(
+        t(overall_centroid),
+        as.matrix(training_data),
+        method = "euclidean"
+      )
+    ),
+    orig_labels = NULL,
+    cluster_assignments = clusts
+  )
+
+  if (length(outlier_idx) > 0) {
+    summary$centroids[outlier_idx, ] <- rep(NA, ncol(summary$centroids))
+    summary$sse_within_total_total[outlier_idx] <- NA
+  }
+
+  # reorder centroids
+  summary$centroids <- summary$centroids[
+    c(outlier_idx, setdiff(1:n_clust, outlier_idx)),
+  ]
+
+  summary
+}
+
 #' @export
 extract_fit_summary.ms <- function(
   object,

R/predict_helpers.R

@@ -251,6 +251,21 @@ make_predictions_w_outliers <- function(x, prefix, n_clusters, labels = NULL) {
   make_predictions_w_outliers(clusters, prefix, n_clusters, labels)
 }
 
+.db_clust_predict_hdbscan <- function(
+  object,
+  new_data,
+  prefix = "Cluster_",
+  labels = NULL
+) {
+  training_data <- attr(object, "training_data")
+  clusters <- object$cluster
+  n_clusters <- length(unique(clusters[clusters != 0])) + 1
+
+  preds <- stats::predict(object, newdata = new_data, data = training_data)
+
+  make_predictions_w_outliers(preds, prefix, n_clusters, labels)
+}
+
 .mean_shift_predict_LPCM <- function(
   object,
   new_data,

R/tunable.R

@@ -100,6 +100,9 @@ tunable.db_clust <- function(x, ...) {
   if (x$engine == "dbscan") {
     res <- add_engine_parameters(res, dbscan_db_clust_engine_args)
   }
+  if (x$engine == "hdbscan") {
+    res <- add_engine_parameters(res, hdbscan_db_clust_engine_args)
+  }
   res
 }
 
@@ -118,6 +121,17 @@ dbscan_db_clust_engine_args <-
     component_id = "engine"
   )
 
+hdbscan_db_clust_engine_args <-
+  tibble::tibble(
+    name = "min_cluster_size",
+    call_info = list(
+      list(pkg = "dials", fun = "min_points")
+    ),
+    source = "cluster_spec",
+    component = "hdbscan",
+    component_id = "engine"
+  )
+
 
 #' @rdname tunable.cluster_spec
 #' @export