Hierarchical clustering and dendrograms in R: grouping sites and reading the tree

community ecology

ordination

A dissimilarity matrix turns into groups through hierarchical clustering. How to build a dendrogram with hclust, pick a linkage, cut it into clusters, and check the result against an NMDS ordination, in R.

Author

Tidy Ecology

Published

2026-06-21

An ordination like NMDS draws a continuous map: sites that sit close together have similar communities, and there are no hard borders. Hierarchical clustering takes the same dissimilarity information and does the opposite, carving the sites into discrete, nested groups and drawing the result as a tree. The two are complementary. Clustering is what you reach for when the question is “how many community types are here, and which site belongs to which?” rather than “how do communities grade into one another?”

This post builds a dendrogram in R from a Bray-Curtis dissimilarity matrix using hclust, walks through choosing a linkage method, cuts the tree into groups, and then checks those groups against an NMDS ordination of the same sites. That last step matters: a clustering you have not cross-checked is just a picture, and the quickest sanity test is whether an independent ordination agrees with it.

The ingredients

Hierarchical agglomerative clustering needs two decisions. The first is the dissimilarity between sites, which for community abundance data is usually Bray-Curtis (Bray and Curtis 1957), computed with vegan::vegdist. The second is the linkage rule, which decides how the distance between two groups of sites is measured once sites start merging. hclust then repeatedly fuses the two closest groups, from individual sites up to one tree.

A worked example

Here is a synthetic survey of 24 sites belonging to three community types, eight sites each. Types A and B share part of their species pool, so they should resemble each other more than either resembles C. The data is synthetic and illustrative.

library(vegan)
library(factoextra)
library(ggplot2)
library(dplyr)

set.seed(5)
types <- rep(c("A", "B", "C"), each = 8)
n_sites <- length(types); Sp <- 26

mean_mat <- matrix(0.4, nrow = n_sites, ncol = Sp)
for (i in seq_len(n_sites)) {
  t <- types[i]
  if (t == "A") mean_mat[i, 1:8]   <- 8     # A: species 1-8
  if (t == "B") mean_mat[i, 5:12]  <- 8     # B: species 5-12 (overlaps A)
  if (t == "C") mean_mat[i, 17:26] <- 8     # C: species 17-26 (distinct)
  mean_mat[i, ] <- mean_mat[i, ] * runif(Sp, 0.7, 1.3)
}
comm <- matrix(rpois(n_sites * Sp, mean_mat), nrow = n_sites)
rownames(comm) <- paste0(types, ave(seq_len(n_sites), types, FUN = seq_along))
colnames(comm) <- paste0("sp", seq_len(Sp))

bc     <- vegdist(comm, method = "bray")
hc_avg <- hclust(bc, method = "average")

The dendrogram below is read vertically. Each leaf is a site; the height at which two branches join is the Bray-Curtis dissimilarity at which those groups fuse. Low joins mean very similar sites; the high join near the top is the most distinct split in the data.

clusters    <- cutree(hc_avg, k = 3)
# one colour per cluster, reused in the ordination below for consistency
clust_cols  <- c("1" = "#c98a2e", "2" = "#3a7ca5", "3" = "#275139")
lr_order    <- unique(clusters[hc_avg$order])   # cluster order, left to right in the tree

fviz_dend(hc_avg, k = 3, k_colors = clust_cols[as.character(lr_order)],
          color_labels_by_k = TRUE, lwd = 0.8, cex = 0.62, rect = FALSE,
          main = "", ylab = "Bray-Curtis dissimilarity",
          ggtheme = theme_minimal(base_size = 12)) +
  scale_x_continuous(breaks = NULL) +
  theme(panel.grid.major.x = element_blank(),
        axis.text.x = element_blank(), axis.ticks.x = element_blank())

Dendrogram with three colored clusters. The A and B clusters join at a moderate height; the C cluster joins them only at the top. — Figure 1

The structure matches the design. The C sites form one cluster that stays separate until the very top, while A and B fuse with each other at a much lower height before joining C. The tree has recovered not just three groups but their relationship: A and B are sibling types within a broader division from C, which is exactly how the species pools were built.

Choosing a linkage

The linkage rule changes the shape of the tree, sometimes dramatically, so it is worth knowing what the common choices do. Single linkage measures group distance by the closest pair of members, which tends to produce straggly, chained clusters. Complete linkage uses the farthest pair, producing tight compact clusters. Average linkage (UPGMA) uses the mean distance between groups and sits between the two. Ward’s method (Ward 1963) merges the groups that least increase within-cluster variance, often giving clean, roughly equal-sized clusters.

There is a quantitative way to compare them: the cophenetic correlation (Sokal and Rohlf 1962), the correlation between the original dissimilarities and the distances implied by the tree. A higher value means the dendrogram distorts the real distances less.

linkages <- c("average", "complete", "single", "ward.D2")
coph <- sapply(linkages, function(m) cor(bc, cophenetic(hclust(bc, method = m))))
round(sort(coph, decreasing = TRUE), 3)

 average  ward.D2   single complete 
   0.986    0.984    0.981    0.976

For this dataset, average linkage preserves the distances best, which is common for Bray-Curtis community data and a reasonable default. It is not a universal rule: cophenetic correlation rewards faithfulness to the raw distances, while Ward’s method is often preferred when you specifically want compact, interpretable groups even at a small cost in fidelity. One R-specific trap is worth flagging: hclust offers both "ward.D" and "ward.D2", and only "ward.D2" implements Ward’s actual variance criterion on a distance matrix. Use "ward.D2" unless you have a particular reason not to.

Cutting the tree

A dendrogram on its own does not assign sites to groups; you have to cut it at some height, or equivalently ask for a fixed number of clusters with cutree. Cutting this tree into three recovers the planted types cleanly.

table(type = types, cluster = cutree(hc_avg, k = 3))

    cluster
type 1 2 3
   A 8 0 0
   B 0 8 0
   C 0 0 8

Choosing the number of clusters is the genuinely subjective part. The height of the joins is one guide: a long vertical gap before the next merge suggests a natural stopping point. Interpretability is another: clusters should correspond to something you can name. And the most reliable check is external, comparing the cut against an ordination.

Cross-checking with NMDS

If the clusters are real structure rather than artefacts of the cutting, they should hold up as coherent regions in an ordination that was computed independently. Here is an NMDS of the same sites, with points colored by their dendrogram cluster and shaped by their true type.

set.seed(5)
nmds <- metaMDS(comm, distance = "bray", k = 2, trymax = 50, trace = 0)

scores_df <- as.data.frame(scores(nmds, display = "sites"))
scores_df$site    <- rownames(scores_df)
scores_df$cluster <- factor(clusters[scores_df$site])
scores_df$type    <- types[match(scores_df$site, rownames(comm))]
hulls <- scores_df %>% group_by(cluster) %>% slice(chull(NMDS1, NMDS2))

ggplot(scores_df, aes(NMDS1, NMDS2, color = cluster)) +
  geom_polygon(data = hulls, aes(fill = cluster, group = cluster), alpha = 0.12, color = NA) +
  geom_point(aes(shape = type), size = 3, stroke = 1) +
  scale_color_manual(values = clust_cols, name = "cluster") +
  scale_fill_manual(values = clust_cols, guide = "none") +
  scale_shape_manual(values = c(A = 16, B = 17, C = 15), name = "true type") +
  annotate("text", x = Inf, y = -Inf, label = paste0("stress = ", round(nmds$stress, 3)),
           hjust = 1.1, vjust = -0.8, size = 3.2, color = "#5d6b61") +
  theme_minimal(base_size = 12) + theme(legend.position = "right")

NMDS ordination with three well-separated point clouds, each matching one dendrogram cluster and one true community type. — Figure 2

The three clusters land as three separate clouds, and within each cloud every point carries the same true-type shape. The dendrogram and the ordination, built from the same dissimilarities but displayed in completely different ways, tell the same story. That agreement is the evidence that the grouping is worth trusting. When clustering and ordination disagree, usually a cluster that the tree splits but the ordination shows as a single smear, it is a sign to revisit the number of clusters or the linkage.

Caveats

Hierarchical clustering will always return clusters, even when the data is a structureless cloud. The method cannot tell you whether groups exist, only how it would carve the data if they did. Treat the dendrogram as a hypothesis and validate it, with the cophenetic correlation, an ordination, or a test like ANOSIM or PERMANOVA on the resulting groups.

The number of clusters is a decision, not an output. Different cut heights give different group counts, and there is rarely one correct answer; report how you chose and check that the conclusion is stable across reasonable alternatives.

The linkage and the dissimilarity both shape the result, so state them. The same sites under single linkage and Ward’s method can produce visibly different trees. And as always with Bray-Curtis on raw abundances, a few hyperabundant species can dominate the dissimilarities; if that is a concern, transform the data (a square-root or a Hellinger transformation) before clustering.

References

Bray, J.R., and Curtis, J.T. (1957). An ordination of the upland forest communities of southern Wisconsin. Ecological Monographs, 27, 325-349.

Legendre, P., and Legendre, L. (2012). Numerical Ecology, 3rd English edition. Elsevier, Amsterdam.

Sokal, R.R., and Rohlf, F.J. (1962). The comparison of dendrograms by objective methods. Taxon, 11, 33-40.

Ward, J.H. Jr. (1963). Hierarchical grouping to optimize an objective function. Journal of the American Statistical Association, 58, 236-244.

--- title: "Hierarchical clustering and dendrograms in R: grouping sites and reading the tree" description: "A dissimilarity matrix turns into groups through hierarchical clustering. How to build a dendrogram with hclust, pick a linkage, cut it into clusters, and check the result against an NMDS ordination, in R." date: "2026-06-21 13:00" categories: [community ecology, ordination, R] image: thumbnail.png image-alt: "A dendrogram of survey sites split into three colored clusters, with two of the clusters joining before the third." --- An ordination like NMDS draws a continuous map: sites that sit close together have similar communities, and there are no hard borders. Hierarchical clustering takes the same dissimilarity information and does the opposite, carving the sites into discrete, nested groups and drawing the result as a tree. The two are complementary. Clustering is what you reach for when the question is "how many community types are here, and which site belongs to which?" rather than "how do communities grade into one another?" This post builds a dendrogram in R from a Bray-Curtis dissimilarity matrix using `hclust`, walks through choosing a linkage method, cuts the tree into groups, and then checks those groups against an NMDS ordination of the same sites. That last step matters: a clustering you have not cross-checked is just a picture, and the quickest sanity test is whether an independent ordination agrees with it. ## The ingredients Hierarchical agglomerative clustering needs two decisions. The first is the dissimilarity between sites, which for community abundance data is usually Bray-Curtis (Bray and Curtis 1957), computed with `vegan::vegdist`. The second is the linkage rule, which decides how the distance between two *groups* of sites is measured once sites start merging. `hclust` then repeatedly fuses the two closest groups, from individual sites up to one tree. ## A worked example Here is a synthetic survey of 24 sites belonging to three community types, eight sites each. Types A and B share part of their species pool, so they should resemble each other more than either resembles C. The data is synthetic and illustrative. ```{r} #| label: data #| message: false #| warning: false library(vegan) library(factoextra) library(ggplot2) library(dplyr) set.seed(5) types <- rep(c("A", "B", "C"), each = 8) n_sites <- length(types); Sp <- 26 mean_mat <- matrix(0.4, nrow = n_sites, ncol = Sp) for (i in seq_len(n_sites)) { t <- types[i] if (t == "A") mean_mat[i, 1:8] <- 8 # A: species 1-8 if (t == "B") mean_mat[i, 5:12] <- 8 # B: species 5-12 (overlaps A) if (t == "C") mean_mat[i, 17:26] <- 8 # C: species 17-26 (distinct) mean_mat[i, ] <- mean_mat[i, ] * runif(Sp, 0.7, 1.3) } comm <- matrix(rpois(n_sites * Sp, mean_mat), nrow = n_sites) rownames(comm) <- paste0(types, ave(seq_len(n_sites), types, FUN = seq_along)) colnames(comm) <- paste0("sp", seq_len(Sp)) bc <- vegdist(comm, method = "bray") hc_avg <- hclust(bc, method = "average") ``` The dendrogram below is read vertically. Each leaf is a site; the height at which two branches join is the Bray-Curtis dissimilarity at which those groups fuse. Low joins mean very similar sites; the high join near the top is the most distinct split in the data. ```{r} #| label: fig-dendrogram #| message: false #| warning: false #| fig-width: 9 #| fig-height: 5 #| fig-alt: "Dendrogram with three colored clusters. The A and B clusters join at a moderate height; the C cluster joins them only at the top." clusters <- cutree(hc_avg, k = 3) # one colour per cluster, reused in the ordination below for consistency clust_cols <- c("1" = "#c98a2e", "2" = "#3a7ca5", "3" = "#275139") lr_order <- unique(clusters[hc_avg$order]) # cluster order, left to right in the tree fviz_dend(hc_avg, k = 3, k_colors = clust_cols[as.character(lr_order)], color_labels_by_k = TRUE, lwd = 0.8, cex = 0.62, rect = FALSE, main = "", ylab = "Bray-Curtis dissimilarity", ggtheme = theme_minimal(base_size = 12)) + scale_x_continuous(breaks = NULL) + theme(panel.grid.major.x = element_blank(), axis.text.x = element_blank(), axis.ticks.x = element_blank()) ``` The structure matches the design. The C sites form one cluster that stays separate until the very top, while A and B fuse with each other at a much lower height before joining C. The tree has recovered not just three groups but their relationship: A and B are sibling types within a broader division from C, which is exactly how the species pools were built. ## Choosing a linkage The linkage rule changes the shape of the tree, sometimes dramatically, so it is worth knowing what the common choices do. Single linkage measures group distance by the closest pair of members, which tends to produce straggly, chained clusters. Complete linkage uses the farthest pair, producing tight compact clusters. Average linkage (UPGMA) uses the mean distance between groups and sits between the two. Ward's method (Ward 1963) merges the groups that least increase within-cluster variance, often giving clean, roughly equal-sized clusters. There is a quantitative way to compare them: the cophenetic correlation (Sokal and Rohlf 1962), the correlation between the original dissimilarities and the distances implied by the tree. A higher value means the dendrogram distorts the real distances less. ```{r} #| label: cophenetic #| message: false #| warning: false linkages <- c("average", "complete", "single", "ward.D2") coph <- sapply(linkages, function(m) cor(bc, cophenetic(hclust(bc, method = m)))) round(sort(coph, decreasing = TRUE), 3) ``` For this dataset, average linkage preserves the distances best, which is common for Bray-Curtis community data and a reasonable default. It is not a universal rule: cophenetic correlation rewards faithfulness to the raw distances, while Ward's method is often preferred when you specifically want compact, interpretable groups even at a small cost in fidelity. One R-specific trap is worth flagging: `hclust` offers both `"ward.D"` and `"ward.D2"`, and only `"ward.D2"` implements Ward's actual variance criterion on a distance matrix. Use `"ward.D2"` unless you have a particular reason not to. ## Cutting the tree A dendrogram on its own does not assign sites to groups; you have to cut it at some height, or equivalently ask for a fixed number of clusters with `cutree`. Cutting this tree into three recovers the planted types cleanly. ```{r} #| label: cutree #| message: false #| warning: false table(type = types, cluster = cutree(hc_avg, k = 3)) ``` Choosing the number of clusters is the genuinely subjective part. The height of the joins is one guide: a long vertical gap before the next merge suggests a natural stopping point. Interpretability is another: clusters should correspond to something you can name. And the most reliable check is external, comparing the cut against an ordination. ## Cross-checking with NMDS If the clusters are real structure rather than artefacts of the cutting, they should hold up as coherent regions in an ordination that was computed independently. Here is an NMDS of the same sites, with points colored by their dendrogram cluster and shaped by their true type. ```{r} #| label: fig-nmds #| message: false #| warning: false #| fig-width: 8.5 #| fig-height: 5 #| fig-alt: "NMDS ordination with three well-separated point clouds, each matching one dendrogram cluster and one true community type." set.seed(5) nmds <- metaMDS(comm, distance = "bray", k = 2, trymax = 50, trace = 0) scores_df <- as.data.frame(scores(nmds, display = "sites")) scores_df$site <- rownames(scores_df) scores_df$cluster <- factor(clusters[scores_df$site]) scores_df$type <- types[match(scores_df$site, rownames(comm))] hulls <- scores_df %>% group_by(cluster) %>% slice(chull(NMDS1, NMDS2)) ggplot(scores_df, aes(NMDS1, NMDS2, color = cluster)) + geom_polygon(data = hulls, aes(fill = cluster, group = cluster), alpha = 0.12, color = NA) + geom_point(aes(shape = type), size = 3, stroke = 1) + scale_color_manual(values = clust_cols, name = "cluster") + scale_fill_manual(values = clust_cols, guide = "none") + scale_shape_manual(values = c(A = 16, B = 17, C = 15), name = "true type") + annotate("text", x = Inf, y = -Inf, label = paste0("stress = ", round(nmds$stress, 3)), hjust = 1.1, vjust = -0.8, size = 3.2, color = "#5d6b61") + theme_minimal(base_size = 12) + theme(legend.position = "right") ``` The three clusters land as three separate clouds, and within each cloud every point carries the same true-type shape. The dendrogram and the ordination, built from the same dissimilarities but displayed in completely different ways, tell the same story. That agreement is the evidence that the grouping is worth trusting. When clustering and ordination disagree, usually a cluster that the tree splits but the ordination shows as a single smear, it is a sign to revisit the number of clusters or the linkage. ## Caveats Hierarchical clustering will always return clusters, even when the data is a structureless cloud. The method cannot tell you whether groups exist, only how it would carve the data if they did. Treat the dendrogram as a hypothesis and validate it, with the cophenetic correlation, an ordination, or a test like ANOSIM or PERMANOVA on the resulting groups. The number of clusters is a decision, not an output. Different cut heights give different group counts, and there is rarely one correct answer; report how you chose and check that the conclusion is stable across reasonable alternatives. The linkage and the dissimilarity both shape the result, so state them. The same sites under single linkage and Ward's method can produce visibly different trees. And as always with Bray-Curtis on raw abundances, a few hyperabundant species can dominate the dissimilarities; if that is a concern, transform the data (a square-root or a Hellinger transformation) before clustering. ## References Bray, J.R., and Curtis, J.T. (1957). An ordination of the upland forest communities of southern Wisconsin. *Ecological Monographs*, 27, 325-349. Legendre, P., and Legendre, L. (2012). *Numerical Ecology*, 3rd English edition. Elsevier, Amsterdam. Sokal, R.R., and Rohlf, F.J. (1962). The comparison of dendrograms by objective methods. *Taxon*, 11, 33-40. Ward, J.H. Jr. (1963). Hierarchical grouping to optimize an objective function. *Journal of the American Statistical Association*, 58, 236-244.