企业官网建设流程全解析

简介：这份压缩包是面向生物信息学与科研论文写作的R语言GSEA多组富集图绘制工具包，针对高通量转录组数据中基因集富集分析结果的SCI级可视化需求，适用于需要批量比较KEGG、GO等多通路并统一输出图表的科研人员。压缩包内共五百一十三份文件，以二百九十三张png富集图、一百一十份xls统计结果、一百零二份html交互式报告为主体，同时配套r脚本、rnk排序文件、gmt基因集、cls样本分组与rpt运行记录，整体大小约四十一点六三MB，目录分区明确，便于按用途直接查找。目前已有四百二十三人学习或下载，适合希望快速复现多组GSEA绘图流程的初学者，也适合需要调整细节的进阶用户。使用者可直接基于输入文件和R代码批量生成富集图，学习fgsea、clusterProfiler等工具的参数调节、富集方向标注、阈值筛选与结果整理方法，还可参考已有png图表样式，对照xls结果核实富集得分和显著性，从而系统掌握多组GSEA可视化的完整分析与绘图链路，既便于论文组图，也方便后续二次开发，大幅缩短SCI图表制作周期。

1. 多GSEA富集图到底是什么：一个标题背后的三种输入与整版出图需求

审稿意见里最常见的补充要求之一，就是“补充多组GSEA富集图”。很多人第一反应是：单条通路的GSEA图我会画，多组对比的怎么排？实际上多GSEA富集图不是一个函数画一张图，而是一个完整流程：从表达矩阵、分组表、基因集这三份输入出发，用R对多个比较组分别跑GSEA，再拼成一版能直接放进SCI论文的图板。这套流程适合刚做完RNA-seq或芯片差异分析、准备写文章的人，也适合需要批量复查多个实验对比结果的老手。先说结论：多数人在这上面翻车，问题不在画图函数，而在输入文件没统一、基因ID类型不一致；把输入理顺，代码反而只占半小时。R语言在这里的主场优势很明显——从富集计算到拼图出PDF，全程不需要离开R环境。

2. 输入文件怎么备：表达矩阵、分组表和GMT基因集的组织方式

2.1 表达矩阵的两种形态：raw counts 与 TPM/FPKM，GSEA要的是哪种

GSEA的输入其实不是一个矩阵，而是一个按指标降序排列的基因列表（geneList）。这个列表的数值最常见来源是组间差异的log2FC，也可以用信号噪声比。所以表达矩阵喂给下游代码之前，先想清楚指标从哪来：用limma做差异得到logFC很常规，直接从表达矩阵算组均值比也常见。

# 假设拿到了一个“基因×样本”的TPM矩阵 expr <- read.delim("expression_matrix.txt", row.names = 1, check.names = FALSE) # 分组表：两列，sample和group group <- read.delim("sample_group.txt", stringsAsFactors = FALSE) # 按分组计算两组均值，取log2倍变化 treat_mean <- rowMeans(expr[, group$sample[group$group == "treat"]]) ctrl_mean <- rowMeans(expr[, group$sample[group$group == "control"]]) log2fc <- log2((treat_mean + 0.01) / (ctrl_mean + 0.01)) geneList <- sort(log2fc, decreasing = TRUE) head(geneList, 20)

这段代码的逻辑：先把样本按分组表归类，算出每个基因在两组里的平均表达量，再转成log2倍变化。+0.01是一个伪计数，避免除零。如果你用limma或DESeq2拿到了正式的logFC，直接替换这里的计算会更稳。sort(..., decreasing=TRUE)不能省，GSEA要求geneList按数值从大到小排列，否则后面富集分数整体错位。GSEA吃的是排序而不是P值，所以用log2FC做排序值是最常见的可靠做法。

2.2 分组信息表的标准组织：列名、样本名与配对设计

分组表是连接表达矩阵和比较组的桥梁。我最常踩的坑是列名不统一：有的叫SampleName，有的叫condition，代码里写group$sample就全报错。建议进代码前统一成两列：第一列sample，第二列group。如果做的是配对样本，比如同一批样本用药前后，分组表里还要有第三列pair_id。

# 推荐的分组表结构 # sample group pair_id # S01_treat treat P01 # S01_ctrl control P01 table(group$group) # 检查每个组的样本量，至少3个，少于3个时GSEA排序不稳定

这一小步看起来无关紧要，但它决定了后面循环里group$group == "treat"能不能直接复用。样本量方面，每组3个是底线。少于3个我不建议做GSEA，排序列表噪声太大，富集结果很难复现，不如改用单基因集富集分析看趋势。

2.3 基因集的获取与GMT格式：Hallmark、KEGG与自定义通路

GSEA的第三个输入是基因集。最常见的是MSigDB Hallmark集合，适合看核心生物学过程；KEGG通路适合看代谢与信号通路。GMT文件每行代表一个基因集：第一列基因集名，第二列描述（可留空），第三列起是基因成员。注意分隔符必须是Tab。

library(clusterProfiler) gmt <- read.gmt("h.all.v2023.1.Hs.symbols.gmt") dim(gmt) head(gmt) # 检查与geneList的匹配率 matched <- length(intersect(gmt$gene, names(geneList))) cat("匹配基因数:", matched, "\n")

read.gmt()会把GMT读成两列data.frame：gs_name和gene_symbol。匹配率过低时不要往下跑，先回去查ID类型。常见问题是表达矩阵用ENTREZID、基因集里是symbol，或者反过来。统一ID这步值得花10分钟，比后面排查快得多。匹配率在30%以上可以继续，低于20%基本是ID类型错了。

3. 用R跑多组GSEA并批量出图：clusterProfiler核心代码与参数调法

3.1 GSEA的排序原理：为什么geneList顺序决定一切

GSEA做的事情是：给一个排好序的基因列表，看某个基因集的成员是集中在列表顶部（富集于上调端）还是底部（富集于下调端），并计算富集分数ES。这个分数本质上是加权随机游走的最大偏离值，所以列表排序方向直接影响富集方向的判断。同一份log2FC向量，只要翻转正负号，再看同一个基因集，富集方向就会完全颠倒。

library(clusterProfiler) library(enrichplot) # 用一个比较组跑GSEA gsea_res <- GSEA( geneList, TERM2GENE = gmt, minGSSize = 10, maxGSSize = 500, pvalueCutoff = 1, pAdjustMethod = "BH", seed = 2024, by = "fgsea" ) head(gsea_res@result[, c("ID", "NES", "pvalue", "p.adjust")])

这里minGSSize和maxGSSize控制参与计算的基因集大小范围，防止过小或过大的集合制造假信号。pvalueCutoff=1表示先不过滤任何通路，让所有计算结果完整返回，绘图时再筛。seed=2024固定置换检验的随机种子，保证多人复跑结果一致。by="fgsea"是更快的富集分数实现，大数据集建议保留。

3.2 多组比较的循环结构：用一个列表装下所有GSEA结果

多GSEA的核心是循环。常见做法是把要比较的组队写成一个列表，一个for循环跑完所有比较组，结果存进列表。

comparisons <- list( c("treat", "control"), c("drugA", "control"), c("drugB", "treat") ) gsea_results <- list() for (i in seq_along(comparisons)) { grp1 <- comparisons[[i]][1] grp2 <- comparisons[[i]][2] fc <- log2(rowMeans(expr[, group$sample[group$group == grp1]]) + 0.01) - log2(rowMeans(expr[, group$sample[group$group == grp2]]) + 0.01) gl <- sort(fc, decreasing = TRUE) gsea_results[[i]] <- GSEA( gl, TERM2GENE = gmt, minGSSize = 10, maxGSSize = 500, pvalueCutoff = 1, pAdjustMethod = "BH", seed = 2024, by = "fgsea" ) } names(gsea_results) <- sapply(comparisons, function(x) paste0(x[1], "_vs_", x[2]))

循环里每次都重新算log2FC，原因是不同比较组的样本划分不一样，排序列表也不同。给结果命名时用组1_vs_组2格式，后面拼图时标题直接引用，避免忘了哪个图是哪个比较。如果你已经有一个多组差异分析的总表，也可以按组别拆出来复用里面的logFC。

3.3 用gseaplot2画单条通路：三段式子图的组成

gseaplot2是enrichplot里最常用的出图函数。一张标准GSEA子图包含三块内容：上方的running enrichment score曲线、中部的hit基因分布竖线、下方的ranked list基因位置。审稿人很看重这个三段式结构，不要只用上面一条曲线去充数。

p1 <- gseaplot2( gsea_results[[1]], geneSetID = 1, title = "treat vs control", base_size = 12 ) print(p1)

geneSetID是结果表里的行号，title显示在子图顶部。想选特定通路时，先在结果表里找到行号：

gsea_results[[1]]@result[1:5, c("ID", "NES", "p.adjust")]

我一般把每个比较组里p.adjust最小的通路挑出来，或者按研究主题固定选几条通路，两类图各留一张。不要只看p值，还要看NES绝对值；NES接近0的通路即使p值小，生物学意义也有限。

3.4 多图拼接与批量保存：patchwork的简洁写法

多个比较组拼在一起，用patchwork最省事。/表示上下堆叠，|表示左右并排。拼图后存PDF或高分辨率PNG。

library(patchwork) p2 <- gseaplot2(gsea_results[[2]], geneSetID = 1, title = "drugA vs control", base_size = 12) p3 <- gseaplot2(gsea_results[[3]], geneSetID = 1, title = "drugB vs treat", base_size = 12) p_combined <- p1 / p2 / p3 ggsave( "multi_GSEA_panel.pdf", p_combined, width = 7, height = 10, units = "cm", dpi = 600 ) ggsave( "multi_GSEA_panel.png", p_combined, width = 14, height = 20, units = "cm", dpi = 600 )

PDF用于投稿，矢量图排版无损缩放；PNG用于给合作者快速过目，dpi设为600避免预审时被判定清晰度不足。宽度和高度按子图数量动态调整：每张子图大约7cm宽、10cm高，三张垂直排布就是7×10。

3.5 参数速查：一组常用的GSEA参数值

4. SCI拼图与排版细节：标题、坐标系与分辨率的收敛

4.1 子图标题与比较组命名：让审稿人一眼看清每个图在比什么

拼图后最容易出现的问题是：三张子图只有通路名没有比较组信息，审稿人得回正文猜。gseaplot2的title参数放一行，我的习惯是写成“treat vs control”这种格式。通路名可以放到subtitle，补一行。

p1 <- gseaplot2( gsea_results[[1]], geneSetID = 1, title = "treat vs control", subtitle = "HALLMARK_GLYCOLYSIS", base_size = 12 )

多组图放在论文里时，图注里也要写清楚每个子图对应哪个比较组，主图内标题和正文术语保持一致。不要一个图里写“treat vs control”、另一个图里写“T vs C”，这种不一致在校对阶段非常耗时间。

4.2 字号与坐标轴刻度：base_size统一全文风格

gseaplot2内部是ggplot对象，全局字号由base_size控制，不是单个theme()。想让整篇文章的字号统一，生成每个子图时传同一个base_size。

p1 <- gseaplot2(gsea_results[[1]], geneSetID = 1, title = "treat vs control", base_size = 12) p2 <- gseaplot2(gsea_results[[2]], geneSetID = 1, title = "drugA vs control", base_size = 12) p3 <- gseaplot2(gsea_results[[3]], geneSetID = 1, title = "drugB vs treat", base_size = 12)

y轴刻度范围在拼图时会被patchwork自动对齐，不需要手动设置。如果某张图的running score绝对值比其他图大很多，拼图后会自动压缩其它图的高度，这时用plot_layout(heights = c(1, 1, 1))强制三张图等高。

4.3 多图排版成两列：超过4张图时的版面控制

超过4张图时，垂直堆叠会让页面过长。常见做法是排成两列，每列两到三张图。wrap_plots按列表自动分列。

plot_list <- list(p1, p2, p3, p4) p_final <- wrap_plots(plot_list, ncol = 2) ggsave( "GSEA_two_columns.pdf", p_final, width = 16, height = 18, units = "cm", dpi = 600 )

用ncol=2时，1和2并排，3和4并排。好处是每张子图的宽度不会因为拼图而过窄，running score曲线的细节仍然可辨。4张图以内我通常单列，超过4张再考虑两列。两列时图题和图注的位置也要跟着调整，图片宽度超过16cm的投稿图一般会被期刊要求单独一个文件上传。

5. 多GSEA富集图的5个常见坑：现象、原因与排查顺序

5.1 坑一：跑完结果为空，一条通路都没富集到

现象：GSEA返回的结果表是空的，nrow(result)为0。 原因：最常见的是基因ID类型不匹配。geneList的names是symbol，但TERM2GENE里是ENTREZID，匹配率接近0。其次是minGSSize设得太大，把所有通路都过滤掉了。 解决：先用intersect(names(geneList), gmt$gene)检查匹配基因数量。匹配率低时统一ID类型；匹配率正常但结果仍为空，把minGSSize临时降到5重新跑，确认是不是过滤阈值的问题。

5.2 坑二：富集曲线方向反过来，原本上调的通路显示成下调

现象：同一组对比，用其它工具看是上调，gseaplot2画出来方向完全相反。 原因：给geneList做了abs()，或者用了负值排序，破坏了符号信息；也可能比较组写反，treat_vs_control写成了control_vs_treat。 解决：检查geneList头部基因是否确实是目标上调基因。如果头部是负的大绝对值，就把比较组的前后顺序对调再重跑，不要手动翻转结果。

5.3 坑三：多组拼图后y轴范围不一致导致视觉误导

现象：三张子图拼出来后，running enrichment score曲线的最高点位置参差不齐，看起来像富集强度差异很大。 原因：每个子图独立计算坐标系范围，曲线波形差异是正常的，但y轴范围不一致会放大视觉差。 解决：拼图时保持默认坐标系看整体趋势。想严格对齐y轴，可以在patchwork拼图后加plot_layout(axes = "collect")坐标轴统一收集。少数情况下直接比较NES数值比肉眼比曲线高度更可靠。

5.4 坑四：GMT文件分隔符写错，基因匹配率骤降

现象：用read.gmt()读文件没有报错，但匹配数量比预期低很多。 原因：GMT是Tab分隔，但用Excel打开重存后分隔符变成了空格，或者第二列描述里有空格导致解析错位。 解决：用文本编辑器打开确认GMT原始格式，不要用Excel改GMT。重新下载原始GMT，或者用read.delim(..., sep="\t", header=FALSE)检查前两行，确认每行拆出来的列数符合预期。

5.5 坑五：pvalueCutoff设成0.05导致多组之间缺项

现象：每组单独跑能出图，拼在一起时某个比较组没有对应通路，版面不对称。 原因：每个比较组在pvalueCutoff=0.05过滤后，保留的通路集合不一致，组间P值波动让某些组一条都留不下。 解决：跑GSEA时一律设pvalueCutoff=1，过滤留到绘图阶段。挑通路时按p.adjust值固定阈值，比如每组挑前3条，这样多组图的结构才能对称。

6. 进阶：用NES气泡图总览多组结果，并固定随机种子复现

6.1 提取NES矩阵，一张图看完全部比较组的富集方向

多GSEA富集图之外，还有一个更常用的总览图：用NES做气泡图，横轴放比较组，纵轴放通路，点的大小和颜色表示富集方向与强度。做法是把每个GSEA结果里的NES和p.adjust提取出来。

nes_mat <- do.call(cbind, lapply(gsea_results, function(x) { res <- x@result setNames(res$NES, res$ID) })) padj_mat <- do.call(cbind, lapply(gsea_results, function(x) { res <- x@result setNames(res$p.adjust, res$ID) })) # 转长表后用ggplot画气泡图 library(ggplot2) df <- expand.grid(pathway = rownames(nes_mat), comparison = colnames(nes_mat)) df$NES <- as.vector(nes_mat) df$padj <- as.vector(padj_mat) df$sig <- ifelse(df$padj < 0.05, "yes", "no") ggplot(df, aes(x = comparison, y = pathway, size = abs(NES), color = NES)) + geom_point(aes(alpha = sig)) + scale_color_gradient2(low = "blue", mid = "white", high = "red") + theme_minimal(base_size = 12)

气泡图不需要额外装饰，保留坐标轴标签即可。这张图适合做论文补充图，可以把十几条通路和五六个比较组的结论压缩到一张图里。

6.2 固定随机种子，让结果可复现

GSEA的置换检验本质是随机过程，不固定种子时同一份数据两次运行的结果在小数位上有差异。这对提交论文是个隐患，审稿人复跑对不上数时会质疑分析可靠性。seed=2024这一行就是给整个流程上保险。

我的习惯是所有R分析脚本第一行就写set.seed(2024)，不只GSEA，凡是涉及重采样的分析全部固定。这件事我是吃过亏才记住的——某次合作者复跑我的富集分析，P值在小数点后第二位对不上，折腾了两天才发现是随机种子没固定。现在我把这点写进模板，后面所有项目照此执行，希望帮到你。

本文还有配套的精品资源，点击获取