企业官网建设流程全解析

1. 前言：一次基因注释任务的复盘

这两年单细胞测序和RNA-seq项目越来越多，手里攒着几百上千个差异基因做功能注释，几乎是每天都要面对的场景。很多人手里拿到一个基因列表，第一步就是往网页工具里一贴，点个提交，然后拿着产出的GO term和KEGG通路列表开始写文章。这样做本身没什么问题，但当你手里的基因数量从几十个变成几千个，或者项目里同时有多个分组的基因列表需要横向对比，在线工具那套玩法的效率就明显跟不上了。

这篇文章要聊的，就是对基因列表进行批量GO和KEGG注释这件事。我会从思路拆解、工具选型、实操流程、常见坑位这几个角度，把整个环节从头到尾过一遍，给出一个可复现、可落地的参考方案。无论你是刚开始接触组学数据分析，还是已经有一定基础、想优化现有流程，这篇文章都应该能提供一些有价值的东西。

2. 批量注释的整体思路与方案选型

2.1 为什么在线工具做批量注释不合适

很多人喜欢先走DAVID、Enrichr这些在线平台，图的是不用装软件、界面友好。但在实际项目中，在线注释有几道迈不过去的坎：

• 一次最多处理几百个基因，上千个基因就得拆成好几批；
• 每次只能处理一个列表，多个分组的列表比对时要手动导出再合并；
• 注释的版本和基因ID类型经常变化，结果无法保证完全可重现；
• 上传的基因列表涉及未发表项目的核心数据，放在第三方服务器上本身就是风险。

所以当项目里真正出现“批量”两个字，离线脚本基本就是唯一合理的选择。

2.2 主流工具链的取舍逻辑

目前做GO和KEGG注释的主流工具链大概有三条路线：

我个人最常用的是R的clusterProfiler包，原因很简单：它的输入输出设计得很统一，输入一个基因向量，输出标准的富集结果对象，不管后续是做可视化、合并多组结果、还是导出成Excel/CSV，都有现成的函数可以调用。而且它在处理物种注释时需要加载对应的OrgDb包，这意味着你用一个物种的时候，整个注释体系自洽性比较一致。

2.3 先搞清楚ID类型再动手

这一步听起来基础，但大多数注释流程跑不动、结果为空，根子上的问题都出在这里。

在做GO和KEGG富集分析时，工具本身是识别基因ID的，而基因列表里的ID类型五花八门：常见的有Symbol、Entrez ID、Ensembl ID、RefSeq ID等。clusterProfiler内部的主键是Entrez ID，这不代表你得手动把所有Symbol都转成Entrez，但你需要让它知道输入的是哪种类型。

实际操作时要先做一次ID映射检查：

• 如果输入是Entrez ID，直接跑就行；
• 如果是Symbol，需要借助bitr()函数转换，同时要留意存在部分基因转换不到的情况，这些基因的比例如果太高，结果就不太可信；
• 如果是Ensembl ID，转换时会遇到前缀问题（比如人源的ENSG开头、鼠源的ENSMUSG开头），需要注意提前处理。

我曾经处理过一批某物种的RNA-seq数据，最开始直接用Symbol做注释，结果通路富集出来的结果几乎全和线粒体相关的通路相关，后来排查才发现输入的Symbol和参考数据库的ID格式不一致，导致大量基因没有匹配上，剩下来的少量基因恰好富集在一些特殊通路上。从那以后，我开始习惯性先统计ID转换的成功率，低于85%就会回头检查输入格式，而不是直接跑结果。

3. 核心环节：从数据准备到富集结果

3.1 环境准备与依赖安装

如果你从零开始搭环境，建议使用R 4.x版本，用BiocManager来安装生物信息学相关的包。基础依赖装起来通常不会出大问题，但有几个包的安装可能需要一点耐心：

if (!require("BiocManager", quietly = TRUE)) install.packages("BiocManager") BiocManager::install("clusterProfiler") BiocManager::install("org.Hs.eg.db") BiocManager::install("org.Mm.eg.db") BiocManager::install("pathview")

其中org.Hs.eg.db对应人类，org.Mm.eg.db对应小鼠。如果做的是其他物种，可以用AnnotationHub来检索对应的注释包，比如斑马鱼、果蝇、拟南芥之类都有对应的OrgDb。

有一个点需要提前提醒：安装这些包时R的版本要和BioConductor的版本匹配，否则容易出现依赖不兼容。遇到这种情况不要硬碰硬去手动装依赖，先检查一下BiocManager::version()输出的版本号，再对照BioConductor官网给出的对应版本要求。

3.2 读取基因列表与ID转换实操

假设你手里有一个差异基因列表文件，第一列是基因名，页面长这样：

gene_id TP53 BRCA1 EGFR MYC ...

读入R并进行ID转换：

library(clusterProfiler) library(org.Hs.eg.db) # 读取基因列表 deg <- read.csv("deg_list.csv", stringsAsFactors = FALSE) genes <- deg$gene_id # Symbol转Entrez entrez_ids <- bitr(genes, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = org.Hs.eg.db) # 查看转换效率 cat("原始基因数：", length(genes), "\n") cat("成功转换数：", nrow(entrez_ids), "\n")

bitr()函数的输出可能包含一对多的情况，即一个Symbol对应多个Entrez ID，或者反过来。这时候默认参数会保留所有对应关系，如果你在意多对一映射的干扰，可以加上unique()去重处理，但这个操作通常影响不大，因为大部分基因是唯一对应的。

从实践来看，转换失败的那部分基因，通常问题集中在：

• Symbol里带了版本号或者额外后缀（比如TP53.1这种）；
• 基因名是旧版本命名，无法匹配当前注释库；
• 基因来自非模式物种，OrgDb里的注释信息不全。

3.3 GO富集分析详解

GO分析的核心逻辑是把输入的基因集合放到GO层次结构里，看哪些功能分类显著富集。GO分成三个子本体：生物学过程（BP）、细胞组分（CC）、分子功能（MF）。默认情况下clusterProfiler会同时跑三个，也可以指定某个子本体单独跑。

# GO富集分析 go_enrich <- enrichGO(gene = entrez_ids$ENTREZID, OrgDb = org.Hs.eg.db, keyType = "ENTREZID", ont = "ALL", # BP/CC/MF或ALL pAdjustMethod = "BH", pvalueCutoff = 0.05, qvalueCutoff = 0.2, readable = TRUE) # 直接将Entrez映射回Symbol # 查看结果 head(as.data.frame(go_enrich))

几个关键参数值得展开说明：

pAdjustMethod选“BH”是Benjamini-Hochberg校正，也就是FDR控制。做多组学项目时，很少有人会用未校正的p值，因为基因数一多，假阳性会非常可观。

qvalueCutoff设为0.2，这个阈值看起来宽松，但富集分析里q值0.05的门槛会让很多生物学上真实存在的弱信号丢失。当然具体阈值得根据自己的项目要求来。

readable = TRUE会让结果里显示的基因从Entrez ID变成容易被读懂的Symbol。需要注意，设置了这个参数后，结果数据框里会出现“geneID”列，这一列是富集到该分类的所有基因用/分隔组成的长字符串，后续如果要统计单个基因的富集情况，需要拆开处理。

GO结果输出的结构是这样的：

GeneRatio是富集到该GO term的输入基因数除以所有输入到分析中的基因数，BgRatio是注释到该term的背景基因数除以全基因组背景基因数。这两个比例一对比，富集倍数这个概念就立起来了，连带着富集因子（fold enrichment）也就出来了。

实操中我很喜欢用enrichplot包里的几个画图函数来快速把握富集结果全貌：

library(enrichplot) # 气泡图 dotplot(go_enrich, showCategory = 15) # 条形图 barplot(go_enrich, showCategory = 15)

dotplot是我用得最多的一个，横轴是GeneRatio，纵轴是GO term，点的大小代表富集的基因数，颜色代表p.adjust值。一张图就能看出哪些term富集程度高、哪些term富集的基因多，做PPT和组会汇报时特别有用。

3.4 KEGG富集分析详解

KEGG思路和GO不太一样，它更侧重代谢和信号通路层面的功能富集。clusterProfiler处理KEGG时，不再使用OrgDb作为注释来源，而是通过在线KEGG API获取最新通路注释数据。

# KEGG富集分析 kegg_enrich <- enrichKEGG(gene = entrez_ids$ENTREZID, organism = "hsa", # 人类；小鼠是mmu keyType = "kegg", pAdjustMethod = "BH", pvalueCutoff = 0.05, qvalueCutoff = 0.2) # 查看结果 head(as.data.frame(kegg_enrich))

organism参数是三字母的物种缩写，人类是hsa、小鼠是mmu、大鼠是rno、斑马鱼是dre。如果你是做其他物种的，需要去KEGG官网查相应缩写，同时在BRITE数据库里确认这个物种有通路注释数据。很多小众物种KEGG注释不完整，强行跑出来结果会很稀疏。

有一个我反复遇到的坑：KEGG数据库更新后，有些学术申请模式以外的访问方式会超时。enrichKEGG()底层在线访问KEGG API，如果运行时出现网络超时，可以尝试用download.kegg = FALSE参数，让它在本地缓存了的部分结果上分析。但如果网络环境确实有问题，最好还是用KEGG官方离线数据库文件配合enrichKEGG()的use_internal_data = TRUE参数来跑。

KEGG结果里除了标准的富集条目，还有一个值得关注的信息——geneID列里每个基因是否映射到了KEGG里的KO（KEGG Orthology）条目。后续如果你想画通路图，比如用pathview把差异基因标到特定通路上，就需要这一层映射关系。实操中，我喜欢在KEGG富集结果里先挑感兴趣的显著通路，然后单独跑一下pathview，把表达变化信息直接上色到通路图上，那种图拿去写文章真的非常加分。

3.5 从富集结果到可视化输出

分析完的富集结果是要拿去写报告、做PPT的，不可能每次都重新去跑R脚本。把结果导出成表格，保存好图表，我一般是这么处理的：

# 导出GO结果 write.csv(as.data.frame(go_enrich), "GO_enrich_results.csv", row.names = FALSE) # 导出KEGG结果 write.csv(as.data.frame(kegg_enrich), "KEGG_enrich_results.csv", row.names = FALSE) # 保存R对象，方便后续重新出图 saveRDS(go_enrich, "go_enrich.rds") saveRDS(kegg_enrich, "kegg_enrich.rds")

另外，同一个项目的多个分组基因列表，我习惯用一个循环去批量跑富集，把结果汇总到一个列表对象里。这样后续可以做组间的GO term对比，比如多个分组之间哪些生物学过程是共同富集的，哪些是各自特有的。

这类多列表对比场景，clusterProfiler提供的compareCluster()函数相当顺手：

# 构建一个包含多组基因的list gene_list <- list( group1 = group1_entrez, group2 = group2_entrez, group3 = group3_entrez ) # 批量比较富集 comp_go <- compareCluster(geneClusters = gene_list, fun = "enrichGO", OrgDb = org.Hs.eg.db, ont = "BP", pvalueCutoff = 0.05) # 绘制分组气泡图 dotplot(comp_go, showCategory = 10)

4. 实操过程：一个完整案例跑通全流程

4.1 案例背景与数据准备

我在某次项目中拿到了一套某物种的转录组差异基因数据，三个比较组各有一份基因列表，每组基因数大约在800到1500之间。需求是每个组分别做GO和KEGG富集，最后汇总成一个Excel表格作为项目报告的补充材料，再者挑选显著通路做通路图展示。

为了演示，这里用一个小规模模拟数据集来说明完整流程。假设目前有这样一个基因列表文件deg_list.csv：

gene_id,log2FoldChange TP53,2.5 BRCA1,1.8 EGFR,-1.3 MYC,3.1 ...

第一步仍然是读入数据，把Symbol转换为Entrez ID。这一步要做成功率的统计，必要时向项目组反馈有多少基因被过滤掉了。

4.2 一个循环跑完所有分组

实际项目通常是多个分组一起处理，写一个循环把所有分组的数据都跑出来比较省事。我这个案例里有三个分组，流程写成函数：

run_enrichment <- function(gene_symbols, orgdb, organism_kegg) { converted <- bitr(gene_symbols, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = orgdb) go_res <- enrichGO(gene = converted$ENTREZID, OrgDb = orgdb, keyType = "ENTREZID", ont = "ALL", pAdjustMethod = "BH", pvalueCutoff = 0.05, qvalueCutoff = 0.2, readable = TRUE) kegg_res <- enrichKEGG(gene = converted$ENTREZID, organism = organism_kegg, pAdjustMethod = "BH", pvalueCutoff = 0.05, qvalueCutoff = 0.2) list(go = go_res, kegg = kegg_res) } # 假设三个分组基因列表分别存储 results <- list( group1 = run_enrichment(group1_genes, org.Hs.eg.db, "hsa"), group2 = run_enrichment(group2_genes, org.Hs.eg.db, "hsa"), group3 = run_enrichment(group3_genes, org.Hs.eg.db, "hsa") )

这个函数可以复用，下次换物种就换OrgDb和KEGG缩写，换数据就换输入基因列表，结构不需要再改。

4.3 结果汇总与表格整理

把每个组的GO和KEGG结果按统一格式整理到一张表里，方便报告使用：

library(dplyr) go_summary <- bind_rows( lapply(names(results), function(grp) { df <- as.data.frame(results[[grp]]$go) df$group <- grp df }) ) kegg_summary <- bind_rows( lapply(names(results), function(grp) { df <- as.data.frame(results[[grp]]$kegg) df$group <- grp df }) ) write.csv(go_summary, "all_groups_GO_results.csv", row.names = FALSE) write.csv(kegg_summary, "all_groups_KEGG_results.csv", row.names = FALSE)

整理后的表格字段非常直白，每一列都有明确含义，非生信背景的合作者也能看懂哪一列代表什么。我在给实验室其他成员交付结果时，会额外生成一个单独列来说明这个通路的生物学意义，配合显著富集的基因列表，这样别人拿过去就能直接用。

4.4 可视化实操与细节处理

气泡图是最常用的富集可视化方式。如果直接用默认参数，有时候出的图在展示给非生信背景的人时会显得拥挤。我的习惯是先筛选出显著度较高的通路，再控制展示数量。

# 挑选每组GO结果中p.adjust最小的前10个term top_terms <- go_summary %>% filter(p.adjust < 0.05) %>% group_by(group, ONTOLOGY) %>% arrange(p.adjust) %>% slice_head(n = 10) %>% ungroup() # 换个方式：直接用compareCluster结果出图 dotplot(comp_go, showCategory = 10, split = "ONTOLOGY") + facet_grid(ONTOLOGY ~ ., scale = "free")

注意，如果用了split = "ONTOLOGY"加上facet_grid，GO的BP、CC、MF会分开展示，图面会更清楚。这也是我在写文章配图时常用的方式，因为审稿人通常很在意三个本体是否混在一起导致读数困难。

KEGG通路图用pathview来画：

library(pathview) # 以hsa04110（细胞周期）为例 pathview(gene.data = deg$log2FoldChange, gene.idtype = "SYMBOL", pathway.id = "04110", species = "hsa", out.suffix = "cell_cycle")

pathview会输出一个hsa04110.cell_cycle.png文件，通路图上的基因会被标记上红色和绿色来代表上调和下调。这张图放进文章配图或汇报材料里，直观程度远胜一大段文字描述。

4.5 背景基因集的选择问题

富集分析的统计学基础是超几何分布检验，它比对的背景集是什么样的，直接决定结果差异。clusterProfiler默认kegg的 background 是整个基因组注释基因，GO也是用OrgDb中注释的全部基因。但有的时候项目里自己设定了表达基因背景，比如只把表达量高于阈值的基因作为背景，这样会引入“选择偏差”，让富集结果更符合项目的实际数据情况。

实际操作中，对于差异基因富集分析，大多数文献里仍然使用全基因组作为背景，这是默认且较稳妥的做法。但是如果做的是单细胞数据里特定细胞类型的标志基因注释，可以考虑用该细胞类型中检测到的所有基因作为背景。这时候clusterProfiler提供了universe参数，可以传入自己的背景基因集合，这个参数不填就是全基因组，填了就是自定义背景。

那什么时候自定义背景有意义呢？我举一个场景：某项目中我只分析了表达量排前25%的高表达基因，如果背景还是全基因组，这些高表达基因的富集结果会全部偏向核糖体、氧化磷酸化这类高表达基因集中的通路，看起来富集程度很高但缺乏生物学新意。这时候把背景改成所有表达的基因，富集结果才是真正在这群高表达基因中相对显著的通路，才是这个项目真正想找的东西。

这个细节，很多在线工具里根本没法调，也是离线流程的一个关键优势。

5. 常见问题与排查技巧实录

5.1 ID转换失败率过高

现象：转换成功的基因数量不足输入基因总数的70%。

排查思路：

• 先看前几个失败的基因长什么样。如果是“TP53.1”这种带小数点后缀的，需要在转换前先做字符串清洗；
• 如果输入的是Ensembl ID（如ENSG开头），要确认是否带了版本号后缀，比如ENSG00000141510.12，这种后缀在转换时也会导致匹配失败；
• 检查自己用的OrgDb物种和基因列表的物种是否一致。人源基因列表用了小鼠OrgDb，转换率会惨不忍睹。

推荐处理办法：在bitr转换后，把成功的基因列表和失败基因列表分别往外输出一份，保留审计痕迹。这样哪怕后来项目被问起来，也能说清楚哪些基因被排除在富集分析以外了。

5.2 富集结果一个显著条目都没有

现象：跑出来的富集结果空空的，或者只有几个q值很大的条目。

原因通常有三个：

• 基因数量太少（少于30个），统计功效不足。这种情况可以考虑放宽pvalueCutoff或者换用秩次类富集方法；
• 基因列表里大部分基因是同一个通路家族的，比如全部是核糖体蛋白基因，那GO富集确实会非常受限；
• 背景基因集设置不恰当，导致所有条目都不显著。

我曾遇到过一个极端案例：输入的是一个筛选后极其严格的基因集，只有不到20个基因，GO跑出来一个显著项都没有，但这批基因在表达层面差异非常显著。后来我用GSEA这类基于基因排序的方法补充分析，立刻就富集出清晰的通路信号来。这个经验提醒我，富集分析的统计显著性和生物学显著性完全是两码事，不能仅仅因为“没有显著富集”就否定这批基因的功能相关性。

5.3 enrichKEGG网络请求失败

KEGG富集经常遇到网络问题，症状是脚本运行卡住后报错，提示无法连接到KEGG REST API。常见解决办法：

# 设置下载超时时间 options(timeout = 300) # 或者直接尝试用本地的KEGG数据 kegg_enrich <- enrichKEGG(gene = entrez_ids$ENTREZID, organism = "hsa", use_internal_data = TRUE)

use_internal_data = TRUE需要提前通过downloadKEGG()下载过对应物种的KEGG数据并缓存。如果网络环境确实不理想，也可以考虑用clusterProfiler的替代方案，直接把enrichKEGG换成enrichPathway函数，后者使用Reactome数据库，网络依赖相对小一些。

5.4 不同版本注释库导致结果不可重现

生物信息学工具和数据库版本更新极快，同一个基因列表，半年前跑和现在跑，富集结果可能会有差异。为了保证项目的可重现性，强烈建议在项目记录里写明这些信息：

• R版本和clusterProfiler版本
• OrgDb版本
• KEGG数据库的release版本
• 输入的基因ID类型和转换时间

我习惯把所有分析脚本、输入文件和环境版本信息打包存一个文件夹，标注好分析日期。这样无论是后续补充实验数据，还是应对审稿时补充材料的需求，都不会手忙脚乱。

5.5 问题排查速查表

6. 一些高级扩展思路

6.1 不只是富集：基因集评分

富集分析解决的是“这些基因是否在某类功能上显著富集”的问题，但在很多项目里，我们真正想回答的是“某些功能通路的活性在样本间是怎么变化的”。这个问题就要用基因集评分来回答。

常见做法有GSVA和ssGSEA。GSVA把基因表达矩阵转换成通路活性矩阵，后续可以像看基因表达一样做差异分析；ssGSEA则给每个样本算一个通路活性得分。这两个方法可以和富集分析互为补充，而且都能用clusterProfiler周边的包实现。

6.2 多物种比较分析

如果你手里的项目涉及跨物种比较，比如比较人和小鼠的同一生物学过程，KEGG的organism参数需要分别指定，之后把两边的富集结果按通路ID对齐再做交集或差异比较。要注意不同物种之间通路ID虽然是共享的KEGG通路ID，但每个通路对应的基因集合不同，所以“同一个通路富集到了”不代表“富集的是同一群基因”。

6.3 自动化报告生成

分析做完以后，把结果整合成一份自动化报告是很有价值的。R Markdown配合参数化报告能让整条流程彻底摆脱手工复制粘贴。每批数据来了，只要更新输入文件，点一下knit，就能输出一份新的报告，包括背景信息、统计结果、图表、结论。团队有多个项目并行时，这个流程节约的时间非常可观。

7. 实操中我总结的几个心得

GO和KEGG注释是整个组学数据分析流程里最容易入门也最容易做错的一环。它的原理不难，统计模型也不复杂，但真正把结果做得能进文章、能过审稿、能经得起复现，细节非常多。

几点个人经验放在最后：

第一，别盲目迷信默认参数。pvalueCutoff = 0.05不是金科玉律，要根据基因数量和数据质量合理调整。如果基因列表只有50个基因，硬拿0.05的阈值跑，结果大概率不显著，这不能说明你的基因没有功能，只能说明统计功效不够。

第二，富集分析前先画个简单的表达分布图。把差异基因的上下调情况、表达倍数分布先看清楚，有时候一眼就能发现基因列表可能有批次效应或者某个通路占主导。

第三，每个项目的注释版本信息一定要留存。这不是可做可不做的事，而是对有可重复性要求的项目来说必须做的事。短期内也许没人问，但如果文章进入审稿流程，审稿人要求补充材料或者复核数据，版本信息就是你唯一的救命稻草。

第四，在线工具与本地流程并不冲突。前期探索性分析用在线工具快速出图可以，但是项目正式分析一定要切到本地流程。给合作者交付结果时，我不会只给一个图片，而是同时交付表格、脚本、环境配置记录，这样别人拿到手里的是一整套完整产出，而不仅是结果图。

组学数据分析这一行里，真正拉开差距的不是跑通了多复杂的算法，而是能不能稳定、高效、可重复地处理常规任务。批量GO和KEGG注释就是这种最常规的任务之一，希望这篇文章能帮你把这条链路扎扎实实搭建起来。