GO与KEGG富集分析实战：从差异基因到功能注释
1. 从差异基因到功能注释：为什么我们需要GO和KEGG？
如果你刚做完差异表达分析，手里拿到了一长串几百甚至上千个“差异基因”的列表，是不是有点懵？心里肯定在想：这些基因是干嘛的？它们凑在一起，到底在细胞里搞什么“大项目”？这感觉就像你拿到了一份零件清单，但完全不知道它们能组装成什么机器。这时候，GO（Gene Ontology，基因本体论） 和 KEGG（Kyoto Encyclopedia of Genes and Genomes，京都基因与基因组百科全书） 富集分析，就是你手里的“说明书”和“装配图”。
我刚开始做生信分析那会儿，也卡在这里很久。只知道运行代码，出来的图花花绿绿，但根本看不懂背后的生物学意义，更别提给导师或者合作者讲明白了。后来踩过不少坑才明白，富集分析的核心，不是跑通代码，而是理解你的基因列表在生物学功能上的“集体倾向性”。
简单来说，GO就像一个给基因功能分门别类的“图书馆目录系统”。它从三个维度描述基因：
- 细胞组分（Cellular Component, CC）：这个基因的产物（通常是蛋白质）在细胞的哪个位置干活？是在线粒体里提供能量，还是在细胞核里调控DNA？
- 分子功能（Molecular Function, MF）：这个基因产物具体能做什么？是能结合DNA，还是具有酶的催化活性，或者是负责运输物质？
- 生物学过程（Biological Process, BP）：这个基因参与了一个什么样的“连续剧”？比如“细胞周期调控”、“炎症反应”或“葡萄糖代谢”。
而KEGG呢，更像是一张张描绘具体“生产线”或“信号通路”的地图。它告诉你，这些差异基因共同参与了哪些已知的、完整的生物学通路。比如“p53信号通路”、“MAPK信号通路”或者“癌症中的转录调控异常”。如果说GO告诉你零件是“螺丝刀”（功能）和“在工具箱里”（位置），那么KEGG就告诉你这些零件是用来组装“一台自行车”（通路）的。
所以，做富集分析的根本目的，是把冷冰冰的基因ID列表，转化为有血有肉的生物学故事。你的差异基因如果显著富集在“细胞外基质组织”和“TGF-β信号通路”上，那很可能暗示你的实验条件影响了纤维化或组织重塑过程。这个从数据到生物学洞察的飞跃，正是富集分析的价值所在。接下来，我就手把手带你走一遍完整的实战流程，并分享一些我总结出来的、能让结果更可靠的技巧和避坑指南。
2. 实战第一步：数据准备与基因ID转换
拿到差异基因列表后，千万别急着往富集分析工具里扔。第一步的预处理，直接决定了后续分析的成功率。这里最常见的坑就是基因ID格式不匹配。
2.1 理解你的输入：基因标识符的“黑话”
差异分析的结果，基因名可能是各种格式：Gene Symbol（如TP53, ACTB）、Ensembl ID（如ENSG00000141510）、Entrez ID（如7157），或者RefSeq ID（如NM_000546）。大多数富集分析工具，尤其是R包clusterProfiler，其核心计算需要的是Entrez ID（一种数字型的稳定ID）。而GO注释数据库本身则关联着各种ID。
我个人的经验是，始终从最稳定的标识符开始转换。Entrez ID或Ensembl ID通常比Gene Symbol更稳定，因为后者可能存在同名（一个Symbol对应多个基因）或更新（旧的Symbol被废止）的问题。在原始文章复现时，经常因为Symbol对不上而导致富集结果为空或很奇怪，第一步就要检查这里。
2.2 手把手进行ID转换
这里我们用R语言和clusterProfiler系列包来操作，这是目前最主流、最强大的工具集。假设你有一个差异基因列表，基因名是Gene Symbol，保存在DEG_list.txt文件中。
# 加载必要的R包 library(clusterProfiler) library(org.Hs.eg.db) # 如果你是人源数据，这里是人的数据库。小鼠用 org.Mm.eg.db library(dplyr) library(readr) # 1. 读入差异基因列表（假设文件只有一列基因名，没有表头） diff_genes <- read.table("DEG_list.txt", header = FALSE, stringsAsFactors = FALSE)[,1] head(diff_genes) # 查看前几个基因，确认是Symbol，例如 "TP53", "BRCA1" # 2. 进行ID转换：从 Symbol 到 Entrez ID # 使用 bitr 函数，这是 clusterProfiler 里的转换神器 gene_df <- bitr(geneID = diff_genes, fromType = "SYMBOL", # 你输入的ID类型 toType = c("ENTREZID", "ENSEMBL"), # 你想转换成的ID类型，至少需要ENTREZID OrgDb = org.Hs.eg.db) # 指定物种注释数据库 # 查看转换结果 head(gene_df) dim(gene_df) # 看看转换成功了多少个基因 运行后，gene_df这个数据框里就有每个基因对应的SYMBOL、ENTREZID和ENSEMBL了。关键点来了：bi
