用R语言打造5维桑吉气泡图从clusterProfiler结果到高级可视化在生物信息学研究中KEGG和GO富集分析几乎是每个转录组项目的标配。传统的富集气泡图虽然能展示通路名称、富集程度、p值和基因数四个维度的信息但关键的基因列表却往往被隐藏在表格中。本文将带你突破这一局限使用ggplot2和ggsankey包直接从clusterProfiler结果对象中提取数据构建包含基因关联信息的5维桑吉气泡图。1. 准备工作与环境配置在开始之前确保你已经安装了必要的R包。如果你使用Bioconductor管理生物信息学相关的包可以用以下命令安装clusterProfilerif (!require(BiocManager, quietly TRUE)) install.packages(BiocManager) BiocManager::install(clusterProfiler)对于可视化部分我们需要ggplot2作为基础以及ggsankey来创建桑吉图元素。ggsankey不是CRAN官方包需要从GitHub安装install.packages(ggplot2) install.packages(remotes) remotes::install_github(davidsjoberg/ggsankey)提示如果你遇到网络问题导致GitHub包安装失败可以尝试先安装devtools包然后使用devtools::install_github()函数。加载所有需要的包library(clusterProfiler) library(ggplot2) library(ggsankey) library(dplyr) library(tidyr)2. 从clusterProfiler结果中提取5维数据clusterProfiler的富集分析结果通常存储在一个特定的对象中如enrichResult。我们需要从中提取五个关键维度的数据Description通路或GO term的名称GeneRatio富集基因比例pvalue统计显著性geneID富集基因列表Count富集基因数量假设我们已经有一个名为kegg_enrich的富集结果对象下面是提取和整理数据的代码# 提取富集结果数据框 enrich_df - as.data.frame(kegg_enrich) # 处理GeneRatio为数值 enrich_df - enrich_df %% separate(GeneRatio, into c(GeneInPathway, TotalGenes), sep /) %% mutate(GeneRatio as.numeric(GeneInPathway)/as.numeric(TotalGenes)) # 展开geneID列为桑吉图准备数据 sankey_data - enrich_df %% select(Description, geneID) %% separate_rows(geneID, sep /) %% rename(gene geneID, pathway Description)这段代码完成了几个关键操作将GeneRatio从a/b格式转换为数值将geneID列中的基因列表用/分隔拆分为多行为桑吉图准备了pathway-gene的对应关系数据3. 构建基础气泡图在添加桑吉图元素之前我们先创建传统的4维气泡图。这将成为我们5维可视化的基础base_plot - ggplot(enrich_df, aes(x GeneRatio, y reorder(Description, GeneRatio), size Count, color -log10(pvalue))) geom_point(alpha 0.8) scale_color_gradient(low blue, high red, name -log10(pvalue)) scale_size(range c(3, 10), name Gene Count) labs(x Gene Ratio, y Pathway) theme_minimal() theme(axis.text.y element_text(size 10), legend.position right)这个气泡图已经包含了四个维度的信息Y轴通路名称按GeneRatio排序X轴GeneRatio富集程度点大小Count富集基因数量点颜色-log10(pvalue)统计显著性4. 添加桑吉图元素展示基因信息现在我们将在气泡图左侧添加桑吉图来展示第五个维度——基因列表。这需要一些数据转换和巧妙的图形组合# 创建桑吉图 sankey_plot - ggplot(sankey_data, aes(x x, next_x next_x, node node, next_node next_node, fill factor(node), label node)) geom_sankey(flow.alpha 0.5, node.color gray30) geom_sankey_label(size 3, color black, fill white) scale_fill_viridis_d(option D, alpha 0.8) theme_sankey(base_size 12) theme(legend.position none, axis.title element_blank(), axis.text element_blank(), axis.ticks element_blank(), panel.grid element_blank()) # 组合两个图形 library(patchwork) final_plot - sankey_plot base_plot plot_layout(widths c(1, 2))注意ggsankey需要特定的数据格式其中x表示当前阶段如gene或pathwaynext_x表示下一阶段node和next_node表示连接的节点。5. 高级定制与美化为了让图形更具发表质量我们可以进行一系列的美化调整# 自定义颜色和主题 final_plot - final_plot theme(plot.margin margin(1, 1, 1, 1, cm), plot.background element_rect(fill white, color NA), panel.background element_rect(fill white, color NA)) # 调整图例和标签 final_plot[[2]] - final_plot[[2]] guides(size guide_legend(override.aes list(color darkgray)), color guide_colorbar(barwidth 1, barheight 10)) labs(title KEGG Pathway Enrichment Analysis, subtitle Bubble size represents gene count, color shows -log10(p-value)) # 调整桑吉图部分的宽度比例 final_plot - final_plot plot_layout(widths c(0.7, 2))对于特别大的基因集桑吉图可能会显得过于拥挤。这时可以考虑以下策略筛选显著通路只展示p值最显著的前20个通路enrich_df - enrich_df %% arrange(pvalue) %% head(20)合并相似通路手动或使用语义相似度方法合并冗余通路基因缩写对于特别长的基因列表可以只显示部分代表性基因6. 常见问题与调试技巧在实际操作中你可能会遇到以下问题问题1桑吉图节点标签重叠解决方案调整图形大小或字体大小sankey_plot - sankey_plot geom_sankey_label(size 2.5)问题2气泡图点的大小范围不合适解决方案调整scale_size的参数scale_size(range c(2, 8)) # 根据你的数据调整这两个数字问题3颜色梯度不明显解决方案使用更鲜明的颜色方案或调整p值的转换方式scale_color_gradient2(low blue, mid yellow, high red, midpoint median(-log10(enrich_df$pvalue)))问题4图形元素不对齐解决方案确保两个图形的y轴顺序一致enrich_df - enrich_df %% arrange(desc(GeneRatio)) sankey_data - sankey_data %% mutate(pathway factor(pathway, levels enrich_df$Description))对于特别复杂的富集结果我通常会先导出数据到CSV在Excel中手动筛选和整理然后再导入R中进行可视化write.csv(enrich_df, kegg_enrichment_results.csv, row.names FALSE) # 手动编辑后 enrich_df - read.csv(kegg_enrichment_results_edited.csv)这种5维桑吉气泡图不仅适用于KEGG通路分析也可以应用于GO富集分析、疾病关联分析等各种场景。关键在于理解数据的结构并灵活调整可视化参数以适应不同的数据集。