企业官网建设流程全解析

用R语言打造高阶KEGG富集桑吉气泡图：从数据到可视化的完整指南

在基因功能富集分析领域，气泡图一直是展示KEGG通路结果的经典选择。但传统气泡图只能呈现通路名称、富集倍数、p值和基因数四个维度，关键的基因列表信息往往被埋没在表格中。本文将带你用R语言的clusterProfiler和ggplot2生态系统，实现基因列表与富集结果的五维联动可视化——这就是桑吉气泡图的魅力所在。

1. 环境准备与数据加载

1.1 安装必要R包

首先确保你的R环境（建议4.0+版本）已安装以下关键包：

if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager") BiocManager::install(c("clusterProfiler", "org.Hs.eg.db", "DOSE")) install.packages(c("ggplot2", "ggrepel", "stringr", "tidyr", "dplyr"))

1.2 加载示例数据集

我们使用clusterProfiler内置的示例数据模拟真实分析场景：

library(clusterProfiler) data(geneList, package = "DOSE") gene <- names(geneList)[abs(geneList) > 2] kk <- enrichKEGG(gene = gene, organism = "hsa") kegg_result <- as.data.frame(kk)

典型KEGG富集结果包含以下关键列：

2. 数据预处理：为桑吉图做准备

2.1 基因列表结构化处理

原始geneID列以"/"分隔基因，需要转换为适合绘图的长格式：

library(tidyr) library(dplyr) kegg_long <- kegg_result %>% separate_rows(geneID, sep = "/") %>% select(Description, GeneRatio, pvalue, Count, geneID)

2.2 添加美学映射参数

为后续绘图添加必要的计算列：

kegg_long <- kegg_long %>% mutate( logP = -log10(pvalue), GeneRatio_num = sapply(strsplit(GeneRatio, "/"), function(x) as.numeric(x[1])/as.numeric(x[2])) )

3. 构建基础气泡图框架

3.1 初始化ggplot对象

library(ggplot2) base_plot <- ggplot(kegg_result, aes(x = GeneRatio_num, y = reorder(Description, GeneRatio_num))) + geom_point(aes(size = Count, color = logP)) + scale_color_gradient(low = "blue", high = "red", name = "-log10(p-value)") + scale_size_continuous(range = c(3, 8), name = "Gene count") + labs(x = "Gene Ratio", y = "") + theme_minimal(base_size = 12)

3.2 优化视觉呈现

添加专业期刊级别的格式调整：

enhanced_plot <- base_plot + theme( panel.grid.major = element_line(color = "grey90"), panel.grid.minor = element_blank(), axis.text.y = element_text(color = "black", size = 10), legend.position = "right", legend.box = "vertical" ) + guides( color = guide_colorbar(barwidth = 1, barheight = 10), size = guide_legend(nrow = 3) )

4. 集成桑吉图元素实现五维展示

4.1 添加基因连线

使用geom_segment创建通路与基因的关联网络：

library(ggrepel) sankey_plot <- enhanced_plot + geom_segment( data = kegg_long, aes(x = 0, xend = -0.05, y = Description, yend = geneID), color = "grey70", linewidth = 0.3 ) + geom_text( data = distinct(kegg_long, geneID), aes(x = -0.07, y = geneID, label = geneID), size = 3, hjust = 1 ) + scale_x_continuous( limits = c(-0.1, max(kegg_result$GeneRatio_num) * 1.1), expand = c(0, 0) )

4.2 高级布局技巧

当基因数量较多时，采用以下策略优化布局：

# 基因标签防重叠处理 sankey_plot <- sankey_plot + geom_text_repel( data = distinct(kegg_long, geneID), aes(x = -0.07, y = geneID, label = geneID), size = 3, hjust = 1, direction = "y", box.padding = 0.1, segment.color = NA ) + annotate("text", x = -0.05, y = Inf, label = "Gene List", hjust = 0.5, vjust = -1, fontface = "bold")

5. 输出与进阶定制

5.1 图形导出最佳实践

ggsave("KEGG_Sankey_Dotplot.pdf", plot = sankey_plot, width = 12, height = 8, device = cairo_pdf)

5.2 高级自定义选项

通过修改以下参数实现个性化效果：

custom_plot <- sankey_plot + scale_color_gradientn( colors = c("#4575b4", "#ffffbf", "#d73027"), values = scales::rescale(c(0, 0.5, 1)), breaks = seq(0, 10, by = 2) ) + geom_point( aes(size = Count, color = logP), shape = 21, fill = "white", stroke = 1 )

5.3 交互式版本实现

使用plotly创建可探索的交互式图表：

library(plotly) ggplotly(sankey_plot) %>% layout( hoverlabel = list( bgcolor = "white", font = list(size = 12) ), margin = list(l = 150) )

6. 实战问题解决指南

6.1 处理大量通路的显示问题

当通路超过20条时，建议：

1. 按p值筛选前N条显著通路

top_kegg <- kegg_result %>% arrange(pvalue) %>% head(20)

2. 使用分面显示

faceted_plot <- sankey_plot + facet_grid(cluster ~ ., scales = "free", space = "free")

6.2 基因名称重叠解决方案

• 对长基因列表实施抽样显示

set.seed(123) sampled_genes <- kegg_long %>% group_by(Description) %>% sample_n(min(5, n()))

• 使用缩写形式

kegg_long <- kegg_long %>% mutate(geneID_short = str_sub(geneID, 1, 8))

6.3 多组比较可视化

对于时间序列或多条件实验，可扩展为分组桑吉气泡图：

multi_group_plot <- ggplot() + geom_point( data = kegg_result, aes(x = GeneRatio_num, y = Description, size = Count, color = logP), position = position_dodge(width = 0.5) ) + facet_grid(. ~ Group)

7. 性能优化技巧

处理大规模数据集时（如>1000个基因关联）：

1. 使用data.table加速数据处理

library(data.table) kegg_dt <- as.data.table(kegg_result) kegg_long <- kegg_dt[, .(geneID = unlist(strsplit(geneID, "/"))), by = .(Description, GeneRatio, pvalue, Count)]

2. 简化图形元素

fast_plot <- ggplot(kegg_result) + geom_point(aes(...)) + geom_segment( data = kegg_long[sample(.N, 100)], aes(...), alpha = 0.3 )

3. 预计算中间结果

saveRDS(kegg_long, "kegg_preprocessed.rds")