R语言聚类分析实战：层次聚类与K-means一网打尽

生态博士的R笔记
22 阅读5分钟

R语言聚类分析实战：层次聚类与K-means一网打尽

一文掌握群落数据分析的核心方法：用R语言实现层次聚类（UPGMA）与K-means聚类的完整流程，包含距离计算、树状图绘制、可视化及结果解读，附完整可运行代码。

一、引言

物种组成相似的样方往往反映了相似的生境条件——聚类分析正是量化这种相似性的核心工具。在生态学中，层次聚类帮助我们理解群落之间的层级关系，K-means则能快速将样方划分为离散的群落类型。本文以一份30样方×20物种的仿真群落数据为例，系统演示从距离计算、层次聚类到K-means聚类的完整流程，帮你掌握群落分类的基本分析方法。

二、数据准备

本文使用仿真生成的群落多度数据，包含30个样方的20个物种：

数据结构：

  • 30行（样方），每行一个样方
  • 20列（物种），每列一个物种
  • 值为物种多度（个体数）
  • 样方名：Plot01 ~ Plot30
  • 物种名：Sp01 ~ Sp20

数据设计为3个群落类型（每种10个样方），每个类型有不同的优势物种组合：

类型样方优势物种特征
类型APlot01~10Sp01~Sp08物种1-8多度较高
类型BPlot11~20Sp09~Sp15物种9-15多度较高
类型CPlot21~30Sp16~Sp20物种16-20多度较高
# 读取数据
df <- read.csv("community_data.csv", header = TRUE, row.names = 1)
head(df[, 1:6])  # 查看前6行前6列

配图1.png

三、层次聚类（Hierarchical Clustering）

第一步：计算距离矩阵

聚类分析的第一步是计算样方之间的距离（相异度）。对于物种多度数据，最常用的是 Jaccard距离（支持多度和有无数据）或 Bray-Curtis距离（基于多度）。

library(vegan)
dist_jac <- vegdist(df, method = "jaccard")

vegdist() 是vegan包的距离计算函数，method = "jaccard" 基于物种有无和相对多度计算相异度。

第二步：UPGMA聚类

有了距离矩阵后，用 hclust() 进行层次聚类：

hc_upgma <- hclust(dist_jac, method = "average")

method = "average" 即UPGMA（非加权组平均法），是最常用的聚类方法之一。其他常用选项：

参数值方法特点
"average"UPGMA（非加权组平均）最常用，平衡各分支
"ward.D2"Ward最小方差法倾向于形成紧凑的球状簇
"complete"完全连接法基于最远点距离
"single"单连接法基于最近点距离，易产生链式效应

第三步：绘制树状图

# 默认树状图（垂下式）
par(mar = c(4, 4, 3, 8))
plot(hc_upgma, hang = -1, main = "Hierarchical Clustering (UPGMA, Jaccard)",
     xlab = "Sample", ylab = "Jaccard Distance", sub = "", cex = 0.7)
rect.hclust(hc_upgma, k = 3, border = c("#E64B35", "#4DBBD5", "#00A087"))

参数解析：

  • hang = -1：将所有分支末端对齐到同一水平线，树状图更整洁
  • rect.hclust(k = 3)：在3类处画矩形框，颜色区分不同簇
  • cex = 0.7：样方标签字体大小，样方多时调小避免重叠

配图2.png

# 水平树状图（适合样方多的情况）
plot(as.dendrogram(hc_upgma), horiz = TRUE,
     main = "Hierarchical Clustering (UPGMA, Jaccard)")

配图3.png

第四步：切割聚类结果

# 在3类处切割
hc_clusters <- cutree(hc_upgma, k = 3)
print(table(hc_clusters))

cutree() 函数将树状图在指定高度处切割，返回每个样方所属的类别编号。

四、K-means 聚类

层次聚类不要求预先指定类别数就能看到层级关系，但当你知道大致有几类时，K-means 能更快速地给出分类结果。

第一步：降维

K-means直接在高维物种数据上表现不佳，建议先降维：

dist_bc <- vegdist(df, method = "bray")
pcoa <- cmdscale(dist_bc, k = 4, eig = TRUE)

用 cmdscale() 对Bray-Curtis距离矩阵做PCoA降维，取前4个轴作为K-means的输入。

第二步：K-means聚类

set.seed(123)
km <- kmeans(pcoa$points, centers = 3, nstart = 50)
km$cluster  # 查看每个样方的分类

参数解析：

  • centers = 3：指定分为3类（可根据先验知识或肘部法则确定）
  • nstart = 50：随机初始中心50次，取最优结果，避免局部最优
  • set.seed(123)：保证结果可复现

第三步：可视化

library(factoextra)
print(fviz_cluster(km, data = pcoa$points,
                   ellipse.type = "norm",
                   palette = c("#E64B35", "#4DBBD5", "#00A087"),
                   ggtheme = theme_bw(),
                   main = "K-means Clustering (k = 3)",
                   labelsize = 0) +
  labs(color = "Cluster", shape = "Cluster"))

配图4.png

fviz_cluster() 自动将数据投影到前两个主坐标轴上展示聚类结果，ellipse.type = "norm" 添加正态置信椭圆。

五、结果解读

从层次聚类树状图和K-means聚类图中可以清楚看到：

1. 三个群落类型区分明显 树状图中三个主要分支清晰可辨（红/蓝/绿色矩形框），对应仿真数据的三类群落结构。这说明基于物种组成的Jaccard距离能有效区分不同生境类型的样方。

2. UPGMA vs K-means 分类结果基本一致 两种方法对30个样方的分类结果高度吻合——层次聚类基于距离矩阵的完整层级信息，K-means基于PCoA降维后的坐标空间，两种方法从不同角度得到了相似的群落划分。

3. K-means可视化直观 fviz_cluster() 的散点图能直观展示三类样方在主坐标空间中的分布，椭圆置信区间帮助判断各类的紧致程度。

4. 聚类分析的应用场景：

  • 根据物种组成划分群落类型
  • 验证预先定义的生境/处理分组是否与物种组成分类一致
  • 为后续的排序分析（NMDS、RDA）提供分类变量

六、完整可运行代码

以下为完整R脚本，复制到RStudio即可运行：

# ============================================================
# 聚类分析完整流程：层次聚类 + K-means
# ============================================================
# 首次运行请先安装包：
# install.packages(c("vegan", "factoextra"))

library(vegan)
library(factoextra)

# --- 读取数据 ---
df <- read.csv("community_data.csv", header = TRUE, row.names = 1)

# --- 层次聚类 ---
dist_jac <- vegdist(df, method = "jaccard")
hc_upgma <- hclust(dist_jac, method = "average")

# 垂下式树状图
pdf("dendrogram.pdf", width = 12, height = 7)
par(mar = c(4, 4, 3, 8))
plot(hc_upgma, hang = -1, main = "Hierarchical Clustering (UPGMA, Jaccard)",
     xlab = "Sample", ylab = "Jaccard Distance", sub = "", cex = 0.7)
rect.hclust(hc_upgma, k = 3, border = c("#E64B35", "#4DBBD5", "#00A087"))
dev.off()

# 水平树状图
pdf("dendrogram_horizontal.pdf", width = 8, height = 10)
par(mar = c(4, 4, 2, 2))
plot(as.dendrogram(hc_upgma), horiz = TRUE,
     main = "Hierarchical Clustering (UPGMA, Jaccard)")
dev.off()

# 切割聚类结果
hc_clusters <- cutree(hc_upgma, k = 3)

# --- K-means 聚类 ---
dist_bc <- vegdist(df, method = "bray")
pcoa <- cmdscale(dist_bc, k = 4, eig = TRUE)
set.seed(123)
km <- kmeans(pcoa$points, centers = 3, nstart = 50)

# K-means可视化
pdf("kmeans_clusters.pdf", width = 8, height = 6)
print(fviz_cluster(km, data = pcoa$points,
                   ellipse.type = "norm",
                   palette = c("#E64B35", "#4DBBD5", "#00A087"),
                   ggtheme = theme_bw(),
                   main = "K-means Clustering (k = 3)",
                   labelsize = 0) +
  labs(color = "Cluster", shape = "Cluster"))
dev.off()

# 聚类结果对比
data.frame(Sample = rownames(df),
           Hierarchical = paste0("G", hc_clusters),
           Kmeans = paste0("G", km$cluster))

替换为自己数据的步骤：

  1. 准备一个CSV文件：行=样方，列=物种，值为多度
  2. 第一列为样方名，第一行为物种名
  3. 将 read.csv("community_data.csv") 替换为你的文件名
  4. 根据你的研究调整 k = 3 中的类别数

七、实用技巧总结

方法选择速查：

想知道样方之间的层级关系？→ 层次聚类（UPGMA） 已有预期类别数，想快速分类？→ K-means 样方数>50？→ K-means更高效，层次聚类树状图太拥挤 数据中零值很多？→ Jaccard距离 优于欧氏距离

距离方法选择：

数据类型推荐距离函数
物种多度（原始值）Bray-Curtisvegdist(method = "bray")
物种有无（0/1）或多度Jaccardvegdist(method = "jaccard")
环境变量（连续）欧氏距离dist(method = "euclidean")
物种多度（经标准化的）欧氏距离dist(method = "euclidean")

注意事项：

  • 层次聚类前对数据做标准化（如Hellinger转化），可减少稀有物种的权重影响
  • K-means对初始中心敏感，务必设置 nstart >= 50
  • 聚类结果建议用两种不同方法交叉验证，结论更可靠
  • 样方命名避免使用纯数字（如"1","2"），加前缀（如"Plot1","S1"）避免R自动转为行号

如果对你有帮助，欢迎 点赞 + 收藏。欢迎评论区交流 🙌

avatar
文章
阅读
粉丝