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

这次我们来看一个R语言生信分析中的实用技巧如何用一套代码同时生成KEGG气泡图和桑基流向图。对于做通路富集分析的研究者来说KEGG气泡图能直观展示通路富集程度而桑基流向图则能清晰呈现基因在不同通路间的流动关系。通常我们需要分别写两套代码、调整两次数据格式才能得到这两张图过程繁琐且容易出错。这个教程的核心价值在于“一套代码两图齐出”。它解决了生信分析中多图联动、数据一致性的痛点。你只需要准备一份标准的KEGG富集分析结果运行同一段R脚本就能同时获得出版级的气泡图和具有动态流向美感的桑基图。这不仅极大提升了绘图效率也保证了数据分析结果在不同可视化图表间的高度统一。本文将带你完整走通这个流程。我们会从最基础的环境准备和R包安装开始一步步讲解数据格式的整理、核心绘图代码的编写与解读直到最终图形的输出与美化。无论你是刚接触R语言的生信新手还是希望优化现有分析流程的老手这套方法都能让你快速上手将富集分析结果转化为更具表现力的科研图表。1. 核心能力速览在深入代码之前我们先通过下表快速了解本方案的核心特性和要求。能力项说明核心功能使用同一套数据与核心代码一次性生成KEGG通路富集气泡图和桑基流向图。主要技术栈R语言依赖ggplot2,ggsankey,dplyr,clusterProfiler(或enrichplot) 等包。输入数据标准KEGG富集分析结果data.frame格式通常包含通路描述、基因数、P值、富集因子、基因ID列表等列。输出成果1. 静态高清KEGG气泡图.png/.pdf2. 静态高清桑基流向图.png/.pdf注可通过调整参数保存为不同格式环境门槛安装R和RStudio或其它IDE具备基础的R语言操作与包安装能力。性能要求对硬件无特殊要求普通笔记本电脑即可运行。图形渲染速度取决于通路和基因数量。适合场景科研论文图表制作、结题报告可视化、生信数据分析流程自动化、结果复现与展示。2. 适用场景与使用边界这套方法并非万能明确其适用边界能帮助你更好地应用它。最适合谁用生信分析初学者希望用一套清晰的代码模板快速实现从数据到高级图表的跨越。科研工作者需要为论文或报告制作规范且美观的通路富集可视化图表提升图表质量。流程优化者希望将重复的绘图工作脚本化、自动化提高分析复现效率。能解决什么问题效率问题告别分别绘制气泡图和桑基图时数据准备、代码编写、样式调整的重复劳动。一致性问题确保两张图基于完全相同的数据子集和筛选标准避免因手动操作导致图表间结果不一致。美观度问题提供可直接微调的绘图代码生成达到出版要求的图表避免从零开始调样的麻烦。不适合什么场景超大规模数据分析如果富集结果包含上千条通路、上万个基因桑基图可能会因为元素过多而显得杂乱此时需先进行严格筛选如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如 hsa04110Description: 通路描述GeneRatio: 富集到该通路的基因数/背景基因数如 50/5000BgRatio: 通路背景基因数/总背景基因数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. 常见问题与排查方法在实际操作中你可能会遇到一些问题。下表列出了常见错误及其解决方法。问题现象可能原因排查方式解决方案安装ggsankey包失败包不在CRAN上或依赖问题。查看错误信息是否提示package ‘ggsankey’ is not available。从GitHub安装devtools::install_github(davidsjoberg/ggsankey)。确保已安装devtools包。桑基图节点标签重叠严重基因或通路节点过多画布空间不足。检查sankey_df的行数如果超过几百行图形会非常拥挤。1. 在数据预处理阶段通过top_n_pathways减少通路数量。2. 在桑基图数据准备阶段过滤掉只出现在一条通路中的基因低频基因。3. 减小标签字体size或增加图片width和height。气泡图的Y轴通路顺序不对Description列没有被正确转换为因子或因子水平顺序错误。使用class(plot_data$Description)检查是否为因子用levels(plot_data$Description)查看顺序。在创建plot_data时确保使用mutate(Description factor(Description, levels rev(unique(Description))))来根据数据出现的顺序通常是按P值排序后的顺序设置因子水平。错误geneID列不存在你的富集结果文件中基因列表列的名称不是geneID。运行colnames(kegg_raw)查看所有列名。在代码中将所有出现geneID的地方替换为你数据中实际的列名例如可能是geneIds、genes或GeneSymbol。桑基图颜色杂乱或不好看scale_fill_viridis_d()的配色方案不适合。尝试不同的option参数如A,B,C,D。可以改用scale_fill_manual()自定义颜色或使用RColorBrewer包提供的调色板例如scale_fill_brewer(palette Set3)。运行脚本时报错找不到对象脚本中的变量名与数据中的列名不匹配或包未加载。仔细检查错误信息指向的行号对比变量名和列名。1. 确保在脚本开头正确加载了所有必要的包。2. 使用head(your_data)确认数据框的列名并据此修改脚本中的列名引用。生成的图片空白或布局异常图片保存尺寸太小或图形设备问题。检查保存的图片文件大小如果异常小几KB可能是保存失败。1. 增加ggsave()中的width和height参数值例如从10增加到15。2. 尝试先直接在RStudio的Plots窗口显示图形print(p)确认无误后再保存。3. 尝试保存为.pdf格式矢量图通常更稳定。9. 最佳实践与使用建议掌握基础操作后遵循以下建议能让你的分析更高效、图表更专业。项目文件管理为每个分析项目建立独立的文件夹里面包含原始数据(/raw_data)、R脚本(/scripts)、中间结果(/intermediate)、最终图表(/figures)。使用相对路径增强代码可移植性。数据备份与版本控制对关键的预处理后数据如plot_data使用saveRDS()和readRDS()进行保存和读取避免重复计算。考虑使用Git进行脚本的版本管理。参数脚本化像上面的完整脚本一样将top_n_pathways、pval_cutoff、颜色方案、图片尺寸等作为脚本开头的变量。修改参数只需改动一处无需在代码中四处寻找。图表美学统一如果你需要为多组数据生成图表请固定颜色方案、字体大小、主题风格。可以创建一个自定义主题函数确保所有图表风格一致。my_theme - function(base_size14) { theme_minimal(base_size base_size) theme(plot.title element_text(hjust0.5, facebold), axis.title element_text(facebold), legend.position bottom) } # 然后在每个ggplot对象后加上 my_theme()桑基图复杂度控制桑基图的美观度与数据量高度相关。“少即是多”。强烈建议将通路数量控制在10-20条并通过过滤低频基因如gene_freq 2来大幅减少连线使图形重点突出、脉络清晰。探索交互式可视化如果你需要在线报告或更灵活的探索可以考虑使用networkD3包生成可交互的桑基图D3.js驱动它支持鼠标悬停查看详细信息。结果解读与标注在论文中使用图表时记得在图表标题或图注中清晰说明通路筛选标准如Top 15 by adj.P-value、颜色和大小代表的含义、以及使用的工具ggplot2, ggsankey。10. 总结与下一步通过本教程你已经掌握了用R语言“一套代码生成KEGG气泡图和桑基流向图”的核心流程。这套方法的核心优势在于效率与一致性你只需要维护一份数据清洗和筛选的代码就能驱动两种不同维度的可视化极大减少了出错的可能。最值得尝试的步骤是将第7部分的完整脚本应用到你自己的KEGG富集结果上。你可能会遇到的第一个挑战是数据列名不匹配请根据错误提示耐心调整。第二个挑战可能是桑基图过于复杂请果断使用top_n_pathways和基因频率过滤来简化它。最容易踩的坑是忽略了数据的预处理。务必确保你的输入数据包含Description、p.adjust、Count和geneID或类似这几列并且geneID列中的基因是以特定分隔符如“/”连接的字符串。这是后续所有操作的基础。掌握了这个基本流程后你可以进行多方面的扩展扩展至GO分析尝试将代码适配到GOGene Ontology的富集结果通常需要处理BP、CC、MF三个分类。组合图表使用patchwork或cowplot包将气泡图和桑基图并排组合在一张图中形成更强大的展示效果。自动化报告结合rmarkdown或quarto将数据读取、分析、绘图和结果解读整合成一个可重复生成的HTML或PDF报告。性能优化如果你的数据量极大可以考虑使用data.table包替代dplyr进行数据操作提升处理速度。建议将本文的核心脚本收藏或保存到你的代码库中。下次当你拿到一份KEGG富集结果需要快速产生高质量可视化时它将成为你的得力工具。