☰
pathview实操:把差异基因画到KEGG通路上,从ID映射到在线渲染
2026/9/29 17:12:49 网站建设 项目流程

做转录组差异分析的人,十有八九都有过这种经历:DEG列表出来了,富集分析也做了,但老板一句“这些基因到底在通路里怎么作用的”还是能把你问住。富集结果告诉你哪个通路显著,却没法告诉你某条通路上哪几个节点被激活、哪几个被抑制。我最常用的解决办法,是把核心基因直接画到 KEGG Pathway Map 上——用 pathview 给感兴趣的基因、化合物按表达变化上色,生成一张能直接放进论文或汇报 PPT 里的通路图。

这篇文章既聊 pathview 的在线渲染逻辑,也聊本地 R 环境里的实操细节。内容包括:为什么要上色而不是只看富集列表、输入数据怎么规范才不出乱子、pathview 核心参数逐项拆解、一个转录组联合代谢组的完整案例、以及我在“在线”这件事上的几种落地方式(RStudio Server、网页版 Pathview Web、Shiny 封装)和一路上踩过的坑。不管你是第一次接触 pathview,还是已经能跑通但总被 ID 映射卡住,这份记录应该都能帮你省下不少时间。

1. 为什么给 KEGG 通路的节点上色,而不是只贴一张富集表

1.1 pathview 解决的是哪一类问题

差异基因筛选完成之后,常规操作是跑一个 KEGG 富集分析,拿到一堆显著通路的名字和 p 值。这类结果的问题在于:它把整条通路当成一个“黑盒”,只能告诉你“这条通路相关”,无法区分是通路入口处的受体被激活,还是下游转录因子被抑制。同一个通路里,上游基因上调、下游基因下调,生物学含义截然不同,但富集气泡图看不出来。

pathview 的思路完全不同。它会读取 KEGG 官方发布的 KGML 文件,这个文件保存了通路图的节点坐标、连线关系、基因 ID 和化合物 ID。pathview 把这些信息解析出来之后,把你传入的基因表达值、差异倍数或者化合物丰度映射到对应节点上,按颜色深浅和色相呈现数值大小。最终生成一张和 KEGG 网站上一模一样布局的通路图,只是局部节点被涂上了颜色。

换句话说,pathview 干的事就是“在原版地图上标出你的数据”。它比富集表更强的地方是空间信息:你能看到某条代谢通路的中间产物节点是什么颜色,某个信号通路的磷酸化级联哪一步变化最大。这在多组学联合分析里尤其有用,转录组的基因表达和代谢组的化合物丰度可以同时映射到同一张图上,一眼看出“基因上调导致代谢物积累”这类上下游关系。

1.2 pathview 和 KEGG 官网上色工具、clusterProfiler 的区别

不少人问过我:KEGG 官网本身也提供 pathway mapping 和 color pathway 功能,为什么不直接用网页版,非要绕一圈用 R?我的看法是,官网工具适合小规模手工操作,比如一次看两三个基因;但当你手里有几百个差异基因、化合物数据还要和基因数据叠加时,网页手动提交非常痛苦,而且结果图不可复现。

clusterProfiler 是另一个常被拿来对比的工具。它能做富集分析,也能简单地把基因映射到通路上,但它的可视化重心在富集条形图、气泡图和分类图,不是 KEGG 原生通路布局。pathview 的价值恰恰在“原生布局”四个字:通路图的长相和 KEGG 官网一致,审稿人、老板、合作者看起来没有任何理解成本。

我用一个表把这几个工具的能力差异列出来,方便你按需求选:

工具是否使用KEGG原生通路布局是否支持化合物节点上色是否支持自定义颜色映射批量处理能力可复现性
pathview是是是,非常灵活强,脚本循环即可高
KEGG官网Color Pathway是不支持弱,只能按固定模式差,手动粘贴低
clusterProfiler部分(可输出简单通路图)不支持一般中等中
手画通路模式图否看你怎么画看你怎么画差低

1.3 “在线渲染”在 pathview 语境下的两种含义

标题里提到“在线渲染”,这个词其实可能对应两种完全不同的使用方式,我在后面会分别展开。第一种是底层在线:pathview 运行时需要实时从 KEGG 服务器拉取通路 XML 和底图数据,即使你在本地 R 里运行,它也在“在线干活”。第二种是运行环境在线:你把 R 脚本跑在 RStudio Server、网页版 Pathview Web 或者封装好的 Shiny 应用里,通过浏览器完成整个渲染流程,不需要在一台本地电脑上装 R 环境。

这两种方式我在第 5 章都会给具体操作方案。先把本地跑通的逻辑讲清楚,再谈怎么把它搬到浏览器里,顺序上更合理。

2. 环境准备与数据规范:先把输入数据整清楚

2.1 安装 pathview 和依赖

pathview 是 Bioconductor 家族的包,安装方式和其他 Bioconductor 包一样,不要直接用install.packages(),会装错源。

if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager") BiocManager::install("pathview")

安装过程会拉取一系列依赖包,包括 KEGGgraph、graph、png、Rgraphviz 等。如果公司网络对 Bioconductor 仓库有访问限制,可以先用options("repos")配置国内镜像,或者下载仓库里的压缩包后本地安装。这里提醒一句:pathview 对 R 版本有要求,R 4.2 以下的老版本可能会出现依赖冲突,建议至少用 R 4.2 以上。

装完之后验证一下:

library(pathview) packageVersion("pathview")

能正常打印版本号,说明环境没问题。后面所有代码都假设你已经跑通了这一步。

2.2 物种代码与基因 ID:最容易错的第一步

KEGG 对每个物种有一个三字母代码,例如人的是hsa、小鼠是mmu、大鼠是rno、斑马鱼是dre、拟南芥是ath、酿酒酵母是sce、大肠杆菌是eco。这个代码必须和你的数据物种严格对应,否则 pathview 会下载到错误的通路文件,结果所有基因都匹配不上。

物种代码还有一个容易踩的坑:同一个物种可能有多个别名。比如人有时候被写成human,鼠被写成mouse,pathview 只认官方三字母代码,写错直接报错或者全灰。常见物种对照表我给你列一份:

物种拉丁名KEGG代码
人Homo sapienshsa
小鼠Mus musculusmmu
大鼠Rattus norvegicusrno
斑马鱼Danio reriodre
果蝇Drosophila melanogasterdme
线虫Caenorhabditis eleganscel
拟南芥Arabidopsis thalianaath
酵母Saccharomyces cerevisiaesce
大肠杆菌Escherichia coli K-12eco

关于基因 ID,这是绝大多数人第一次跑 pathview 失败的根本原因。KEGG 通路节点里的基因使用 Entrez Gene ID,比如人类的 TP53 在 KEGG 里对应hsa:7157。如果你手里的数据是基因符号(Symbol)或者 Ensembl ID,必须先把它们转换成 Entrez ID,再传给 pathview。否则 pathview 拿 symbol 去匹配 KEGG 节点 ID,自然会全部落空。

转换用org.Hs.eg.db(人或相应物种的 org.db 包)或者 AnnotationDbi 的mapIds函数就行。我一般在读取差异基因表格之后立刻做这一步,不拖到 pathview 调用时再处理:

library(org.Hs.eg.db) library(AnnotationDbi) deg$entrez <- mapIds(org.Hs.eg.db, keys = deg$symbol, keytype = "SYMBOL", column = "ENTREZID", multiVals = "first")

这里需要注意multiVals = "first",当某个 symbol 对应多个 Entrez ID 时只取第一个。这个设定在大部分情况下够用,但如果你想严谨一点,可以先把多对一的基因单独拎出来人工确认,一般数量很少,不影响整体通路渲染。

2.3 化合物 ID:C 开头不代表格式对了

KEGG 化合物 ID 统一是C加五位数字,比如水是C00001,ATP 是C00002。但现实中代谢组学数据最常见的输入格式是 HMDB ID、ChEBI ID 或者商品化代谢物数据库自己的编号。直接把数据框里写着HMDB0000054的 ID 传给 pathview,它不能识别。

解决办法是用 KEGGREST 包进行 ID 转换。KEGGREST 是 Kyoto Encyclopedia of Genes and Genomes 官方 REST API 的 R 客户端,keggConv 函数专门做这种 ID 对照:

library(KEGGREST) keggConv("compound", "hmdb:HMDB0000054")

返回结果类似cpd:C00005,你只需要把cpd:前缀去掉,剩下的C00005就是 pathview 能认的化合物 ID。如果是一次性转换几十上百个化合物,写个循环或者用vapply批量处理即可。这里有个实际操作经验:有些 HMDB ID 在 KEGG 里没有对应条目,返回结果是 NA,不要强行剔除,先记录下来,后面看这些化合物在通路里是不是本来就该缺失。通常 KEGG 没有收录的代谢物比例不高,不影响整体图面。

3. pathview() 函数核心参数拆解

3.1 基因数据的颜色映射逻辑

pathview 的核心调用就是pathview()这个函数,它有一套默认规则:当你传入的基因数据是像 log2FC 这样含负值的向量时,正值默认显示为红色,负值显示为绿色,绝对值越大颜色越深。灰色表示该基因在数据里没有对应值,或者表达变化接近零。

这里有一个新手经常搞混的点:pathview 的颜色深浅不是自动正则化到数据范围的,而是由limit和bins两个参数控制的。limit决定颜色映射的数值边界,bins决定把数值范围分成多少档。比如:

limit = list(gene = 2, cpd = 2) bins = list(gene = 10, cpd = 10)

含义是:log2FC 绝对值超过 2 的,都按最高饱和度上色;从 0 到 2 之间切 10 档,颜色逐级变深。如果你不设limit,pathview 会取数据的最大绝对值作为边界,这会导致极端离群值把所有其他节点都压成很浅的颜色,图面上看起来只有一两个点有颜色。我建议在跑正式分析前先看数据的分布:

summary(abs(deg$log2FC))

中位数、90 分位数能帮你判断合理的limit应该设多少。一般转录组 log2FC 的 cutoff 是 1,我会把 limit 设在 1 或 2,而不是让离群值主导色阶。

另一个重要参数是both.dirs。当数据里有正有负时,默认both.dirs = list(gene = TRUE, cpd = TRUE),表示正负需要分开映射。如果你传的数据是只有正值的东西,比如 P 值的负对数,那就应该设both.dirs = FALSE,让所有节点都在同一条色阶上。

3.2 化合物数据和基因数据同时渲染

pathview 最吸引人的一点是可以把gene.data和cpd.data同时传入。这在代谢通路、比如中心碳代谢(hsa01100)里特别有用:基因的转录水平变化和代谢物的丰度变化同时上色,能直接判断转录调控和代谢表型是否一致。

gene_data <- c("7157" = 1.8, "1956" = -1.2, "207" = 0.9) cpd_data <- c("C00001" = 0.5, "C00002" = -1.5) pathview(gene.data = gene_data, cpd.data = cpd_data, pathway.id = "hsa01100", species = "hsa", same.layer = TRUE)

same.layer参数控制基因和化合物是否画在同一张图里。TRUE时输出一张综合图,方便看上下游;FALSE时会生成两张独立的图,一张只有基因、一张只有化合物。实际使用中我的建议是:先跑same.layer = TRUE看整体,如果图面太拥挤、节点标签重叠严重,再分别出图。

化合物 ID 匹配不上时,pathview 并不会直接报错,而是默默跳过未匹配的节点,导致你看到的图里有一部分化合物没有上色。这个时候要回头检查第 2.3 节讲过的 ID 转换环节,确认 KEGG compound ID 格式正确。

3.3 图形与布局参数:从能出图到出好图

pathview 默认输出png、svg和xml三种格式。png 是位图,方便直接插入 Word 和 PPT;svg 是矢量图,可以用 Inkscape 或 Adobe Illustrator 继续编辑,期刊投稿一般要求矢量图。文件名格式是hsaXXXXX.pathview.png,其中XXXXX是通路编号。

out.suffix参数用来给输出文件加后缀。例如:

out.suffix = "DEG_treatment"

生成的文件就是hsa04668.pathview.DEG_treatment.png。这个参数在跑多组比较时特别重要:不加后缀的话,第二次调用会直接覆盖第一次的结果图。

节点大小和标签位置也有参数可以调,但我的经验是默认值基本够用,真正需要调的是key.pos。图例默认放在图的右侧,个别通路节点比较密的时候,图例会挡住节点,此时设key.pos = "topright"或者key.pos = "bottomleft"可以缓解。如果你需要把通路图拼成大图,key.pos = "topright"通常更方便对齐。

还有一个细节:node.sum参数。当多个基因映射到同一个 KEGG 节点时,pathview 默认对这些基因的数值求和(sum)。这意味着如果节点上有三个基因,一个上调 2、一个下调 1、一个不变,求和结果是 1,看起来是轻度上调。换用mean会取平均,得到 0.33,更接近整体趋势。到底用哪个,取决于你想强调总量还是均值,我习惯用mean,不容易被基因数量误导。

4. 一个完整的实战案例

4.1 数据读取与 ID 转换

假设你手上有这样一张差异基因表deg_table.csv,包含两列:symbol和log2FC。我们想要把这张表渲染到 TNF 信号通路上。

library(pathview) library(org.Hs.eg.db) library(AnnotationDbi) deg <- read.csv("deg_table.csv", stringsAsFactors = FALSE) head(deg) deg$entrez <- mapIds(org.Hs.eg.db, keys = deg$symbol, keytype = "SYMBOL", column = "ENTREZID", multiVals = "first") cat("转换后丢失的基因数:", sum(is.na(deg$entrez)), "\n") deg <- deg[!is.na(deg$entrez), ] gene_data <- deg$log2FC names(gene_data) <- deg$entrez gene_data <- gene_data[!duplicated(names(gene_data))]

这段代码里有几个容易忽略的细节。首先是cat输出丢失基因数,如果丢失比例超过 10%,说明你下载的注释包可能跟数据物种不匹配,先停下来检查,不要一路往下跑。其次是!duplicated(names(gene_data)),多个 symbol 转到同一个 Entrez ID 时,如果不做去重,pathview 会因为重复名字报错或者警告。

4.2 单条通路渲染与结果文件解读

拿到规范的命名向量之后,直接调pathview():

pv <- pathview(gene.data = gene_data, pathway.id = "hsa04668", species = "hsa", out.suffix = "TNF_deg", limit = list(gene = 1, cpd = 1), bins = list(gene = 10, cpd = 10), both.dirs = list(gene = TRUE, cpd = TRUE))

运行之后,工作目录会多出这些文件:hsa04668.pathview.TNF_deg.png、hsa04668.pathview.TNF_deg.svg、hsa04668.pathview.TNF_deg.xml,可能还有hsa04668.png和hsa04668.xml这两个中间文件,是 pathview 下载的原始 KEGG 底图数据。

打开 png 你会看到,TNF 通路图上的大部分基因节点都有颜色。红色越深的节点表示该基因在实验组里上调越显著,绿色越深表示下调越显著,灰白色节点表示样本里没有检测到这个基因。观察这类图有一个技巧:先看通路入口处的配体和受体节点,再看中间信号分子,最后看末端转录因子和效应分子。如果红色集中在入口处而末端是绿色,说明信号传导被阻断;如果全线都是红色,说明整个通路被激活。

中间文件千万别急着删。hsa04668.xml就是 KGML 数据,后面排查 ID 匹配、做缓存、二次编辑图都要用到它。

4.3 批量渲染多条通路

实际项目里很少只看一条通路。你通常是从富集分析结果里拿到一串显著通路 ID,然后想一次把所有这些通路都跑一遍。推荐做法是先构造通路 ID 向量,再循环调用 pathview:

pathways <- c("hsa04668", "hsa04010", "hsa04210", "hsa04630") for (pid in pathways) { pathview(gene.data = gene_data, pathway.id = pid, species = "hsa", out.suffix = "TNF_deg", limit = list(gene = 1, cpd = 1), bins = list(gene = 10, cpd = 10), both.dirs = list(gene = TRUE, cpd = TRUE)) }

循环里加一个进度提示会更友好,比如每跑完一条通路打印一次当前文件名。还有一点实操建议:如果通路数量多,同一组颜色参数保持统一,最后拼图时才不会因为色阶不一致被审稿人挑毛病。

批量跑完之后,可以用系统命令把同一后缀的 png 文件归到子目录,方便后续统一查看。如果你用的是 RStudio Server,RStudio 的 Files 面板可以直接预览 png,这就已经是最轻量的“在线查看通路图”方式了。

5. 在线渲染的几种落地方式

5.1 RStudio Server:浏览器里的完整 R 环境

如果你的工作流里需要频繁更新数据、反复跑通路渲染,最实用的在线方案是 RStudio Server。它在服务器上运行完整的 RStudio 界面,你本机只需要一个浏览器,不需要安装 R。

部署方式不展开讲,核心逻辑就是把 R 脚本和差异基因表放到服务器上,浏览器访问 RStudio Server 的地址,然后在网页里运行脚本。RStudio 的 Viewer 面板能直接预览 png,Files 面板能快速查看所有输出文件。对于团队协作场景,我给同事的流程是:他们只改deg_table.csv,我负责的渲染脚本里把gene_data的构造写成了可复用函数,跑一条命令就能更新全部通路图。

这里有个小坑:RStudio Server 默认工作目录是启动 R 时所在的目录,新人经常遇到“文件明明上传了但读不到”的问题。建议在脚本开头专治各种不服:

setwd("/path/to/your/project")

5.2 Pathview Web 网页工具的操作流程

如果你完全不想碰代码,只想偶尔上一次色,可以直接用 Pathview Web 这个网页版工具。它的界面是传统的表单式提交:

  1. 选择物种代码,比如 hsa。
  2. 填写或粘贴基因列表,每行一个基因,可以带数值列(表达量或 log2FC)。
  3. 选择目标通路,按通路 ID 或者浏览分类目录选择。
  4. 提交后等待服务器渲染,页面会直接展示着色后的通路图,并提供下载链接。

网页版对输入格式比较严格,基因 ID 建议直接用 Entrez Gene ID,避免它在线上帮你做 ID 映射时因为符号歧义而失败。数据量方面,提交几千个基因没有问题,但如果同时开启多通路批量渲染,服务器响应时间会明显变长,耐心等即可。

这种方式的局限也很明显:个性化参数少,不能调颜色阈值、不能叠加化合物数据、批量渲染能力弱。所以我的建议是把它当成“应急工具”,正式项目还是用 R 脚本,方便复现。

5.3 把 pathview 封装成 Shiny 小工具

如果团队里有人不写 R,但你希望他也能自助渲染通路图,可以考虑用 Shiny 把 pathview 包一层网页界面。这个方案在生信支撑类的场景下特别实用,分析人员把差异基因上传到网页,选择物种和通路,点一下按钮就出图。

一个最简版本的 Shiny app,核心逻辑就是文件上传控件 + pathview 调用 + 图片展示。服务端代码大致是这个样子:

library(shiny) library(pathview) server <- function(input, output, session) { observeEvent(input$run, { req(input$degFile) deg <- read.csv(input$degFile$datapath, stringsAsFactors = FALSE) gene_data <- deg$log2FC names(gene_data) <- deg$entrez pathview(gene.data = gene_data, pathway.id = input$pathway, species = "hsa", out.suffix = "shiny_user") output$pathPlot <- renderImage({ list(src = paste0(input$pathway, ".pathview.shiny_user.png"), contentType = "image/png") }, deleteFile = FALSE) }) }

界面端就放一个文件上传框、一个通路 ID 输入框、一个“渲染”按钮和一个图片输出区域。这套东西看起来简单,实际用起来能省掉大量“帮我看一下这个基因在哪个通路有变化”的重复请求。我在项目里甚至加了一个下拉框让用户选择颜色主题,这样不同组的人出图风格也能统一。

6. 我踩过的几个坑

6.1 KEGG 数据下载失败与缓存策略

pathview 第一次跑某条通路时,需要从 KEGG 服务器下载该通路的 XML 和 png 底图。如果服务器网络不稳定,或者访问外网受限,会直接报出类似 “cannot open URL” 的错误。

我的处理方法是提前手动下载缓存文件,把网络请求放到可控的时间点做。只要你能访问 KEGG 公共数据库,可以先用 KEGGREST 把 KGML 文件拉下来:

library(KEGGREST) kgml_file <- keggGet("hsa04668", "kgml") writeLines(kgml_file, "hsa04668.xml")

png 底图可以用download.file()从 KEGG 的网站地址拉取。把这两个文件放在工作目录之后,再运行 pathview 时它就能直接使用本地文件,不再依赖网络。这样既保证了离线环境下也能出图,也避免在正式分析跑到一半时因为网络抖动中断。

关于“在线渲染”这个词,这里正好说清楚:pathview 的核心逻辑可以完全本地跑,但首次获取 KEGG 数据这一步确实依赖网络。所以最好的实践是“在线拉数据、离线出图”,而不是反复在线。

6.2 基因节点全是灰色:ID 映射失败的排查链路

这是 pathview 新手最常遇到的问题,没有之一。现象是通路图完整,但所有基因节点都是灰白色,一个颜色都没有。

我的排查顺序一定是这样的:

第一步,检查names(gene_data)的前几个值:

head(names(gene_data))

如果看到的是 “TP53”“EGFR” 这种 symbol,问题就找到了——pathview 默认只认 Entrez Gene ID。需要回到第 2.2 节做 mapIds 转换。

第二步,如果已经是数字 ID,打开hsa04668.xml或者用readLines()看一下 KEGG 通路里实际存在的 gene entry ID 是不是hsa:7157这种格式,检查是不是混入了不该有的物种前缀。names(gene_data)里应该是7157,而不是hsa:7157。

第三步,检查是否因为multiVals = "first"导致了大量 NA。如果 mapIds 之后丢失的基因数量异常,考虑是不是注释包装错了。用人的数据却装了小鼠的 org.Mm.eg.db,这种低级错误我也犯过。

另外一个隐藏坑是:gene_data的数值全是 NA。如果 read.csv 之后没仔细看数据,log2FC列里可能混入了 "NA" 字符串,导致向量里的元素都是 NA。pathview 拿到全 NA 数据也不会报错,结果就是全灰。跑之前一定要加一句:

stopifnot(!all(is.na(gene_data)))

6.3 化合物和基因图层错位

同时传入gene.data和cpd.data时,一个常见的现象是基因上色正常,但化合物节点没颜色,或者化合物节点上有个奇怪的灰色标签。这里的原因多半出在化合物 ID 的格式上。

前面说过,KEGG 化合物的标准 ID 是C00001这种。如果你的names(cpd_data)里写成了带cpd:前缀的格式,pathview 在解析的时候不会自动去掉前缀,就会全部匹配失败。我自己写过一个很隐蔽的 bug:用keggConv转换 HMDB ID 后,直接拿返回结果去当names,忘了去掉cpd:,然后排查了整整一个下午。

正确做法:

converted <- keggConv("compound", "hmdb:HMDB0000054") cpd_id <- sub("^cpd:", "", converted)

批量转换时也保持同样的处理,确保最终传给 pathview 的 names 是裸的C加五位数字。

6.4 多通路批量渲染时的内存与输出管理

批量跑二三十条通路的时候,pathview 每次都会把 KGML 解析成图形对象,内存占用上升很快。如果机器内存不够,跑到一半会莫名中断,或者图片输出不完整。

我的常规做法是在循环里用显式垃圾回收:

for (pid in pathways) { cat("Processing", pid, "\n") pathview(gene.data = gene_data, pathway.id = pid, species = "hsa", out.suffix = "batch") gc() }

另外,批量跑之前先确认输出目录有写入权限。服务器上如果路径不对,pathview 会把图片写到“你以为的目录”和“实际的工作目录”不一致的地方,最后整理图片时找不到文件,也是常见情况。脚本开头用getwd()打印一下当前目录,至少能少一个定位问题的维度。

最后说几句个人体会

pathview 我用到现在,最大的感受是它的学习曲线其实很短。一次性把 ID 规范、参数含义和输出文件搞明白之后,剩下的就是从“能出图”到“出好图”的打磨过程。我给自己的标准是:每张图拿到手里,你要能说清楚为什么红色集中在上游、为什么某个化合物没有上色、为什么这组参数配这个通路最合适。能解释清楚这三个“为什么”,这个工具才算真正属于你了。

如果你刚开始接触,建议先挑一条你熟悉信号通路,比如 TNF 或者 MAPK,用差异基因列表跑通一次全流程,再逐步叠加化合物数据、批量渲染和在线部署。这条路我已经走过一遍,坑也都替你踩完了,按这篇文章的顺序做,出图是迟早的事。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询