R语言一键生成KEGG气泡图与桑基图：生信结果可视化完整工作流

如果你在生物信息学分析中，已经完成了差异表达基因的筛选和KEGG通路富集分析，面对一长串富集结果，下一步该怎么做？

直接贴一个满是P值和基因数的表格到论文里？审稿人大概率会皱眉头。用Excel手动做个柱状图？不仅费时费力，而且信息呈现单一，难以同时展示通路的统计学显著性和基因的流向关系。

这正是很多生信初学者和科研工作者遇到的真实瓶颈：分析做完了，但不知道如何高效、专业地将结果可视化，做出既能上得了台面，又能讲清楚故事的图表。

今天要解决的问题就是：如何用同一套R语言代码和数据，一键生成两种科研级图表——直观展示富集显著性的KEGG气泡图，和清晰揭示基因-通路关联的桑基图（Sankey Diagram）。

本文不会停留在“调用某个函数出图”的层面。我们将深入一套完整的工作流：从富集结果表格的整理开始，到使用ggplot2绘制可高度定制化的气泡图，再到利用networkD3包构建交互式桑基图，最后实现两图联动，用数据讲述一个完整的生物学故事。无论你是刚接触R语言的生信新手，还是想优化分析流程的老手，这套方法都能直接提升你的成果展示效率与专业度。

1. 核心价值：为什么你需要掌握“一码两图”？

在生信数据分析中，可视化不是分析的终点，而是沟通的起点。KEGG富集分析的结果至少包含两个维度的信息：通路的显著性（如P值、Q值、富集因子）和基因与通路的归属关系。单一图表往往只能强调其中一个维度。

• KEGG气泡图：擅长展示“哪些通路最显著”。通过点的大小（基因数）和颜色（P值），一眼就能锁定关键通路，适用于结果摘要和初步结论展示。
• 桑基图：擅长展示“基因具体是如何分配到各个通路中的”。它能清晰呈现基因与通路之间的多对多关系，直观显示哪些基因是“多面手”（参与多个通路），哪些通路共享关键基因，非常适合深入机制探讨。

传统做法是分别用不同工具（如R的ggplot2画气泡图，再用Origin或在线工具画桑基图）处理，流程割裂，且当数据更新时需要重复劳动。而“一码两图”的核心优势在于：

1. 效率提升：基于同一份清洁数据，运行一次脚本，同时获得两种视角的图表，避免重复数据整理。
2. 一致性保证：两图来源于同一数据源，确保在汇报或论文中数据口径绝对一致。
3. 深度挖掘：气泡图帮你筛选出TOP通路，桑基图则帮你深入解读这些通路背后的基因联系，形成分析闭环。
4. 灵活性高：全部流程在R中完成，从数据预处理到图形美化，每一步都可控、可复现。

接下来，我们将从零开始，拆解整个流程。

2. 环境准备：构建你的R生信绘图工作流

工欲善其事，必先利其器。本节将确保你拥有一个可稳定运行的环境。

2.1 R与RStudio的安装与配置

首先，你需要安装R语言和RStudio集成开发环境（IDE）。R是引擎，RStudio是方向盘和仪表盘，能极大提升编码效率。

• 安装R：访问R官网（https://www.r-project.org/），选择与操作系统对应的CRAN镜像下载安装。建议选择较新的稳定版本（如4.3.x）。
• 安装RStudio：访问RStudio官网（https://posit.co/download/rstudio-desktop/）下载免费的Desktop版本并安装。

安装后，打开RStudio，界面通常分为四个窗格：脚本编辑器、控制台、环境/历史、文件/图/帮助。

2.2 必需R包的安装与加载

我们将主要依赖以下三个包：

• ggplot2：绘图语法之王，用于绘制静态气泡图。
• dplyr/tidyr：数据整理神器，用于清洗和转换富集结果。
• networkD3：用于创建交互式桑基图（D3.js的R接口）。

在RStudio的控制台（Console）中，一次性安装并加载它们：

# 安装包（如果尚未安装） install.packages(c("ggplot2", "dplyr", "tidyr", "networkD3")) # 加载包到当前会话 library(ggplot2) library(dplyr) library(tidyr) library(networkD3)

注意：networkD3包生成的桑基图是HTML交互式图表，可在浏览器中查看，支持鼠标悬停查看详细信息，这比静态图更具探索性。

2.3 准备示例数据

为了演示，我们模拟一份典型的KEGG富集分析结果。在实际项目中，你只需将你的结果文件（通常是.csv或.txt格式）读入即可。

# 模拟创建一份富集结果数据框 set.seed(123) # 确保结果可重复 kegg_result <- data.frame( Pathway = c("Metabolic pathways", "Biosynthesis of secondary metabolites", "Microbial metabolism in diverse environments", "Carbon metabolism", "Biosynthesis of amino acids", "2-Oxocarboxylic acid metabolism", "Fatty acid metabolism", "Degradation of aromatic compounds", "ABC transporters", "Quorum sensing"), GeneRatio = c(120/2500, 85/2500, 78/2500, 65/2500, 54/2500, 48/2500, 42/2500, 38/2500, 35/2500, 30/2500), BgRatio = c(1500/8000, 1200/8000, 1100/8000, 900/8000, 800/8000, 700/8000, 600/8000, 500/8000, 450/8000, 400/8000), pvalue = c(1.2e-12, 3.5e-09, 2.1e-07, 8.7e-06, 4.5e-05, 0.00012, 0.00045, 0.0012, 0.0033, 0.0088), p.adjust = c(1.2e-10, 1.8e-07, 7.0e-06, 2.2e-04, 9.0e-04, 0.0020, 0.0060, 0.0150, 0.0330, 0.0660), qvalue = c(1.0e-10, 1.5e-07, 5.8e-06, 1.8e-04, 7.5e-04, 0.0017, 0.0050, 0.0125, 0.0275, 0.0550), geneID = c("geneA,geneB,geneC,geneD,geneE,geneF,geneG", "geneH,geneI,geneJ,geneK,geneL", "geneM,geneN,geneO,geneP,geneQ,geneR", "geneS,geneT,geneU,geneV", "geneW,geneX,geneY,geneZ,geneAA", "geneB,geneF,geneL,geneP,geneV", "geneC,geneI,geneO,geneU,geneAA,geneBB", "geneD,geneJ,geneP,geneV,geneCC", "geneE,geneK,geneQ,geneW,geneDD", "geneG,geneR,geneX,geneBB,geneCC,geneDD") ) # 查看数据结构 head(kegg_result)

这份数据包含了通路名称、富集基因比例、背景基因比例、P值、校正P值、Q值以及富集到的基因ID列表（以逗号分隔）。geneID列是连接气泡图和桑基图的关键。

3. 数据预处理：从原始结果到绘图就绪数据

原始数据通常不能直接用于绘图，尤其是geneID列。我们需要将其转换为适合ggplot2和networkD3的格式。

3.1 清洗与筛选富集结果

首先，我们通常只关注最显著的一些通路，比如选择校正P值（p.adjust）小于0.05的，并按P值排序。

# 筛选显著通路，并按p.adjust升序排序 kegg_sig <- kegg_result %>% filter(p.adjust < 0.05) %>% arrange(p.adjust) # 为了演示，我们选取前8条最显著的通路 top_pathways <- head(kegg_sig, 8) print(top_pathways)

3.2 关键步骤：拆分基因列表，构建“基因-通路”关联表

这是生成桑基图的核心准备步骤。我们需要将geneID列（一个包含多个基因的字符串）拆分成多行，每行代表一个基因与其所属通路的关系。

# 拆分geneID列，构建长格式的基因-通路关联表 gene_pathway_df <- top_pathways %>% select(Pathway, geneID) %>% # 选择需要的列 separate_rows(geneID, sep = ",") %>% # 按逗号拆分，一行变多行 rename(Gene = geneID) %>% # 重命名列 mutate(Gene = trimws(Gene)) # 去除基因名两端的空格 # 查看关联表的前几行 head(gene_pathway_df)

现在，gene_pathway_df数据框的每一行都是一个明确的Gene属于某个Pathway的记录。这个格式完美契合桑基图对“源-目标”链接数据的要求。

3.3 为绘图准备衍生数据

对于气泡图，我们可能希望用-log10(p.adjust)来表示显著性，使得值越大点颜色越深（越显著），同时计算富集因子（Enrichment Factor）。

# 为气泡图准备数据：计算 -log10(p.adjust) 和富集因子 bubble_data <- top_pathways %>% mutate( `-log10(p.adjust)` = -log10(p.adjust), EnrichmentFactor = (GeneRatio) / (BgRatio) # 简化计算，实际需注意格式转换 ) %>% # 重新调整GeneRatio的格式，便于理解 mutate(GeneCount = as.numeric(sub("/.*", "", GeneRatio))) # 提取基因数 head(bubble_data)

至此，我们得到了两个核心数据对象：用于气泡图的bubble_data和用于桑基图的gene_pathway_df。它们同源，但形态各异。

4. 核心图表绘制：ggplot2气泡图实战

我们将使用ggplot2的图层语法，逐步构建一个出版级的气泡图。

4.1 基础气泡图绘制

最基本的映射关系是：X轴=富集因子（或GeneRatio），Y轴=通路名称，点大小=基因数，点颜色=显著性。

# 基础气泡图 p_bubble <- ggplot(bubble_data, aes(x = EnrichmentFactor, y = reorder(Pathway, `-log10(p.adjust)`), # 按显著性排序通路 size = GeneCount, color = `-log10(p.adjust)`)) + geom_point(alpha = 0.8) + # 添加点图层，设置透明度 scale_size_area(name = "Gene Count", max_size = 12) + # 控制点大小范围 scale_color_gradient(low = "blue", high = "red", name = "-log10(p.adjust)") + # 设置颜色渐变 labs(x = "Enrichment Factor", y = "Pathway", title = "KEGG Pathway Enrichment Analysis") + theme_bw(base_size = 14) + # 使用黑白主题，设置基础字体大小 theme(axis.text.y = element_text(size = 10, color = "black"), plot.title = element_text(hjust = 0.5, face = "bold")) # 标题居中加粗 print(p_bubble)

这段代码会生成一个可用的气泡图。但我们可以做得更好。

4.2 高级美化与定制

科研图表讲究清晰、准确、美观。以下是一些常见的优化技巧：

# 高级美化版气泡图 p_bubble_enhanced <- p_bubble + # 1. 优化图例 guides(size = guide_legend(order = 1), color = guide_colorbar(order = 2)) + # 2. 调整坐标轴和网格线 theme( panel.grid.major.y = element_line(linetype = "dashed", color = "grey90"), # 横向虚线网格 panel.grid.major.x = element_blank(), # 去除纵向主网格线 panel.grid.minor = element_blank(), # 去除次要网格线 axis.line.x = element_line(color = "black"), # X轴线 axis.ticks.y = element_blank() # 去除Y轴刻度线 ) + # 3. 扩展颜色标度，让极值更突出 scale_color_gradientn( colours = c("#4393C3", "#FFD700", "#D73027"), name = "-log10(p.adjust)" ) + # 4. 添加数值标签（可选，在点上显示基因数） geom_text(aes(label = GeneCount), color = "white", size = 3, show.legend = FALSE) print(p_bubble_enhanced)

通过theme()函数，你可以精细控制图表的每一个元素。现在，你的气泡图已经具备了投稿期刊的潜力。

5. 核心图表绘制：networkD3桑基图实战

桑基图需要一种特定的数据格式：一个包含“链接”（Links）和“节点”（Nodes）的数据框。

5.1 构建桑基图数据格式

“链接”数据框需要三列：source（源节点索引）、target（目标节点索引）、value（链接权重，通常为1）。“节点”数据框需要一列name，按顺序列出所有唯一的节点名。

# 1. 准备节点列表：包含所有唯一的基因和通路 # 注意：节点顺序至关重要，它将决定索引号。 all_genes <- unique(gene_pathway_df$Gene) all_pathways <- unique(gene_pathway_df$Pathway) node_names <- c(all_genes, all_pathways) # 通常将“源”（基因）放在前面，“目标”（通路）放在后面 # 创建节点数据框 nodes <- data.frame(name = node_names, stringsAsFactors = FALSE) # 2. 准备链接数据框 # 为每个链接找到对应的源节点和目标节点索引 links <- gene_pathway_df %>% mutate( source = match(Gene, node_names) - 1, # networkD3索引从0开始 target = match(Pathway, node_names) - 1, value = 1 # 每个链接的权重设为1 ) %>% select(source, target, value) # 查看链接数据前几行 head(links)

关键点：networkD3要求索引从0开始。match()函数返回的是在node_names向量中的位置（R索引从1开始），所以需要减1。

5.2 绘制交互式桑基图

使用sankeyNetwork()函数，传入链接和节点数据。

# 绘制基础桑基图 sankey_plot <- sankeyNetwork(Links = links, Nodes = nodes, Source = "source", Target = "target", Value = "value", NodeID = "name", units = "genes", # 链接的单位 fontSize = 12, nodeWidth = 20, height = 600, width = 900) # 在RStudio的Viewer窗格中显示（交互式） sankey_plot # 保存为独立的HTML文件，可在浏览器中打开并交互 saveNetwork(sankey_plot, file = "KEGG_Sankey_Diagram.html")

运行后，你会在RStudio的Viewer窗口看到一个可交互的流程图。鼠标悬停在节点（基因或通路）上，会高亮显示所有与之相连的流；悬停在链接上，会显示详细信息。这极大地便利了数据探索。

5.3 桑基图的美化与问题处理

默认的桑基图可能颜色单一、节点拥挤。我们可以进行优化：

# 美化桑基图：为基因和通路节点设置不同颜色 # 假设我们想用蓝色系表示基因，橙色系表示通路 node_colors <- c(rep("steelblue", length(all_genes)), # 基因节点颜色 rep("darkorange", length(all_pathways))) # 通路节点颜色 sankey_plot_enhanced <- sankeyNetwork(Links = links, Nodes = nodes, Source = "source", Target = "target", Value = "value", NodeID = "name", units = "genes", fontSize = 14, nodeWidth = 25, nodePadding = 15, # 增加节点间距 height = 700, width = 1000, colourScale = JS('d3.scaleOrdinal().range(["#999"])'), # 链接颜色 NodeGroup = "name", # 用于分组的列，这里我们用名字，但配合自定义颜色 LinkGroup = "source", # 链接按源节点分组着色 sinksRight = TRUE) # 右侧节点对齐 # 直接修改HTML对象的样式是复杂的，更简单的办法是在保存后编辑HTML，或使用更高级的包。 # 一个实用技巧：通过调整nodeWidth, nodePadding, height, width来改善布局。

如果基因和通路数量过多，桑基图会变得非常复杂难以阅读。最佳实践是：先用气泡图筛选出最显著的少数几个通路（如5-10个），再用这些通路对应的基因子集来绘制桑基图。这正是“一码两图”工作流的精髓：气泡图用于筛选，桑基图用于深挖。

6. 流程整合与自动化脚本

将以上步骤整合到一个R脚本或函数中，即可实现“一键出两图”。

# 文件名：kegg_dual_plot.R # 功能：输入KEGG富集结果数据框，输出气泡图和桑基图 generate_kegg_plots <- function(enrichment_df, p_adjust_cutoff = 0.05, top_n = 8) { # 加载必要库（如果在函数外未加载） library(ggplot2) library(dplyr) library(tidyr) library(networkD3) # 1. 数据筛选与排序 sig_data <- enrichment_df %>% filter(p.adjust < p_adjust_cutoff) %>% arrange(p.adjust) %>% head(top_n) if (nrow(sig_data) == 0) { stop("No significant pathways found with the given cutoff.") } # 2. 准备气泡图数据 bubble_data <- sig_data %>% mutate(`-log10(p.adjust)` = -log10(p.adjust), GeneCount = as.numeric(sub("/.*", "", GeneRatio)), EnrichmentFactor = GeneCount / as.numeric(sub(".*/", "", BgRatio))) # 3. 绘制气泡图 p_bubble <- ggplot(bubble_data, aes(x = EnrichmentFactor, y = reorder(Pathway, `-log10(p.adjust)`), size = GeneCount, color = `-log10(p.adjust)`)) + geom_point(alpha = 0.7) + scale_size_area(max_size = 10, name = "Gene Count") + scale_color_gradient(low = "lightblue", high = "red", name = "-log10(p.adjust)") + labs(x = "Enrichment Factor", y = NULL, title = "Top KEGG Enriched Pathways") + theme_minimal(base_size = 12) + theme(axis.text.y = element_text(size = 10), plot.title = element_text(hjust = 0.5, face = "bold"), legend.position = "right") # 4. 准备桑基图数据 gene_pathway_long <- sig_data %>% select(Pathway, geneID) %>% separate_rows(geneID, sep = ",") %>% mutate(geneID = trimws(geneID)) %>% rename(Gene = geneID) all_genes <- unique(gene_pathway_long$Gene) all_pathways <- unique(gene_pathway_long$Pathway) node_names <- c(all_genes, all_pathways) nodes <- data.frame(name = node_names) links <- gene_pathway_long %>% mutate(source = match(Gene, node_names) - 1, target = match(Pathway, node_names) - 1, value = 1) %>% select(source, target, value) # 5. 绘制桑基图 sankey_plot <- sankeyNetwork(Links = links, Nodes = nodes, Source = "source", Target = "target", Value = "value", NodeID = "name", fontSize = 10, nodeWidth = 20, height = 500, width = 800, sinksRight = TRUE) # 6. 返回结果 return(list(bubble_plot = p_bubble, sankey_plot = sankey_plot, bubble_data = bubble_data, sankey_data = list(links = links, nodes = nodes))) } # 使用函数 # 假设你的富集结果在 `my_kegg_results` 数据框中 plots <- generate_kegg_plots(my_kegg_results, p_adjust_cutoff = 0.05, top_n = 6) # 查看气泡图 print(plots$bubble_plot) # 查看并保存桑基图 plots$sankey_plot saveNetwork(plots$sankey_plot, file = "My_Analysis_Sankey.html")

这个函数封装了核心流程，你只需要提供自己的enrichment_df，调整p_adjust_cutoff和top_n参数，即可快速生成图表。

7. 常见问题与排查指南

在实际操作中，你可能会遇到以下问题：

最重要的建议：始终从一个小规模的、可控的测试数据集开始（比如只选3个通路），确保每一步的代码都按预期工作，然后再应用到全数据集上。

8. 最佳实践与进阶技巧

掌握了基础流程后，这些技巧能让你的分析更上一层楼：

1. 数据溯源与可复现性：

   • 在脚本开头使用set.seed()保证随机过程可重复。
   • 使用sessionInfo()记录R和包的版本。
   • 将原始数据、处理脚本和最终图表放在同一个项目目录下。

2. 桑基图性能优化：

   • 当基因数过多时，桑基图会变得极其复杂。考虑在桑基图中只展示差异最显著的基因（如logFC绝对值最大的Top 50），而不是所有富集基因。
   • 使用networkD3的nodePadding和margin参数调整布局，避免节点挤在一起。

3. 气泡图的美学定制：

   • 使用scale_color_gradient2()可以设置中间色（如白色），让颜色对比更柔和。
   • 通过theme(legend.position = “bottom”)将图例放在底部，节省纵向空间。
   • 使用ggsave(“bubble_plot.png”, width=8, height=6, dpi=300)导出高清图片用于投稿。

4. 结果的生物学解读：

   • 气泡图帮你找到“什么”通路重要。
   • 桑基图帮你回答“为什么”重要——通过展示哪些核心基因同时参与了多个关键通路，提示潜在的调控枢纽。
   • 将两图并列放在报告或论文中，并配文说明：左图展示了显著性排名靠前的通路，右图揭示了这些通路之间通过共享基因形成的功能网络。

5. 扩展到其他富集分析：

   • 本流程不仅适用于KEGG，稍作修改即可用于GO、Reactome、MSigDB等任何提供基因列表的富集分析结果。关键在于结果表中需包含geneID这类列。

从混乱的富集结果表格，到直观的气泡图和揭示内在联系的桑基图，你不仅完成了一次可视化升级，更构建了一个从宏观显著性筛选到微观基因网络探查的完整分析叙事。这套基于R的“一码两图”工作流，其价值在于将固定的分析模式转化为可复用的自动化脚本。下次当你拿到新的测序数据并完成富集分析后，只需将结果文件路径指向这个脚本，几分钟内就能获得可用于组会、报告或论文插图的专业图表。