简介这份资源面向生物信息学入门学习者与R语言数据分析实践者围绕TCGA-BRCA乳腺癌基因表达数据展开聚类分析教学案例。内容涵盖利用层次聚类按基因表达水平对病人分型、选择average距离度量、绘制heatmap可视化聚类效果并实现PCA降维、依据主成分数目重新聚类再结合临床数据中的ER_Status_nature2012标签评估聚类是否符合预期帮助读者理解无监督学习在肿瘤分型中的应用逻辑。资源包共22个文件约10.91MB包含8个png与8个pdf图示、2个md说明、2个txt数据文件、1个license及1个R源码脚本图示与文档便于对照结果与复现流程。目前已有687人学习下载适合需要完整赛题式分析方案、可运行代码与可视化输出参考的读者。1. 从一份 TCGA-BRCA 聚类作业说起R 软件怎么把基因表达矩阵变成分型图如果你手头有一份基因表达矩阵行是基因、列是病人导师让你“按表达水平把病人分个类”你大概率会先愣一下从哪下手这份「生物信息学概论——聚类分析TCGA-BRCA数据.zip」就是干这个的。它把一套完整的 R 分析流程、示例数据、出图脚本和结果图打包在一起核心任务是用层次聚类和 PCA 对 TCGA-BRCA 的基因表达数据做无监督分型再用临床里的 ER 状态做外部验证。适合正在做生物信息课程设计、想跑通“表达矩阵 → 聚类 → 热图 → PCA 降维 → 再聚类”这条链的人。我拆完的第一感觉是代码不长但坑全在数据预处理和参数选择上下面按能复现的顺序讲。2. 数据拆包与矩阵对齐GeneMatrix 和 clinical_data 怎么接上拿到压缩包先别急着跑cluster.R把文件结构看清楚。包里有GeneMatrix.txt、clinical_data.txt、cluster.R、README.md以及Figures-PNG和Figures-PDF两个结果目录。图里包含cluster、pca_cluster、pca_screeplot、heatmap、ER_heatmap、pca_cumulative、pca_heatmap、ER_cluster这几类基本覆盖了作业要求的全部输出。2.1 两个核心文件的格式与读入方式GeneMatrix.txt是表达矩阵一行一个基因一列一个病人带行名和列名。clinical_data.txt是临床信息一行一个病人病人编号与表达矩阵的列名对应但只覆盖其中一部分。这意味着你不能直接merge完事得先取交集。# 读入表达矩阵check.namesFALSE 防止列名里的连字符被改掉 expr - read.table(GeneMatrix.txt, header TRUE, row.names 1, sep \t, check.names FALSE) # 读入临床信息 clin - read.table(clinical_data.txt, header TRUE, row.names 1, sep \t, check.names FALSE) # 取病人交集顺序以表达矩阵列为准 common_samples - intersect(colnames(expr), rownames(clin)) expr - expr[, common_samples] clin - clin[common_samples, , drop FALSE] # 确认没有全零或全 NA 的基因 expr - expr[rowSums(is.na(expr)) 0, ] expr - expr[rowSums(expr) 0, ]逻辑说明check.names FALSE很关键TCGA 的样本编号常带-R 默认会把连字符转成点号后面按编号取子集就会对不上。取交集后必须让clin的行顺序跟expr的列顺序一致否则后面 ER 状态和聚类标签会错位这种错位不会报错但热图上的注释条会全乱属于典型的“玄学翻车”。参数上sep \t对应制表符分隔如果你的文件是逗号分隔就改sep ,。row.names 1表示第一列是行名若第一列是数字索引而非基因名需要改成row.names NULL再手动赋值。2.2 表达值分布检查与标准化选择聚类对量纲敏感。如果矩阵里既有原始 counts 又有 TPM或者不同样本测序深度差异大直接算距离会被高表达基因主导。常见做法是先做log2(x 1)转换再按行做 z-score 标准化。# 判断是否需要 log 转换若最大值远大于 50 且分布右偏先 log summary(as.numeric(expr[1, ])) # log2 转换加 1 避免 log(0) expr_log - log2(expr 1) # 按基因行做 z-score让每个基因在样本间均值为 0、标准差为 1 expr_scaled - t(scale(t(expr_log))) # 检查是否有 NaN标准差为 0 的基因会产生 NaN sum(is.nan(expr_scaled)) expr_scaled - expr_scaled[!is.nan(rowSums(expr_scaled)), ]这里有个选择理由作业要求“按基因表达水平聚类病人”距离用 average即 UPGMA。如果不对基因做标准化高表达基因会主导距离矩阵聚类结果基本只反映那几个高变基因。z-score 之后每个基因权重接近聚类更能反映整体表达模式。注意scale默认按列所以要用t()转置两次这是 R 里最容易写反的一步。提示如果sum(is.nan(expr_scaled))不为 0说明有基因在所有样本中表达值完全相同标准化后标准差为 0必须剔除否则dist()会报 NA 错误。3. 层次聚类实操average 距离、树切分与热图输出数据对齐后进入核心分析。作业明确要求用层次聚类、average 距离并输出 heatmap。这一章把hclust、cutree和pheatmap的配合讲清楚。3.1 距离矩阵与 hclust 参数# 计算样本间距离用欧氏距离 dist_mat - dist(t(expr_scaled), method euclidean) # 层次聚类average 即 UPGMA hc - hclust(dist_mat, method average) # 查看聚类树 plot(hc, labels FALSE, main Hierarchical Clustering (average))逻辑说明dist作用于t(expr_scaled)因为我们要聚类的是病人列不是基因行。method average对应作业要求它的合并准则是对两类间所有点对距离取平均比 single 更稳健比 complete 更不容易被离群点拉偏。hclust的method还可选ward.D2如果你想让类内方差最小化ward 通常更紧凑但作业指定 average 就按 average 来。参数上dist的method可选euclidean、manhattan、correlation。对表达数据有人偏好correlation因为它关注表达模式形状而非绝对幅度。但作业没指定用欧氏距离配合 z-score 是常见组合。3.2 切树数目与热图注释切几类作业没给死数但临床里有 ER 状态可以拿它做参照。常见做法是看树高突变点或者直接切 2 到 4 类对比。# 切成 2 类 clusters - cutree(hc, k 2) # 与 ER 状态做交叉表 table(clusters, clin$ER_Status_nature2012) # 用 pheatmap 出图带 ER 注释条 library(pheatmap) annotation_col - data.frame(ER_Status clin$ER_Status_nature2012, Cluster factor(clusters)) rownames(annotation_col) - colnames(expr_scaled) pheatmap(expr_scaled, annotation_col annotation_col, show_rownames FALSE, show_colnames FALSE, clustering_distance_cols euclidean, clustering_method average, filename Figures-PNG/heatmap.png)逻辑说明cutree的k是切分簇数也可以给h按树高切。交叉表能直观看出聚类标签和 ER 状态是否吻合。pheatmap的annotation_col必须是有行名的数据框行名与矩阵列名一致否则注释条对不上。clustering_distance_cols和clustering_method要与前面hclust保持一致不然热图上的树和你的聚类结果不是一回事。参数上show_rownames FALSE是因为基因数通常上万显示会糊成一片。filename直接输出 PNG省去手动保存。包里的ER_heatmap.PNG应该就是加了 ER 注释的版本heatmap.PNG是不带注释的版本。注意pheatmap默认会对行做聚类如果你只想展示样本聚类把cluster_rows设为FALSE否则行树会占用大量绘图空间。4. PCA 降维与再聚类主成分数目怎么定、和层次聚类比什么作业的第二问是 PCA 降维后再聚类并与第一问比较。这一步的坑在于主成分数目选择选少了信息不够选多了等于没降维。4.1 PCA 实现与碎石图、累积方差图# 转置后做 PCAprcomp 默认中心化 pca_res - prcomp(t(expr_scaled), center TRUE, scale. FALSE) # 碎石图 png(Figures-PNG/pca_screeplot.png, width 800, height 600) plot(pca_res, type l, main Scree Plot) dev.off() # 累积方差贡献 var_explained - summary(pca_res)$importance[3, ] png(Figures-PNG/pca_cumulative.png, width 800, height 600) plot(var_explained, type b, xlab PC, ylab Cumulative Variance, main Cumulative Variance Explained) abline(h 0.8, col red, lty 2) dev.off()逻辑说明prcomp的scale. FALSE是因为前面已经做过 z-score再 scale 会重复标准化。碎石图看拐点累积方差图看达到 80% 需要几个主成分。常见做法是选累积方差到 80% 或碎石图拐点处的主成分数。包里的pca_screeplot.PNG和pca_cumulative.PNG就是这两张图。参数上center TRUE是必须的PCA 对中心化敏感。scale.如果设为TRUE会对每个基因再除标准差若前面没做 z-score这里设TRUE也能达到类似效果但两者叠加会过度标准化。4.2 用前几个主成分再聚类并对比# 假设选前 3 个主成分根据累积方差图调整 pc_use - 3 pca_data - pca_res$x[, 1:pc_use] # 对主成分得分做层次聚类 dist_pca - dist(pca_data, method euclidean) hc_pca - hclust(dist_pca, method average) clusters_pca - cutree(hc_pca, k 2) # 与原始聚类比较 table(Original clusters, PCA clusters_pca) # PCA 聚类热图 pheatmap(t(pca_data), annotation_col data.frame(Cluster factor(clusters_pca)), show_rownames TRUE, show_colnames FALSE, filename Figures-PNG/pca_heatmap.png)逻辑说明pca_res$x是样本在主成分上的得分取前pc_use列作为新特征。再聚类后与原始clusters做交叉表看一致性。如果 PCA 聚类结果与原始层次聚类高度一致说明前几个主成分已经抓住了主要变异如果差异大可能是主成分选少了或者数据本身没有明显分群。参数上pc_use的选择要有依据不能拍脑袋。我一般会看累积方差图选达到 80% 的最小主成分数同时参考碎石图拐点。包里的pca_cluster.PNG应该是 PCA 得分的散点图按聚类标签着色。提示PCA 聚类和原始聚类比较时簇的编号是任意的交叉表要允许标签置换不能只看对角线。5. 避坑与排查从数据错位到聚类不收敛的 5 个血泪经验这一章记录我复现这类作业时踩过的坑每条按现象、原因、解决写。现象一热图注释条颜色和聚类树完全对不上。原因clin取子集后行顺序没跟expr列顺序对齐或者annotation_col的行名与矩阵列名不一致。解决取交集后强制clin - clin[colnames(expr), , drop FALSE]并在pheatmap前用all(rownames(annotation_col) colnames(expr_scaled))检查返回FALSE就停下来查。现象二hclust报错“NA/NaN/Inf in foreign function call”。原因表达矩阵里有 NA、Inf或者 z-score 后出现 NaN。解决读入后先expr - expr[complete.cases(expr), ]再检查sum(is.nan(expr_scaled))把标准差为 0 的基因剔除。别用na.omit直接删样本会丢病人。现象三聚类结果每次跑都不一样。原因层次聚类本身是确定性的但如果前面用了随机初始化比如 k-means或者数据里有重复列导致距离矩阵不稳定。解决确认用的是hclust而非kmeans检查duplicated(t(expr_scaled))有重复样本就合并或剔除。现象四PCA 累积方差图显示前两个主成分就到 95%但聚类效果很差。原因数据没做标准化少数高表达基因主导了方差前几个主成分只反映那几个基因。解决回到 2.2 节做 z-score再重新 PCA。标准化前后各跑一次对比累积方差曲线你会看到明显差别。现象五ER 状态交叉表里 NA 一大堆。原因clinical_data.txt里 ER 状态列有缺失值或者病人编号匹配时因为连字符/点号问题没对上。解决读入时check.names FALSE匹配前用gsub(\\., -, colnames(expr))统一格式对 ER 缺失的样本在交叉表里单独列为 NA不要直接删否则样本量对不上。6. 进阶技巧用轮廓系数和 ER 一致性给聚类结果打分跑完图不算完作业要求“评估聚类是否符合预期”。除了看 ER 交叉表还可以算轮廓系数和调整兰德指数。轮廓系数衡量簇内紧密度和簇间分离度值越接近 1 越好调整兰德指数衡量聚类标签与 ER 状态的一致性越接近 1 越吻合。library(cluster) library(mclust) # 轮廓系数需要距离矩阵 sil - silhouette(clusters, dist_mat) png(Figures-PNG/silhouette.png, width 800, height 600) plot(sil, main Silhouette Plot) dev.off() mean(sil[, 3]) # 调整兰德指数比较聚类与 ER 状态 eri - adjustedRandIndex(clusters, clin$ER_Status_nature2012) print(paste(Adjusted Rand Index:, round(eri, 3)))逻辑说明silhouette的第一个参数是簇标签第二个是距离矩阵必须与聚类时用的距离一致。mean(sil[, 3])是平均轮廓系数大于 0.5 说明结构清晰0.2 到 0.5 算一般低于 0.2 就要怀疑是否真的存在分群。adjustedRandIndex对标签置换不敏感比直接看交叉表对角线更客观。参数上silhouette对大规模数据计算慢如果病人数上千可以先抽样。mclust包提供adjustedRandIndex如果不想装额外包也可以手写交叉表算 Rand 指数但调整版更严谨。我一般会把这几个指标和 ER 交叉表放一起看如果轮廓系数低但 ARI 高说明聚类结构弱但恰好和 ER 吻合如果两者都低说明这批数据在表达层面没有明显分型或者预处理有问题。从那以后我每次跑完聚类都强制走一遍“轮廓系数 ARI 交叉表”三件套不再只看热图颜色好不好看。希望帮到你。本文还有配套的精品资源点击获取 R语言聚类分析实战：TCGA-BRCA基因表达数据分型与PCA降维
简介这份资源面向生物信息学入门学习者与R语言数据分析实践者围绕TCGA-BRCA乳腺癌基因表达数据展开聚类分析教学案例。内容涵盖利用层次聚类按基因表达水平对病人分型、选择average距离度量、绘制heatmap可视化聚类效果并实现PCA降维、依据主成分数目重新聚类再结合临床数据中的ER_Status_nature2012标签评估聚类是否符合预期帮助读者理解无监督学习在肿瘤分型中的应用逻辑。资源包共22个文件约10.91MB包含8个png与8个pdf图示、2个md说明、2个txt数据文件、1个license及1个R源码脚本图示与文档便于对照结果与复现流程。目前已有687人学习下载适合需要完整赛题式分析方案、可运行代码与可视化输出参考的读者。1. 从一份 TCGA-BRCA 聚类作业说起R 软件怎么把基因表达矩阵变成分型图如果你手头有一份基因表达矩阵行是基因、列是病人导师让你“按表达水平把病人分个类”你大概率会先愣一下从哪下手这份「生物信息学概论——聚类分析TCGA-BRCA数据.zip」就是干这个的。它把一套完整的 R 分析流程、示例数据、出图脚本和结果图打包在一起核心任务是用层次聚类和 PCA 对 TCGA-BRCA 的基因表达数据做无监督分型再用临床里的 ER 状态做外部验证。适合正在做生物信息课程设计、想跑通“表达矩阵 → 聚类 → 热图 → PCA 降维 → 再聚类”这条链的人。我拆完的第一感觉是代码不长但坑全在数据预处理和参数选择上下面按能复现的顺序讲。2. 数据拆包与矩阵对齐GeneMatrix 和 clinical_data 怎么接上拿到压缩包先别急着跑cluster.R把文件结构看清楚。包里有GeneMatrix.txt、clinical_data.txt、cluster.R、README.md以及Figures-PNG和Figures-PDF两个结果目录。图里包含cluster、pca_cluster、pca_screeplot、heatmap、ER_heatmap、pca_cumulative、pca_heatmap、ER_cluster这几类基本覆盖了作业要求的全部输出。2.1 两个核心文件的格式与读入方式GeneMatrix.txt是表达矩阵一行一个基因一列一个病人带行名和列名。clinical_data.txt是临床信息一行一个病人病人编号与表达矩阵的列名对应但只覆盖其中一部分。这意味着你不能直接merge完事得先取交集。# 读入表达矩阵check.namesFALSE 防止列名里的连字符被改掉 expr - read.table(GeneMatrix.txt, header TRUE, row.names 1, sep \t, check.names FALSE) # 读入临床信息 clin - read.table(clinical_data.txt, header TRUE, row.names 1, sep \t, check.names FALSE) # 取病人交集顺序以表达矩阵列为准 common_samples - intersect(colnames(expr), rownames(clin)) expr - expr[, common_samples] clin - clin[common_samples, , drop FALSE] # 确认没有全零或全 NA 的基因 expr - expr[rowSums(is.na(expr)) 0, ] expr - expr[rowSums(expr) 0, ]逻辑说明check.names FALSE很关键TCGA 的样本编号常带-R 默认会把连字符转成点号后面按编号取子集就会对不上。取交集后必须让clin的行顺序跟expr的列顺序一致否则后面 ER 状态和聚类标签会错位这种错位不会报错但热图上的注释条会全乱属于典型的“玄学翻车”。参数上sep \t对应制表符分隔如果你的文件是逗号分隔就改sep ,。row.names 1表示第一列是行名若第一列是数字索引而非基因名需要改成row.names NULL再手动赋值。2.2 表达值分布检查与标准化选择聚类对量纲敏感。如果矩阵里既有原始 counts 又有 TPM或者不同样本测序深度差异大直接算距离会被高表达基因主导。常见做法是先做log2(x 1)转换再按行做 z-score 标准化。# 判断是否需要 log 转换若最大值远大于 50 且分布右偏先 log summary(as.numeric(expr[1, ])) # log2 转换加 1 避免 log(0) expr_log - log2(expr 1) # 按基因行做 z-score让每个基因在样本间均值为 0、标准差为 1 expr_scaled - t(scale(t(expr_log))) # 检查是否有 NaN标准差为 0 的基因会产生 NaN sum(is.nan(expr_scaled)) expr_scaled - expr_scaled[!is.nan(rowSums(expr_scaled)), ]这里有个选择理由作业要求“按基因表达水平聚类病人”距离用 average即 UPGMA。如果不对基因做标准化高表达基因会主导距离矩阵聚类结果基本只反映那几个高变基因。z-score 之后每个基因权重接近聚类更能反映整体表达模式。注意scale默认按列所以要用t()转置两次这是 R 里最容易写反的一步。提示如果sum(is.nan(expr_scaled))不为 0说明有基因在所有样本中表达值完全相同标准化后标准差为 0必须剔除否则dist()会报 NA 错误。3. 层次聚类实操average 距离、树切分与热图输出数据对齐后进入核心分析。作业明确要求用层次聚类、average 距离并输出 heatmap。这一章把hclust、cutree和pheatmap的配合讲清楚。3.1 距离矩阵与 hclust 参数# 计算样本间距离用欧氏距离 dist_mat - dist(t(expr_scaled), method euclidean) # 层次聚类average 即 UPGMA hc - hclust(dist_mat, method average) # 查看聚类树 plot(hc, labels FALSE, main Hierarchical Clustering (average))逻辑说明dist作用于t(expr_scaled)因为我们要聚类的是病人列不是基因行。method average对应作业要求它的合并准则是对两类间所有点对距离取平均比 single 更稳健比 complete 更不容易被离群点拉偏。hclust的method还可选ward.D2如果你想让类内方差最小化ward 通常更紧凑但作业指定 average 就按 average 来。参数上dist的method可选euclidean、manhattan、correlation。对表达数据有人偏好correlation因为它关注表达模式形状而非绝对幅度。但作业没指定用欧氏距离配合 z-score 是常见组合。3.2 切树数目与热图注释切几类作业没给死数但临床里有 ER 状态可以拿它做参照。常见做法是看树高突变点或者直接切 2 到 4 类对比。# 切成 2 类 clusters - cutree(hc, k 2) # 与 ER 状态做交叉表 table(clusters, clin$ER_Status_nature2012) # 用 pheatmap 出图带 ER 注释条 library(pheatmap) annotation_col - data.frame(ER_Status clin$ER_Status_nature2012, Cluster factor(clusters)) rownames(annotation_col) - colnames(expr_scaled) pheatmap(expr_scaled, annotation_col annotation_col, show_rownames FALSE, show_colnames FALSE, clustering_distance_cols euclidean, clustering_method average, filename Figures-PNG/heatmap.png)逻辑说明cutree的k是切分簇数也可以给h按树高切。交叉表能直观看出聚类标签和 ER 状态是否吻合。pheatmap的annotation_col必须是有行名的数据框行名与矩阵列名一致否则注释条对不上。clustering_distance_cols和clustering_method要与前面hclust保持一致不然热图上的树和你的聚类结果不是一回事。参数上show_rownames FALSE是因为基因数通常上万显示会糊成一片。filename直接输出 PNG省去手动保存。包里的ER_heatmap.PNG应该就是加了 ER 注释的版本heatmap.PNG是不带注释的版本。注意pheatmap默认会对行做聚类如果你只想展示样本聚类把cluster_rows设为FALSE否则行树会占用大量绘图空间。4. PCA 降维与再聚类主成分数目怎么定、和层次聚类比什么作业的第二问是 PCA 降维后再聚类并与第一问比较。这一步的坑在于主成分数目选择选少了信息不够选多了等于没降维。4.1 PCA 实现与碎石图、累积方差图# 转置后做 PCAprcomp 默认中心化 pca_res - prcomp(t(expr_scaled), center TRUE, scale. FALSE) # 碎石图 png(Figures-PNG/pca_screeplot.png, width 800, height 600) plot(pca_res, type l, main Scree Plot) dev.off() # 累积方差贡献 var_explained - summary(pca_res)$importance[3, ] png(Figures-PNG/pca_cumulative.png, width 800, height 600) plot(var_explained, type b, xlab PC, ylab Cumulative Variance, main Cumulative Variance Explained) abline(h 0.8, col red, lty 2) dev.off()逻辑说明prcomp的scale. FALSE是因为前面已经做过 z-score再 scale 会重复标准化。碎石图看拐点累积方差图看达到 80% 需要几个主成分。常见做法是选累积方差到 80% 或碎石图拐点处的主成分数。包里的pca_screeplot.PNG和pca_cumulative.PNG就是这两张图。参数上center TRUE是必须的PCA 对中心化敏感。scale.如果设为TRUE会对每个基因再除标准差若前面没做 z-score这里设TRUE也能达到类似效果但两者叠加会过度标准化。4.2 用前几个主成分再聚类并对比# 假设选前 3 个主成分根据累积方差图调整 pc_use - 3 pca_data - pca_res$x[, 1:pc_use] # 对主成分得分做层次聚类 dist_pca - dist(pca_data, method euclidean) hc_pca - hclust(dist_pca, method average) clusters_pca - cutree(hc_pca, k 2) # 与原始聚类比较 table(Original clusters, PCA clusters_pca) # PCA 聚类热图 pheatmap(t(pca_data), annotation_col data.frame(Cluster factor(clusters_pca)), show_rownames TRUE, show_colnames FALSE, filename Figures-PNG/pca_heatmap.png)逻辑说明pca_res$x是样本在主成分上的得分取前pc_use列作为新特征。再聚类后与原始clusters做交叉表看一致性。如果 PCA 聚类结果与原始层次聚类高度一致说明前几个主成分已经抓住了主要变异如果差异大可能是主成分选少了或者数据本身没有明显分群。参数上pc_use的选择要有依据不能拍脑袋。我一般会看累积方差图选达到 80% 的最小主成分数同时参考碎石图拐点。包里的pca_cluster.PNG应该是 PCA 得分的散点图按聚类标签着色。提示PCA 聚类和原始聚类比较时簇的编号是任意的交叉表要允许标签置换不能只看对角线。5. 避坑与排查从数据错位到聚类不收敛的 5 个血泪经验这一章记录我复现这类作业时踩过的坑每条按现象、原因、解决写。现象一热图注释条颜色和聚类树完全对不上。原因clin取子集后行顺序没跟expr列顺序对齐或者annotation_col的行名与矩阵列名不一致。解决取交集后强制clin - clin[colnames(expr), , drop FALSE]并在pheatmap前用all(rownames(annotation_col) colnames(expr_scaled))检查返回FALSE就停下来查。现象二hclust报错“NA/NaN/Inf in foreign function call”。原因表达矩阵里有 NA、Inf或者 z-score 后出现 NaN。解决读入后先expr - expr[complete.cases(expr), ]再检查sum(is.nan(expr_scaled))把标准差为 0 的基因剔除。别用na.omit直接删样本会丢病人。现象三聚类结果每次跑都不一样。原因层次聚类本身是确定性的但如果前面用了随机初始化比如 k-means或者数据里有重复列导致距离矩阵不稳定。解决确认用的是hclust而非kmeans检查duplicated(t(expr_scaled))有重复样本就合并或剔除。现象四PCA 累积方差图显示前两个主成分就到 95%但聚类效果很差。原因数据没做标准化少数高表达基因主导了方差前几个主成分只反映那几个基因。解决回到 2.2 节做 z-score再重新 PCA。标准化前后各跑一次对比累积方差曲线你会看到明显差别。现象五ER 状态交叉表里 NA 一大堆。原因clinical_data.txt里 ER 状态列有缺失值或者病人编号匹配时因为连字符/点号问题没对上。解决读入时check.names FALSE匹配前用gsub(\\., -, colnames(expr))统一格式对 ER 缺失的样本在交叉表里单独列为 NA不要直接删否则样本量对不上。6. 进阶技巧用轮廓系数和 ER 一致性给聚类结果打分跑完图不算完作业要求“评估聚类是否符合预期”。除了看 ER 交叉表还可以算轮廓系数和调整兰德指数。轮廓系数衡量簇内紧密度和簇间分离度值越接近 1 越好调整兰德指数衡量聚类标签与 ER 状态的一致性越接近 1 越吻合。library(cluster) library(mclust) # 轮廓系数需要距离矩阵 sil - silhouette(clusters, dist_mat) png(Figures-PNG/silhouette.png, width 800, height 600) plot(sil, main Silhouette Plot) dev.off() mean(sil[, 3]) # 调整兰德指数比较聚类与 ER 状态 eri - adjustedRandIndex(clusters, clin$ER_Status_nature2012) print(paste(Adjusted Rand Index:, round(eri, 3)))逻辑说明silhouette的第一个参数是簇标签第二个是距离矩阵必须与聚类时用的距离一致。mean(sil[, 3])是平均轮廓系数大于 0.5 说明结构清晰0.2 到 0.5 算一般低于 0.2 就要怀疑是否真的存在分群。adjustedRandIndex对标签置换不敏感比直接看交叉表对角线更客观。参数上silhouette对大规模数据计算慢如果病人数上千可以先抽样。mclust包提供adjustedRandIndex如果不想装额外包也可以手写交叉表算 Rand 指数但调整版更严谨。我一般会把这几个指标和 ER 交叉表放一起看如果轮廓系数低但 ARI 高说明聚类结构弱但恰好和 ER 吻合如果两者都低说明这批数据在表达层面没有明显分型或者预处理有问题。从那以后我每次跑完聚类都强制走一遍“轮廓系数 ARI 交叉表”三件套不再只看热图颜色好不好看。希望帮到你。本文还有配套的精品资源点击获取