Calculating {tidyclust} metrics using alternative distance functions

[bvancil]

I've created an example below that uses {tidyclust} and {tidymodels} to calculate resampling metrics. My hope was to use an l^\infty /supremum/Chebyshev norm rather than the conventional Euclidean one, but you can see by running the code that the metric values are too high to be consistent with the method="chebyshev" option for {philentropy}'s distance functions. I don't know how to specify a dist_fun in {tidyclust} other than within the model, and it doesn't seem to be used for calculating metrics. Thanks in advance for your time!

library(philentropy)
library(rlang)
library(tibble)
library(tidyclust)
library(tidymodels)

get_data <- function(dimension = 100, number_of_points = 30) {
  matrix(
    rnorm(number_of_points * dimension), 
    nrow = number_of_points, 
    ncol = dimension
  ) |>
  as_tibble(.name_repair = "unique") |> 
  set_names(paste0("x_", seq_len(dimension)))
}

get_metrics <- function(dataset, distance_method = "euclidean") {
  distance_function1 <- function(x) philentropy::distance(x, method = distance_method)
  distance_function2 <- function(x, y) philentropy::dist_many_many(x, y, method = distance_method)

  a_cv <- vfold_cv(dataset, v = 5)

  a_model <- hier_clust(
      num_clusters = tune(),
      linkage_method = "complete",
      dist_fun = distance_function1
    ) |> 
    set_engine('stats') |> 
    set_mode('partition')
  a_model

  a_recipe <- recipe(dataset) |> 
    update_role(starts_with("x_"), new_role = "predictor")
  summary(a_recipe)

  a_workflow <- workflow() |> 
    add_model(a_model) |> 
    add_recipe(a_recipe)

  n_clusters_grid <- grid_regular(num_clusters(range = c(10L, 25L)), levels = 16)

  a_resample <- tune_cluster(
    a_workflow,
    resamples = a_cv,
    grid = n_clusters_grid,
    control = control_grid(save_pred = TRUE, extract = identity),
    metrics = cluster_metric_set(
      sse_total,
      sse_within_total,
      sse_ratio,
      silhouette_avg
    )
  )

  a_resample_metrics <- a_resample |> 
    collect_metrics()

  a_resample_metrics
}

dataset <- get_data()
e1 <- get_metrics(dataset, distance_method = "euclidean")
# These are too high to be based on the $l^\infty$/supremum/Chebyshev norm.
c1 <- get_metrics(dataset, distance_method = "chebyshev")

[Emilhvitfeldt]

Thanks for the detailed write-up, this surfaced a couple of real things in tidyclust.

The key point is that dist_fun set on the model object isn't used when metrics are computed. Metrics calculate their own distances, independently of how the model was fit, and default to Euclidean. That's why your silhouette values looked Euclidean rather than Chebyshev.

To make a metric use an alternative distance, author a custom metric with new_cluster_metric() and pass it to cluster_metric_set(). A plain wrapper won't work on its own, since cluster_metric_set() requires each input to carry the cluster_metric class that new_cluster_metric() adds.

library(tidyclust)

cheby <- function(x) philentropy::distance(x, method = "chebyshev")

cheby_silhouette_avg <- new_cluster_metric(
  function(object, new_data = NULL, ...) {
    silhouette_avg(object, new_data = new_data, dist_fun = cheby, ...)
  },
  direction = "maximize"
)

my_metrics <- cluster_metric_set(cheby_silhouette_avg, sse_ratio)

This propagates through tune_cluster() too: pass the metric set as the metrics argument and your custom distance is used during resampling. I verified this end-to-end, and the Chebyshev results differ from Euclidean as expected.

I also found a bug along the way: two metrics wrapping the same built-in (e.g. Chebyshev and Euclidean silhouette) previously collided under the same .metric label and were silently averaged. That's now fixed in the development version, where each metric is labeled with the name you gave it.

See tidymodels/tidyclust#254 and #257 for details. Both are on the dev version:

pak::pak("tidymodels/tidyclust")

[bvancil]

Thanks, Emil! I'm posting the code for anyone else who wants to replicate the solution.

library(philentropy)
library(rlang)
library(tibble)
library(tidyclust)
library(tidymodels)

get_data <- function(dimension = 100, number_of_points = 30) {
  matrix(
    rnorm(number_of_points * dimension), 
    nrow = number_of_points, 
    ncol = dimension
  ) |>
  as_tibble(.name_repair = "unique") |> 
  set_names(paste0("x_", seq_len(dimension)))
}

get_metrics <- function(dataset, distance_method = "euclidean") {
  distance_function1 <- function(x) philentropy::distance(x, method = distance_method)
  distance_function2 <- function(x, y) philentropy::dist_many_many(x, y, method = distance_method)

  custom_silhouette_avg <- new_cluster_metric(
    function(object, new_data = NULL, ...) {
      silhouette_avg(object, new_data = new_data, dist_fun = distance_function1, ...)
    },
    direction = "maximize"
  )

  custom_sse_total <- new_cluster_metric(
    function(object, new_data = NULL, ...) {
      sse_total(object, new_data = new_data, dist_fun = distance_function2, ...)
    },
    direction = "zero"
  )

  custom_sse_within_total <- new_cluster_metric(
    function(object, new_data = NULL, ...) {
      sse_within_total(object, new_data = new_data, dist_fun = distance_function2, ...)
    },
    direction = "zero"
  )

  custom_sse_ratio <- new_cluster_metric(
    function(object, new_data = NULL, ...) {
      sse_ratio(object, new_data = new_data, dist_fun = distance_function2, ...)
    },
    direction = "zero"
  )
  
  a_metric_set <- cluster_metric_set(
    custom_sse_total,
    custom_sse_within_total,
    custom_sse_ratio,
    custom_silhouette_avg
  )

  a_cv <- vfold_cv(dataset, v = 5)

  a_model <- hier_clust(
      num_clusters = tune(),
      linkage_method = "complete",
      dist_fun = distance_function1
    ) |> 
    set_engine('stats') |> 
    set_mode('partition')
  a_model

  a_recipe <- recipe(dataset) |> 
    update_role(starts_with("x_"), new_role = "predictor")
  summary(a_recipe)

  a_workflow <- workflow() |> 
    add_model(a_model) |> 
    add_recipe(a_recipe)

  n_clusters_grid <- grid_regular(num_clusters(range = c(10L, 25L)), levels = 16)

  a_resample <- tune_cluster(
    a_workflow,
    resamples = a_cv,
    grid = n_clusters_grid,
    control = control_grid(save_pred = TRUE, extract = identity),
    metrics = a_metric_set
  )

  a_resample_metrics <- a_resample |> 
    collect_metrics()

  a_resample_metrics
}

dataset <- get_data()
e1 <- get_metrics(dataset, distance_method = "euclidean")
# These look right be based on the $l^\infty$/supremum/Chebyshev norm.
c1 <- get_metrics(dataset, distance_method = "chebyshev")