代码分享】“基因集单通路的泛癌GSEA富集分析

【代码分享]基因集单通路的泛癌GSEA富集分析#资料 如图

最近在整理TCGA多组学数据时，发现不少小伙伴对通路活性评估有需求。今天分享一个快速实现泛癌GSEA分析的方法，特别适合需要观察某个特定通路在多个癌症类型中激活状态的情况。这个方法不需要复杂的编程基础，咱们用R语言手把手走流程。

先看核心代码结构（完整脚本见文末github仓库）：

library(clusterProfiler) library(org.Hs.eg.db) # 示例数据预处理 pancan_data <- read.csv("tcga_pancan_expression.csv", row.names=1) gene_list <- pancan_data$log2FC names(gene_list) <- rownames(pancan_data) gene_list <- sort(gene_list, decreasing=TRUE) # 执行GSEA分析 gsea_res <- gseGO(geneList = gene_list, ont = "BP", OrgDb = org.Hs.eg.db, keyType = "SYMBOL", pvalueCutoff = 0.05)

这段代码的玄机在geneList参数的构造。注意这里用的是排序后的基因列表（我习惯用log2FC降序排列），这样GSEA算法才能正确计算富集得分。遇到过新手直接把差异基因扔进去跑，结果出不来富集通路的情况，多半是这个排序步骤被忽略了。

可视化部分推荐用阶梯图展示核心富集区域：

library(enrichplot) gseaplot2(gsea_res, geneSetID = "GO:0007049", title = "Cell Cycle Pathway", color = "#E64B35", base_size = 12)

这里geneSetID需要替换成自己关注的通路编号。有个小技巧：用color参数指定颜色时，建议避开默认的红色系，特别是在需要区分上下调通路时，用对比色更直观。

【代码分享]基因集单通路的泛癌GSEA富集分析#资料 如图

实际操作时发现不同癌症类型的数据结构差异较大，建议在循环分析前统一格式：

cancer_types <- unique(meta_data$cancer_type) for (ctype in cancer_types) { subtype_data <- filter(pancan_data, cancer_type == ctype) # 后续分析步骤... saveRDS(gsea_res, paste0("gsea_",ctype,".rds")) }

这个循环框架可以扩展到33种TCGA癌症类型。注意保存中间结果时建议用.rds格式，比.csv节省空间且保留数据结构完整。

遇到过最坑的报错是基因ID转换问题。当使用不同的注释数据库时，一定要确认keyType参数与实际数据匹配。比如从ENSEMBL转SYMBOL时经常会有版本号问题，这时需要：

library(AnnotationDbi) gene_symbol <- mapIds(org.Hs.eg.db, keys=rownames(expr_matrix), column="SYMBOL", keytype="ENSEMBL")

这种ID映射的细节处理直接影响分析结果的可靠性，建议每次转换后都做缺失值检查。

最后说下结果解读。GSEA的NES（标准化富集分数）绝对值大于1通常认为有生物学意义，但具体阈值需要结合样本量调整。个人经验是泛癌分析中NES>1.5的通路更值得关注。附上一个结果速查表制作方法：

signif_res <- gsea_res[gsea_res$qvalues < 0.1 & abs(gsea_res$NES) >1.5,] DT::datatable(signif_res[,c("ID","Description","NES","qvalue")], options=list(pageLength=10))

用交互式表格呈现结果，点击表头还能按列排序，比静态表格方便得多。这套流程跑完，基本就能快速锁定目标通路在不同癌种中的激活模式了。完整代码和测试数据已上传到github（见评论区），遇到报错欢迎带图踢我。