Comparing soft and hierarchical clusterings with element-centric similarity

In this vignette we will illustrate how element-centric similarity can be used to compare different kinds of clustering results: flat disjoint clusterings, flat overlapping clusterings, and hierarchical clusterings.

library(ClustAssess)
library(e1071)
library(dbscan)
library(ggplot2)
suppressPackageStartupMessages(library(dendextend))

# we will use the Iris dataset for this vignette
data = as.matrix(iris[,1:4])
df.iris = as.data.frame(prcomp(data)$x)
df.iris$species = iris$Species
ggplot(df.iris, aes(x=PC1, y=PC2, color=species)) + 
  geom_point() + 
  labs(title='Iris PCA')

# a flat, disjoint clustering with DBscan
db.res = dbscan(data, eps=1)$cluster
df.iris$db.cluster = as.factor(db.res)
ggplot(df.iris, aes(x=PC1, y=PC2, color=db.cluster)) + 
  geom_point() + 
  labs(title='DBScan clustering')

disjoint.clustering = create_clustering(db.res)

# an overlapping clustering with soft k-means
cmeans.res = cmeans(data, centers=6)$membership
# get the strongest cluster assignment for each observation to plot
df.iris$cmeans.cluster = as.factor(apply(cmeans.res, 1, which.max))
ggplot(df.iris, aes(x=PC1, y=PC2, color=cmeans.cluster)) + 
  geom_point() + 
  labs(title='c-means clustering')

# but we will still use the soft cluster assignments to calculate ECS
soft.clustering = create_clustering(cmeans.res)

# a hierarchical clustering
distances = dist(data, method='euclidean')
hc.res = hclust(distances, method='complete')
# plot the resulting dendrogram
node.colors = c('blue', 'red', 'green')
hc.res %>% as.dendrogram() %>% 
  set("leaves_pch", 19) %>% 
  set("leaves_cex", 0.5) %>% 
  set("leaves_col", node.colors[df.iris$species]) %>% 
  plot(main='Complete linkage hierarchical clustering', leaflab='none')

hierarchical.clustering = create_clustering(hc.res, r=0.8)

Now, we will compare the clustering results using element-centric similarity. ECS allows us to compare different kinds of clustering results, including the overlapping clustering and the hierarchical clustering we just calculated. The results tell us how similarly each observation was clustered by the two methods; the higher the ECS, the more similar the clustering of that observation.

# which observations are clustered more similarly?
# first compare flat disjoint and flat soft clusterings
df.iris$dbscan.cmeans.ecs = element_sim_elscore(disjoint.clustering, 
                                                soft.clustering)
ggplot(df.iris, aes(x=PC1, y=PC2, color=dbscan.cmeans.ecs)) + 
  geom_point() + 
  labs(title='DBScan vs c-means similarity')

mean(df.iris$dbscan.cmeans.ecs)
#> [1] 0.4675256

# next compare flat disjoint and hierarchical disjoint clusterings
df.iris$dbscan.hc.ecs = element_sim_elscore(disjoint.clustering, 
                                    hierarchical.clustering)
ggplot(df.iris, aes(x=PC1, y=PC2, color=dbscan.hc.ecs)) + 
  geom_point() + 
  labs(title='DBScan vs hclust similarity')

mean(df.iris$dbscan.hc.ecs)
#> [1] 0.5935388

# finally compare flat soft and hierarchical disjoint clusterings
df.iris$cmeans.hc.ecs = element_sim_elscore(soft.clustering, 
                                            hierarchical.clustering)
ggplot(df.iris, aes(x=PC1, y=PC2, color=cmeans.hc.ecs)) + 
  geom_point() + 
  labs(title='c-means vs hclust similarity')

mean(df.iris$cmeans.hc.ecs)
#> [1] 0.631482

Session info

sessionInfo()
#> R version 4.0.2 (2020-06-22)
#> Platform: x86_64-apple-darwin17.0 (64-bit)
#> Running under: macOS Catalina 10.15.7
#> 
#> Matrix products: default
#> BLAS:   /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
#> 
#> locale:
#> [1] C/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] dendextend_1.14.0 ggplot2_3.3.2     dbscan_1.1-5      e1071_1.7-4      
#> [5] ClustAssess_0.1.1 Seurat_3.1.4     
#> 
#> loaded via a namespace (and not attached):
#>   [1] TH.data_1.0-10      Rtsne_0.15          colorspace_1.4-1   
#>   [4] ellipsis_0.3.1      class_7.3-17        ggridges_0.5.3     
#>   [7] leiden_0.3.7        listenv_0.8.0       farver_2.0.3       
#>  [10] ggrepel_0.9.1       mvtnorm_1.1-1       mathjaxr_1.4-0     
#>  [13] codetools_0.2-16    splines_4.0.2       mnormt_2.0.2       
#>  [16] knitr_1.29          TFisher_0.2.0       jsonlite_1.7.1     
#>  [19] ica_1.0-2           cluster_2.1.0       png_0.1-7          
#>  [22] uwot_0.1.10         sctransform_0.3.2   compiler_4.0.2     
#>  [25] httr_1.4.2          Matrix_1.2-18       lazyeval_0.2.2     
#>  [28] htmltools_0.5.1.1   tools_4.0.2         rsvd_1.0.3         
#>  [31] igraph_1.2.5        gtable_0.3.0        glue_1.4.1         
#>  [34] RANN_2.6.1          reshape2_1.4.4      dplyr_1.0.0        
#>  [37] Rcpp_1.0.5          Biobase_2.50.0      vctrs_0.3.1        
#>  [40] multtest_2.46.0     ape_5.4-1           nlme_3.1-148       
#>  [43] lmtest_0.9-38       xfun_0.15           fastcluster_1.1.25 
#>  [46] stringr_1.4.0       globals_0.14.0      rbibutils_2.0      
#>  [49] lifecycle_0.2.0     irlba_2.3.3         future_1.21.0      
#>  [52] MASS_7.3-51.6       zoo_1.8-9           scales_1.1.1       
#>  [55] parallel_4.0.2      sandwich_3.0-0      RColorBrewer_1.1-2 
#>  [58] yaml_2.2.1          reticulate_1.18     pbapply_1.4-3      
#>  [61] gridExtra_2.3       stringi_1.4.6       mutoss_0.1-12      
#>  [64] plotrix_3.8-1       BiocGenerics_0.34.0 Rdpack_2.1.1       
#>  [67] rlang_0.4.10        pkgconfig_2.0.3     matrixStats_0.58.0 
#>  [70] evaluate_0.14       lattice_0.20-41     ROCR_1.0-11        
#>  [73] purrr_0.3.4         patchwork_1.0.1     htmlwidgets_1.5.1  
#>  [76] labeling_0.3        cowplot_1.1.1       tidyselect_1.1.0   
#>  [79] parallelly_1.24.0   RcppAnnoy_0.0.18    plyr_1.8.6         
#>  [82] magrittr_1.5        R6_2.4.1            generics_0.0.2     
#>  [85] multcomp_1.4-16     withr_2.2.0         pillar_1.4.6       
#>  [88] sn_1.6-2            fitdistrplus_1.1-3  survival_3.1-12    
#>  [91] tibble_3.0.3        future.apply_1.7.0  tsne_0.1-3         
#>  [94] crayon_1.3.4        KernSmooth_2.23-17  tmvnsim_1.0-2      
#>  [97] plotly_4.9.2.1      rmarkdown_2.3       viridis_0.5.1      
#> [100] grid_4.0.2          data.table_1.13.0   metap_1.4          
#> [103] digest_0.6.25       tidyr_1.1.0         numDeriv_2016.8-1.1
#> [106] stats4_4.0.2        munsell_0.5.0       viridisLite_0.3.0