从R的clusterProfiler到Python的gseapy：手把手教你做基因富集分析（ORA）及结果解读避坑指南

weixin_30402343

298人浏览 · 2026-06-04 11:10:58

weixin_30402343 · 2026-06-04 11:10:58 发布

从R的clusterProfiler到Python的gseapy：基因富集分析实战与科学解读指南

在生物信息学研究中，基因富集分析（Gene Set Enrichment Analysis）已成为解读高通量组学数据的标准流程。其中，过表达分析（Over-Representation Analysis, ORA）因其直观性和易解释性，成为最常用的方法之一。本文将聚焦两大主流工具——R语言的clusterProfiler和Python的gseapy，通过对比演示帮助研究者快速掌握从数据准备到结果解读的全流程。

1. 环境准备与数据加载

1.1 安装必要工具包

对于R用户，clusterProfiler是Bioconductor生态系统的重要组成部分：

if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("clusterProfiler")

Python用户则可以通过pip安装gseapy：

pip install gseapy

1.2 输入数据规范

两种工具都接受基因列表作为输入，但格式要求略有不同：

工具	输入格式要求	示例数据结构
clusterProfiler	字符型向量或数据框	`c("GeneA", "GeneB")`
gseapy	列表或pandas Series	`["GeneA", "GeneB"]`

关键点 ：确保基因标识符（如Symbol或Entrez ID）与注释数据库保持一致。常见的标识符问题会导致30%的分析失败案例。

2. 核心分析流程对比

2.1 R/clusterProfiler完整工作流

以下是使用GO数据库进行富集分析的典型代码：

library(clusterProfiler)

# 准备基因列表
gene_list <- c("TP53", "BRCA1", "CDK2", "MYC", "EGFR")

# 执行富集分析
ego <- enrichGO(gene = gene_list,
                OrgDb = "org.Hs.eg.db",
                keyType = "SYMBOL",
                ont = "BP",
                pvalueCutoff = 0.05,
                qvalueCutoff = 0.2)

# 结果可视化
dotplot(ego, showCategory=15)

2.2 Python/gseapy实现方案

gseapy提供了更灵活的数据库选择：

import gseapy as gp

# 基因列表准备
gene_list = ["TP53", "BRCA1", "CDK2", "MYC", "EGFR"]

# 执行分析
enrich_results = gp.enrichr(gene_list=gene_list,
                           gene_sets=['GO_Biological_Process_2021'],
                           organism='human')

# 结果导出
enrich_results.results.to_csv("enrichment_results.csv", index=False)

注意：gseapy默认使用Enrichr在线数据库，当处理敏感数据时，建议通过 no_plot=True 参数禁用自动联网绘图功能。

3. 关键参数解析与优化

3.1 统计检验方法选择

两种工具都支持多种校正方法：

校正方法	clusterProfiler参数	gseapy参数
Benjamini-Hochberg	`pAdjustMethod = "BH"`	`adjust_method='fdr'`
Bonferroni	`pAdjustMethod = "bonferroni"`	`adjust_method='bonferroni'`
不校正	`pAdjustMethod = "none"`	`adjust_method=None`

经验建议 ：对于初步筛选，推荐使用FDR校正；当需要严格控制假阳性时，可考虑Bonferroni方法。

3.2 背景基因集设定

背景基因的选择直接影响结果可靠性：

clusterProfiler ：通过 universe 参数指定
gseapy ：通过 background 参数指定

# R中设置背景基因
background_genes <- read.csv("all_expressed_genes.csv")$Gene
ego <- enrichGO(gene = gene_list, 
                universe = background_genes,
                ...)

常见误区 ：使用全基因组作为背景虽然保守，但会降低检测灵敏度。最佳实践是使用实验检测到的所有基因作为背景。

4. 结果解读与可视化

4.1 核心结果字段解析

典型富集分析结果包含以下关键指标：

P-value ：原始显著性水平
Adjusted P-value ：多重检验校正后的显著性
Odds Ratio ：富集程度量化指标
Gene Ratio ：目标基因集中显著基因比例
Count ：显著基因数量

4.2 可视化技巧对比

clusterProfiler 内置丰富的绘图函数：

# 气泡图
dotplot(ego, showCategory=20)

# 网络图
cnetplot(ego, categorySize="pvalue", foldChange=gene_fc)

gseapy 则依赖matplotlib进行自定义：

import matplotlib.pyplot as plt

top_terms = enrich_results.results.head(10)
plt.figure(figsize=(10,6))
plt.barh(top_terms['Term'], -np.log10(top_terms['Adjusted P-value']))
plt.xlabel('-log10(Adj.P-value)')
plt.title('Top Enriched Terms')

专业提示：当展示多个比较组时，考虑使用热图或火山图展示跨条件的富集模式。

5. 实战中的陷阱与解决方案

5.1 基因标识符转换问题

约40%的分析错误源于基因ID不一致。推荐使用专业转换工具：

R ： clusterProfiler::bitr
Python ： mygene 包

from mygene import MyGeneInfo
mg = MyGeneInfo()
result = mg.querymany(["TP53", "BRCA1"], scopes="symbol", fields="entrezgene", species="human")

5.2 数据库版本控制

不同版本的GO/KEGG数据库可能导致结果差异。建议：

记录分析时使用的数据库版本
对于长期项目，固定数据库版本
使用 sessionInfo() (R)或 gseapy.__version__ (Python)保存环境信息

5.3 生物学解释的层次性

避免简单地根据p值排序选择结果。有效的解读策略应包括：

检查富集项间的层级关系
识别核心生物学主题（如通过词云分析）
结合实验假设进行针对性筛选

6. 高级应用场景

6.1 跨物种分析策略

当研究非模式生物时：

# 使用AnnotationHub获取注释
library(AnnotationHub)
ah <- AnnotationHub()
query(ah, c("OrgDb", "Arabidopsis"))

6.2 时间序列数据分析

对于多时间点数据，可进行动态富集分析：

# 使用gseapy分析时间序列
time_points = ["0h", "6h", "12h"]
for tp in time_points:
    genes = get_deg_list(tp)
    gp.enrichr(gene_list=genes, 
              gene_sets=['KEGG_2021_Human'],
              outdir=f"results/{tp}")

6.3 自定义基因集分析

创建个性化基因集：

# 构建GMT格式文件
custom_gmt <- rbind(
  c("Pathway1", "Description", "GeneA\tGeneB\tGeneC"),
  c("Pathway2", "Description", "GeneD\tGeneE")
)
writeLines(apply(custom_gmt, 1, paste, collapse="\t"), "custom.gmt")

7. 性能优化技巧

7.1 并行计算加速

对于大规模分析：

R 方案：

library(BiocParallel)
register(SnowParam(4))  # 使用4个核心
ego <- enrichGO(..., BPPARAM = SnowParam())

Python 方案：

from multiprocessing import Pool

def run_enrich(gene_set):
    return gp.enrichr(gene_list=gene_set, ...)

with Pool(4) as p:
    results = p.map(run_enrich, gene_sets_list)