Clustering Methods
Clustering is an unsupervised method that groups similar observations together based on their feature values, without using predefined labels. The goal is to find structure in data where none has been imposed externally.
Two broad families of methods are covered here: hierarchical clustering, which builds a tree of nested groups, and k-means, which partitions data into a fixed number of groups. Both rely on a distance matrix, so the choice of distance metric is as important as the choice of algorithm.
Distance Metrics
Before clustering, you must define how dissimilarity between observations is measured. The right choice depends on your data type: Euclidean for scaled continuous variables, Bray-Curtis for species counts and OTU tables, correlation distance for expression profiles. A full reference with formulas and R code is on the Distance Metrics page.
Always scale continuous data before clustering
Variables on different scales will dominate distance calculations. Use scale() unless you are working with a distance metric designed for count data such as Bray-Curtis.
Hierarchical Clustering
Hierarchical clustering builds a tree (dendrogram) by iteratively merging the two most similar clusters. Starting with each observation as its own cluster, it continues until all observations are joined into one.
Linkage methods
The linkage method defines how distance between clusters is calculated as they grow:
| Method | Distance used | Tends to produce |
|---|---|---|
| Complete | Maximum between any two points | Compact, spherical clusters |
| Single | Minimum between any two points | Elongated chains |
| Average | Mean of all pairwise distances | Compromise |
| Ward.D2 | Minimises within-cluster variance | Equal-sized, compact clusters |
Ward's method is recommended for most biological applications. It merges the two clusters whose fusion produces the smallest increase in total within-cluster variance (error sum of squares), rather than defining inter-cluster distance directly from pairwise distances as single, complete, or average linkage do. It tends to produce compact, relatively balanced clusters and often yields dendrograms with pronounced height differences, making it easier to choose a cut level.
Because both Ward's method and k-means minimise within-cluster variance, they often produce similar solutions when applied to Euclidean data; disagreement between them can indicate weak or ambiguous cluster structure.
A key limitation is that Ward's method is derived for Euclidean distances. When working with non-Euclidean dissimilarities such as Bray-Curtis, average linkage is typically preferred.
Implementation
library(factoextra)
# Scale and compute distances
df_scaled <- scale(iris[, 1:4])
dist_matrix <- dist(df_scaled, method = "euclidean")
# Cluster with Ward's method
hc <- hclust(dist_matrix, method = "ward.D2")
# Plot dendrogram
plot(hc, main = "Hierarchical Clustering", xlab = "", sub = "")
# Cut into k clusters
clusters_hc <- cutree(hc, k = 3)
rect.hclust(hc, k = 3, border = 2:4)
Reading a dendrogram
- Height on the y-axis is the distance at which two clusters merge. Greater height means greater dissimilarity.
- Long branches before a merge indicate well-separated clusters.
- Where to cut the tree determines the number of clusters. Look for large jumps in height.
K-means Clustering
K-means partitions data into exactly k groups by minimising within-cluster variance. It assigns each observation to the nearest cluster centroid, recalculates centroids, and repeats until convergence.
K-means is faster than hierarchical clustering and scales well to large datasets, but requires you to specify k in advance and assumes clusters are roughly spherical.
# Always set a seed for reproducibility
set.seed(123)
km <- kmeans(df_scaled, centers = 3, nstart = 25)
# nstart = 25 runs the algorithm from 25 random starts
# and returns the best solution — always use nstart > 1
# Visualise
fviz_cluster(km, data = df_scaled,
palette = c("#E7B800", "#00AFBB", "#FC4E07"),
ggtheme = theme_minimal())
Limitations
K-means works poorly when clusters are elongated or irregularly shaped, when cluster sizes differ greatly, or when outliers are present. For these cases, hierarchical clustering with complete or Ward linkage is more robust, or consider PAM (see Going Further).
Choosing the Number of Clusters
No single method reliably identifies the true number of clusters. Use at least two approaches and consider biological plausibility alongside statistical criteria.
Elbow plot shows within-cluster sum of squares against k. Look for a bend where adding more clusters gives diminishing returns:
fviz_nbclust(df_scaled, kmeans, method = "wss", k.max = 10)
Silhouette analysis measures how similar each point is to its own cluster compared to neighbouring clusters. Values range from -1 to +1; higher is better. Choose k with the highest average silhouette width:
fviz_nbclust(df_scaled, kmeans, method = "silhouette", k.max = 10)
Gap statistic compares observed within-cluster dispersion to a null reference distribution. Choose the smallest k where the gap statistic is within one standard error of its maximum:
library(cluster)
gap_stat <- clusGap(df_scaled, FUN = kmeans, K.max = 10, B = 50)
fviz_gap_stat(gap_stat)
Practical advice
The three methods will not always agree. A biologically interpretable solution with k = 3 is often more useful than a statistically optimal solution with k = 7 that cannot be explained.
Validation
Clustering will always find groups, even in random data. Validation checks whether the clusters are meaningful.
Silhouette score is the most practical internal measure. Values above 0.5 indicate well-separated clusters; values below 0.25 suggest poor structure:
library(cluster)
sil <- silhouette(km$cluster, dist_matrix)
fviz_silhouette(sil)
mean(sil[, "sil_width"])
Comparison between methods
The Adjusted Rand Index mclust::adjustedRandIndex measures the similarity between two cluster assignments, where 1 represents perfect agreement and 0 indicates random chance. If hierarchical and k-means clustering broadly agree, the structure is likely real:
library(mclust)
clusters_km <- km$cluster
adjustedRandIndex(clusters_km, clusters_hc)
PERMANOVA tests whether clusters differ significantly in multivariate space:
library(vegan)
adonis2(dist_matrix ~ km$cluster, permutations = 999)
The output partitions total variation into Model (explained by cluster assignments) and Residual (unexplained).
- R² = 0.169 means the two clusters explain 16.9% of total variation. This is your effect size.
- F = 9.77, p = 0.001 means the separation is highly significant. Always report R² and p together.
One caveat: the clusters were defined on the same data being tested. This is circular validation! Significance confirms the clusters are distinct but says nothing about biological relevance.
External validation
If true labels are available, compare cluster assignments directly against known groups:
table(iris$Species, km$cluster)
adjustedRandIndex(iris$Species, km$cluster)
# 1.0 = perfect agreement, 0.0 = random
The table() shows which species land in which cluster. The Adjusted Rand Index (ARI) quantifies overall agreement on a scale from 0 to 1, where 1 is perfect agreement and 0 is no better than random. An ARI above 0.8 is generally considered strong agreement.
Note that this external validation is only possible when true labels are known. In most real clustering applications no ground truth exists, which is why internal validation metrics such as silhouette width are needed instead.
Silhouette width measures how well each sample fits its assigned cluster compared to neighbouring clusters.
For each sample it asks two questions: how similar is this sample to the other samples in its own cluster, and how similar is it to samples in the nearest other cluster. The silhouette width is the difference between those two numbers, scaled to a range from -1 to 1.
library(cluster)
sil <- silhouette(km$cluster, dist_matrix)
summary(sil)
plot(sil)
Interpreting the value:
Close to 1 sample fits well in its cluster
Close to 0 sample sits on the boundary between two clusters
Close to -1 sample probably belongs in the neighbouring cluster
Silhouette width answers the question you always want to ask after clustering: did the algorithm actually find meaningful groups or did it just divide the data arbitrarily? It is the most widely used internal validation metric precisely because it does not require any external labels.
Common Pitfalls
Not scaling data. Variables with large ranges dominate Euclidean distance.
Using Euclidean distance for OTU counts. Use Bray-Curtis instead.
Setting nstart = 1 in k-means. The algorithm may converge to a local minimum. Always use nstart = 25 or higher.
Over-interpreting clusters. Clustering is exploratory. Clusters found in the data are not necessarily biologically real — validate them.
Exercise
Using the USArrests dataset:
- Scale the data and compute Euclidean distances
- Cluster with Ward's method and plot the dendrogram
- Cut the tree into 4 clusters and identify which states group together
- Confirm with k-means using the same
k - Compare the two solutions with the Adjusted Rand Index (
mclust::adjustedRandIndex)
Solution
library(factoextra)
library(mclust)
## 1. Scale and compute distances
df_scaled <- scale(USArrests)
dist_matrix <- dist(df_scaled)
## 2. Hierarchical clustering
hc <- hclust(dist_matrix, method = "ward.D2")
plot(hc, main = "US Arrests")
## 3. Cut and inspect
clusters_hc <- cutree(hc, k = 4)
rect.hclust(hc, k = 4, border = 2:5)
lapply(1:4, function(i) names(clusters_hc[clusters_hc == i]))
## 4. K-means
set.seed(123)
km <- kmeans(df_scaled, centers = 4, nstart = 25)
## 5. Compare
adjustedRandIndex(clusters_hc, km$cluster)
# This is a metric used to compare two cluster assignments. It tells you how much the two methods agree on which data points belong together, corrected for chance [random (0) - perfect (1)].
# A score of ~0.82 indicates a very strong agreement between the two methods.
# This suggests that despite using different algorithms, both Hierarchical Clustering and K-Means identified very similar structures in your data. This is often a good sign that the clusters are robust and not just an artifact of a specific algorithm.