R语言生信分析：一套代码同时生成KEGG气泡图与桑基流向图

这次我们来看一个R语言生信分析中的实用技巧：如何用一套代码同时生成KEGG气泡图和桑基流向图。对于做通路富集分析的研究者来说，KEGG气泡图能直观展示通路富集程度，而桑基流向图则能清晰呈现基因在不同通路间的流动关系。通常，我们需要分别写两套代码、调整两次数据格式才能得到这两张图，过程繁琐且容易出错。

这个教程的核心价值在于“一套代码，两图齐出”。它解决了生信分析中多图联动、数据一致性的痛点。你只需要准备一份标准的KEGG富集分析结果，运行同一段R脚本，就能同时获得出版级的气泡图和具有动态流向美感的桑基图。这不仅极大提升了绘图效率，也保证了数据分析结果在不同可视化图表间的高度统一。

本文将带你完整走通这个流程。我们会从最基础的环境准备和R包安装开始，一步步讲解数据格式的整理、核心绘图代码的编写与解读，直到最终图形的输出与美化。无论你是刚接触R语言的生信新手，还是希望优化现有分析流程的老手，这套方法都能让你快速上手，将富集分析结果转化为更具表现力的科研图表。

1. 核心能力速览

在深入代码之前，我们先通过下表快速了解本方案的核心特性和要求。

2. 适用场景与使用边界

这套方法并非万能，明确其适用边界能帮助你更好地应用它。

最适合谁用？

• 生信分析初学者：希望用一套清晰的代码模板，快速实现从数据到高级图表的跨越。
• 科研工作者：需要为论文或报告制作规范且美观的通路富集可视化图表，提升图表质量。
• 流程优化者：希望将重复的绘图工作脚本化、自动化，提高分析复现效率。

能解决什么问题？

1. 效率问题：告别分别绘制气泡图和桑基图时，数据准备、代码编写、样式调整的重复劳动。
2. 一致性问题：确保两张图基于完全相同的数据子集和筛选标准，避免因手动操作导致图表间结果不一致。
3. 美观度问题：提供可直接微调的绘图代码，生成达到出版要求的图表，避免从零开始调样的麻烦。

不适合什么场景？

• 超大规模数据分析：如果富集结果包含上千条通路、上万个基因，桑基图可能会因为元素过多而显得杂乱，此时需先进行严格筛选（如Top 20通路）。
• 非KEGG富集分析：本教程代码逻辑针对KEGG富集结果的数据结构设计。如果是GO、Reactome等其他数据库的富集结果，需要调整数据列的名称和含义。
• 完全交互式可视化：本文生成的是静态图片。如果需要可交互的、能悬停查看基因详情的高级桑基图，需采用networkD3等专门用于Web交互的包。

合规与伦理提醒：

• 本方法仅用于科研数据的可视化展示，所有分析应基于合法、合规获取的生物学数据。
• 在论文中使用图表时，应遵循目标期刊的图表格式要求，并注明所使用的工具和包（如ggplot2, ggsankey）。
• 确保对使用的基因数据集拥有相应的使用权限，特别是在涉及人类遗传数据时，需严格遵守相关的伦理审查和数据安全规定。

3. 环境准备与前置条件

工欲善其事，必先利其器。开始之前，请确保你的电脑已经搭建好基础的R语言工作环境。

3.1 软件安装

• R语言：访问 R语言官网 下载并安装最新版本的R。安装过程通常很简单，一路点击“下一步”即可。
• RStudio (推荐)：这是一个强大的R语言集成开发环境（IDE），能极大提升编码效率。从 RStudio官网 下载免费的Desktop版本并安装。

3.2 必备R包安装与加载我们将使用以下几个核心R包，请在新开的R脚本中运行以下命令进行安装。如果已安装，则跳过安装步骤直接加载。

# 安装CRAN上的包 install.packages(c("ggplot2", "dplyr", "tidyr", "stringr", "ggsankey", "RColorBrewer")) # 安装Bioconductor上的包（用于富集分析） if (!require("BiocManager", quietly = TRUE)) install.packages("BiocManager") BiocManager::install("clusterProfiler") BiocManager::install("org.Hs.eg.db") # 以人类为例，请根据你的物种替换，如 org.Mm.eg.db（小鼠）

安装完成后，在每次分析开始时，加载这些包：

library(ggplot2) library(dplyr) library(tidyr) library(stringr) library(ggsankey) library(RColorBrewer) library(clusterProfiler) library(org.Hs.eg.db) # 替换为你的物种数据库

3.3 工作目录与数据准备建议为本次分析创建一个专属文件夹，并将你的数据文件放入其中。

# 设置工作目录（请将路径替换为你自己的文件夹路径） setwd("/你的/项目/路径/KEGG_双图教程") # 检查工作目录 getwd()

你的数据应该是一个KEGG富集分析的结果文件。通常，你可以通过clusterProfiler::enrichKEGG()函数得到结果，并将其保存为.csv或.txt文件。假设我们有一个名为kegg_enrichment_result.csv的文件。

4. 数据加载与预处理

一套代码出两图的前提，是有一份格式规范、信息完整的输入数据。这一步至关重要。

4.1 读取数据

# 读取KEGG富集结果 kegg_df <- read.csv("kegg_enrichment_result.csv", stringsAsFactors = FALSE) # 查看数据前几行和结构，确保关键列存在 head(kegg_df) str(kegg_df)

一份典型的富集结果data.frame应包含以下关键列（列名可能因软件而异，后续代码需要对应调整）：

• ID: KEGG通路ID（如 hsa04110）
• Description: 通路描述
• GeneRatio: 富集到该通路的基因数/背景基因数（如 50/5000）
• BgRatio: 通路背景基因数/总背景基因数
• pvalue/p.adjust/qvalue: P值、校正P值
• Count: 富集到该通路的基因数目（重要）
• geneID: 富集到该通路的基因ID列表，以/分隔（最重要！用于桑基图）

4.2 数据清洗与格式化我们需要从原始结果中提取并计算绘图所需的信息。

# 假设我们的数据框列名如下，请根据你的实际数据调整 # ID, Description, p.adjust, Count, geneID # 1. 筛选显著富集通路（例如校正P值<0.05） sig_kegg <- kegg_df %>% filter(p.adjust < 0.05) # 2. 按富集因子（Enrichment Factor）或P值排序，并选择Top N条通路用于绘图（避免过多） # 计算富集因子 (Enrichment Factor) = (Count / Gene总数) / (通路背景基因数/总背景基因数) # 这里简化处理，我们常按Count或-p.adjust排序 top_n <- 15 # 选择前15条通路 plot_data <- sig_kegg %>% arrange(p.adjust) %>% # 按校正P值升序排列（越小越显著） slice_head(n = top_n) # 取前top_n行 # 3. 为气泡图准备数据：计算 -log10(p.adjust)，并确保Description顺序 plot_data <- plot_data %>% mutate(`-log10(p.adjust)` = -log10(p.adjust), Description = factor(Description, levels = rev(unique(Description)))) # 反转因子顺序，使图表从上到下与数据顺序一致 # 查看处理后的数据 print(plot_data[, c("Description", "Count", "p.adjust", "-log10(p.adjust)", "geneID")])

5. 功能一：绘制KEGG气泡图

气泡图是展示富集分析结果最直观的方式之一，X轴通常代表富集因子或基因比例，Y轴是通路名称，气泡大小代表基因数量，颜色代表显著性。

5.1 基础气泡图绘制

# 使用ggplot2绘制气泡图 bubble_plot <- ggplot(plot_data, aes(x = Count / sum(plot_data$Count), # X轴可用基因比例，也可用富集因子 y = Description, size = Count, color = `-log10(p.adjust)`)) + geom_point(alpha = 0.8) + # 绘制点，alpha控制透明度 scale_size(range = c(3, 8)) + # 控制气泡大小范围 scale_color_gradient(low = "blue", high = "red") + # 设置颜色梯度，从蓝（不显著）到红（显著） labs(x = "Gene Ratio", y = "KEGG Pathway", size = "Gene Count", color = "-log10(adj.P)", title = "KEGG Pathway Enrichment Analysis") + theme_bw(base_size = 12) + # 使用黑白主题，设置基础字体大小 theme(axis.text.y = element_text(color = "black"), plot.title = element_text(hjust = 0.5, face = "bold")) # 标题居中加粗 # 显示图形 print(bubble_plot)

5.2 高级美化与自定义出版级的图表需要更精细的调整。

bubble_plot_fancy <- bubble_plot + scale_color_gradientn(colors = rev(brewer.pal(11, "Spectral")), # 使用ColorBrewer的Spectral色带 name = "-log10(adj.P)") + scale_size_continuous(name = "Gene\nCount", # 图例名称换行 breaks = c(min(plot_data$Count), median(plot_data$Count), max(plot_data$Count))) + # 自定义图例断点 theme( panel.grid.major = element_line(linetype = "dashed", color = "grey90"), panel.grid.minor = element_blank(), axis.title = element_text(size = 14, face = "bold"), axis.text.x = element_text(size = 11), axis.text.y = element_text(size = 11, lineheight = 0.9), legend.title = element_text(size = 11, face = "bold"), legend.text = element_text(size = 10), plot.title = element_text(size = 16, hjust = 0.5, margin = margin(b = 15)) ) print(bubble_plot_fancy)

6. 功能二：绘制桑基流向图

桑基图的核心是展示“基因”如何流向不同的“通路”。因此，我们需要将geneID列（包含多个基因ID的字符串）进行“长格式”转换。

6.1 数据转换：从“通路-基因列表”到“基因-通路”关系对

# 准备桑基图数据：将geneID列拆分成多行 sankey_data <- plot_data %>% select(Description, geneID) %>% # 选择需要的列 mutate(geneID = str_split(geneID, "/")) %>% # 将geneID字符串按“/”分割成列表 unnest_longer(geneID) %>% # 将列表展开为多行，关键步骤！ rename(source = geneID, target = Description) # 重命名列，符合ggsankey要求：基因是source，通路是target # 查看转换后的数据格式：每一行代表一个基因属于一个通路 head(sankey_data)

6.2 基础桑基图绘制ggsankey包提供了geom_sankey()和geom_sankey_label()函数来方便地绘图。

# 计算每个节点（基因和通路）的频次，用于排序和标签 node_summary <- sankey_data %>% make_long(source, target) %>% # ggsankey的关键函数，将数据转换为绘图所需的长格式 dplyr::group_by(node) %>% dplyr::summarise(value = n()) # 绘制桑基图 sankey_plot <- ggplot(sankey_data, aes(x = x, # 阶段（自动分配） next_x = next_x, # 下一阶段 node = node, # 当前节点 next_node = next_node, # 下一个节点 fill = node, # 按节点填充颜色 label = node)) + # 节点标签 geom_sankey(flow.alpha = 0.6, # 流向透明度 node.color = "black", # 节点边框色 show.legend = FALSE) + # 不显示图例（节点过多） geom_sankey_label(size = 3.2, color = "black", fill = "white") + # 添加节点标签 scale_fill_viridis_d(option = "C", alpha = 0.8) + # 使用viridis离散色标 theme_sankey(base_size = 14) + # 桑基图主题 labs(x = NULL, y = NULL, title = "Gene - KEGG Pathway Sankey Diagram") + theme(axis.text.x = element_text(face = "bold"), plot.title = element_text(hjust = 0.5, size = 16, face = "bold"), axis.text.y = element_blank(), axis.ticks.y = element_blank(), panel.background = element_blank()) print(sankey_plot)

6.3 处理重叠与美化当基因或通路过多时，标签会重叠。我们可以通过筛选高频基因或调整布局来优化。

# 策略1：只显示连接数较多（如>2）的基因节点，减少杂乱 gene_freq <- table(sankey_data$source) genes_to_keep <- names(gene_freq[gene_freq > 2]) # 只保留出现在多于2条通路中的基因 sankey_data_filtered <- sankey_data %>% filter(source %in% genes_to_keep) # 使用过滤后的数据重新绘图（代码同上，将sankey_data替换为sankey_data_filtered） # ... (绘图代码省略，与6.2节类似) # 策略2：调整标签大小、角度和位置 sankey_plot_adjusted <- sankey_plot + geom_sankey_label(size = 2.8, color = "black", fill = "white", hjust = 0.5, # 水平居中 position = position_nudge(y = 0.05)) + # 垂直方向微调 theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) # X轴标签倾斜45度 print(sankey_plot_adjusted)

7. 一套代码整合与批量输出

现在，我们将两部分代码整合到一个脚本中，实现“一套代码，两图齐出”，并自动保存高清图片。

7.1 完整脚本示例将以下代码保存为一个R脚本文件（如draw_kegg_dual_plots.R）。

#!/usr/bin/env Rscript # 文件名：draw_kegg_dual_plots.R # 功能：读取KEGG富集结果，同时生成气泡图和桑基图 # 1. 加载包 library(ggplot2) library(dplyr) library(tidyr) library(stringr) library(ggsankey) library(RColorBrewer) # 2. 参数设置（用户可根据需要修改） input_file <- "kegg_enrichment_result.csv" # 输入文件名 output_prefix <- "my_kegg_analysis" # 输出图片前缀 top_n_pathways <- 15 # 绘制Top多少条通路 pval_cutoff <- 0.05 # 显著性阈值 width <- 12 # 图片宽度（英寸） height <- 10 # 图片高度（英寸） dpi <- 300 # 图片分辨率 # 3. 数据加载与预处理 kegg_raw <- read.csv(input_file, stringsAsFactors = FALSE) sig_kegg <- kegg_raw %>% filter(p.adjust < pval_cutoff) plot_data <- sig_kegg %>% arrange(p.adjust) %>% slice_head(n = top_n_pathways) %>% mutate(`-log10(p.adjust)` = -log10(p.adjust), Description = factor(Description, levels = rev(unique(Description)))) # 4. 绘制并保存气泡图 cat("正在生成气泡图...\n") bubble_p <- ggplot(plot_data, aes(x = Count / sum(Count), y = Description, size = Count, color = `-log10(p.adjust)`)) + geom_point(alpha = 0.7) + scale_size_continuous(range = c(4, 10), name = "Gene Count") + scale_color_gradientn(colors = brewer.pal(9, "YlOrRd"), name = "-log10(adj.P)") + labs(x = "Gene Ratio", y = NULL, title = "Top Enriched KEGG Pathways") + theme_minimal(base_size = 14) + theme(plot.title = element_text(hjust = 0.5, face = "bold", size = 16), axis.text.y = element_text(color = "black", size = 12), legend.position = "right") bubble_file <- paste0(output_prefix, "_bubble.png") ggsave(bubble_file, plot = bubble_p, width = width, height = height, dpi = dpi) cat("气泡图已保存至：", bubble_file, "\n") # 5. 准备桑基图数据并绘图 cat("正在生成桑基图...\n") sankey_df <- plot_data %>% select(Description, geneID) %>% mutate(geneID = str_split(geneID, "/")) %>% unnest_longer(geneID) %>% rename(source = geneID, target = Description) # 可选：过滤低频基因以简化图形 # gene_counts <- table(sankey_df$source) # keep_genes <- names(gene_counts[gene_counts > 1]) # sankey_df <- sankey_df %>% filter(source %in% keep_genes) sankey_p <- ggplot(sankey_df, aes(x = x, next_x = next_x, node = node, next_node = next_node, fill = node, label = node)) + geom_sankey(flow.alpha = 0.5, show.legend = FALSE) + geom_sankey_label(size = 3, fill = "white", color = "black") + scale_fill_viridis_d(option = "plasma") + theme_sankey(base_size = 16) + labs(title = "Gene to Pathway Flow (Sankey Diagram)") + theme(plot.title = element_text(hjust = 0.5, face = "bold"), axis.text.y = element_blank(), axis.ticks.y = element_blank()) sankey_file <- paste0(output_prefix, "_sankey.png") ggsave(sankey_file, plot = sankey_p, width = width, height = height, dpi = dpi) cat("桑基图已保存至：", sankey_file, "\n") cat("绘图完成！\n")

7.2 如何运行你可以在RStudio中打开这个脚本，全选并运行。或者，在终端（命令行）中，使用以下命令：

Rscript draw_kegg_dual_plots.R

脚本会自动读取当前目录下的kegg_enrichment_result.csv文件，生成my_kegg_analysis_bubble.png和my_kegg_analysis_sankey.png两张高清图片。

8. 常见问题与排查方法

在实际操作中，你可能会遇到一些问题。下表列出了常见错误及其解决方法。

9. 最佳实践与使用建议

掌握基础操作后，遵循以下建议能让你的分析更高效、图表更专业。

1. 项目文件管理：为每个分析项目建立独立的文件夹，里面包含：原始数据(/raw_data)、R脚本(/scripts)、中间结果(/intermediate)、最终图表(/figures)。使用相对路径，增强代码可移植性。
2. 数据备份与版本控制：对关键的预处理后数据（如plot_data）使用saveRDS()和readRDS()进行保存和读取，避免重复计算。考虑使用Git进行脚本的版本管理。
3. 参数脚本化：像上面的完整脚本一样，将top_n_pathways、pval_cutoff、颜色方案、图片尺寸等作为脚本开头的变量。修改参数只需改动一处，无需在代码中四处寻找。
4. 图表美学统一：如果你需要为多组数据生成图表，请固定颜色方案、字体大小、主题风格。可以创建一个自定义主题函数，确保所有图表风格一致。

   my_theme <- function(base_size=14) { theme_minimal(base_size = base_size) + theme(plot.title = element_text(hjust=0.5, face="bold"), axis.title = element_text(face="bold"), legend.position = "bottom") } # 然后在每个ggplot对象后加上 `+ my_theme()`

5. 桑基图复杂度控制：桑基图的美观度与数据量高度相关。“少即是多”。强烈建议将通路数量控制在10-20条，并通过过滤低频基因（如gene_freq > 2）来大幅减少连线，使图形重点突出、脉络清晰。
6. 探索交互式可视化：如果你需要在线报告或更灵活的探索，可以考虑使用networkD3包生成可交互的桑基图（D3.js驱动），它支持鼠标悬停查看详细信息。
7. 结果解读与标注：在论文中使用图表时，记得在图表标题或图注中清晰说明：通路筛选标准（如Top 15 by adj.P-value）、颜色和大小代表的含义、以及使用的工具（ggplot2, ggsankey）。

10. 总结与下一步

通过本教程，你已经掌握了用R语言“一套代码生成KEGG气泡图和桑基流向图”的核心流程。这套方法的核心优势在于效率与一致性：你只需要维护一份数据清洗和筛选的代码，就能驱动两种不同维度的可视化，极大减少了出错的可能。

最值得尝试的步骤，是将第7部分的完整脚本应用到你自己的KEGG富集结果上。你可能会遇到的第一个挑战是数据列名不匹配，请根据错误提示耐心调整。第二个挑战可能是桑基图过于复杂，请果断使用top_n_pathways和基因频率过滤来简化它。

最容易踩的坑是忽略了数据的预处理。务必确保你的输入数据包含Description、p.adjust、Count和geneID（或类似）这几列，并且geneID列中的基因是以特定分隔符（如“/”）连接的字符串。这是后续所有操作的基础。

掌握了这个基本流程后，你可以进行多方面的扩展：

• 扩展至GO分析：尝试将代码适配到GO（Gene Ontology）的富集结果，通常需要处理BP、CC、MF三个分类。
• 组合图表：使用patchwork或cowplot包将气泡图和桑基图并排组合在一张图中，形成更强大的展示效果。
• 自动化报告：结合rmarkdown或quarto，将数据读取、分析、绘图和结果解读整合成一个可重复生成的HTML或PDF报告。
• 性能优化：如果你的数据量极大，可以考虑使用data.table包替代dplyr进行数据操作，提升处理速度。

建议将本文的核心脚本收藏或保存到你的代码库中。下次当你拿到一份KEGG富集结果，需要快速产生高质量可视化时，它将成为你的得力工具。