
做代谢组学的朋友应该都有过这种体验LC-MS上机跑完peak table导出来差异代谢物也筛出来了论文的第一张热图也画完了但往下一做就卡壳。你要解释这批差异代谢物到底是从哪来的哪些代谢物在生物学上是协同变化的一个功能单元结果要么回到HMDB一条条人工翻“Sources”字段要么拿起Cytoscape手动连线再要么就是复制一堆表格到Excel里熬到凌晨。TidyMass2这个工作我是在Nature Communications上看到的它把这件让很多人头疼的事也就是代谢物溯源和功能模块分析做成了一条完整且可复现的R流水线。它最打动我的地方不是又造了一个新的组学算法而是用tidy data的数据哲学把代谢组学数据从预处理到溯源、再到共丰度网络和功能注释全都统一在一个框架里。所以这篇就想从实际使用的角度聊聊TidyMass2解决了哪些真实痛点也把我实际操作过程中踩过的坑、调过的参数、怎么和下游统计对接的经验一并整理出来给正在处理代谢组数据但卡在“溯源”和“模块”这两个环节的人一个可以直接参考的路线。如果你是做临床代谢组、微生物-宿主共代谢、暴露组或者只是想摆脱散装工具链、把分析做得更规范的人这篇应该对你有用。1. 代谢组学跑完流程之后溯源和功能模块为什么让人头疼1.1 代谢物溯源难在“来源”不是一个简单字段所谓代谢物溯源听起来简单就是回答“这个代谢物从哪来”嘛。但实际做起来你会发现“从哪来”有很多层含义。最常见的是区分内源性和外源性比如一个代谢物是宿主自身合成的还是来自饮食、肠道微生物、甚至是环境污染暴露。这个区分在机制研究里特别重要不然你发现一个差异代谢物根本说不清楚它是机体变化的因还是果。麻烦在于同一个代谢物可能同时具备多个来源。比如某些胆汁酸既可以是宿主肝脏合成的初级胆汁酸又可以在肠道里被微生物转化溯源结果天然就不是非此即彼。HMDB里确实每个化合物都有来源描述但那是自由文本写着“microbiota”“dietary”这种非结构化的词你没法直接拿来做统计更没法回答“这个分子有多少可能来自肠道菌群”这种问题。再加上HMDB ID、KEGG ID、PubChem CID之间的ID体系不统一很多代谢物名字里还带着加合物、异构体标记做一次人工追溯几十个代谢物能翻一整天。另一种溯源需求是空间或组织来源。比如外周血里检测到的某类代谢物到底来自肝脏还是肠道还是脂肪组织单靠血液样本很难判断通常需要结合组织特异性的标志物或同位素示踪。这类需求TidyMass2并不能完全自动化解决但它的可扩展设计可以让你把已经验证过的组织来源信息作为额外的注释层加进去。也就是说它能帮你在同一套数据结构里管理“数据库匹配来源”和“人工实验来源”而不是让你分开维护两份Excel。1.2 功能模块分析难在“单个代谢物说了不算”做差异代谢物筛完之后你手上可能有一百多个显著变化的代谢物按字母排一列每个看起来都很有道理但放一起讲不出故事。常规做法是拿去做KEGG通路富集但代谢组学和转录组有个明显区别许多代谢物并不在KEGG通路里尤其是一些脂质、植物来源成分和环境暴露物富集结果往往只覆盖一小部分差异代谢物剩下的就成了“孤儿”写论文时只能选择性忽略。功能模块分析想解决的问题就是不依赖通路注释直接从数据里找“哪些代谢物总是协同变化”。如果一组代谢物的丰度在样本间高度相关它们很可能来自同一个生物过程、同一个来源、受到同一类调控。把这些代谢物聚成一个模块再用模块的整体丰度去和表型关联比单个代谢物做几百次统计检验更稳也更符合生物学直觉。道理我都懂可过去想做这件事路径特别曲折。你要先自己写脚本算相关矩阵然后导出成网络文件再拿去Cytoscape里跑MCODE或者Louvain聚类聚完类再回到R里计算模块特征值最后逐个模块做富集注释。每一步之间都要来回导表格、改格式不同软件对ID字段的要求还不一样稍微哪个环节清洗不到位结果就接不上。WGCNA虽然能做类似的事但它毕竟是转录组数据的逻辑在代谢物上需要手动调一堆参数做得人很头大。1.3 最关键的问题其实是格式割裂源头溯要在HMDB里翻模块网络要在Cytoscape里画富集注释又要去MetaboAnalyst里跑统计关联还得自己在R里写数据格式完全割裂。TidyMass2的切入点我认为特别准它没有花大力气设计新的数学方法而是先把整个分析过程中的数据结构统一成tidy的tibble和统一的对象类型。这一步看似基础实际上解决了代谢组学分析里最多的时间黑洞。数据结构统一之后溯源结果可以直接用于模块分组模块分组又可以无缝对接富集注释和表型关联因为所有中间结果都是同一套对象派生出来的不再需要手工bridge。2. TidyMass2到底做了什么设计让这两个难题能串起来2.1 Tidy data“每个变量一列、每个观测一行”为什么救命TidyMass2的设计核心是tidy data原则。这个概念来自tidyverse说穿了就三条每个变量是一列每个观测是一行每个表是一个主题。听起来像废话但代谢组学数据常年违反这个规则。比如仪器导出的原始表格通常是宽表结构一行一个feature后面几百列是样本的峰面积。这其实既不是“一个观测一行”也不好做后续分析。你如果想按分组筛选样本、按时间点重排就得不停地把宽表转长表、把长表再转宽表。TidyMass2把数据按三个实体分开管理一个是feature级别的变量信息包含mz、rt、注释的HMDB ID和KEGG ID一个是sample级别的样本信息包含分组、批次、临床变量一个是feature-sample对级别的定量表也就是峰面积。主对象把这几个表聚合成一个整体所有分析都围绕这个对象展开。这个设计带来的实际好处是你可以大胆执行filter、mutate、select这类dplyr操作而不用担心把feature表和sample表对应乱也不会出现“我删掉了三个样本但后面所有统计还是把它们算进去了”这种低级错误。我自己的感受是用了这套数据结构之后代码的可复用性明显提升换一个数据集时只需要重新准备三个表后面的溯源和模块流程基本上原样跑通。2.2 溯源模块把“来源”从文本变成了可计算的结构化标签TidyMass2对溯源的实现本质上是建立了一个代谢物与来源类别之间的可扩展映射关系。它先把来自HMDB、KEGG、文献报道的来源信息整理成结构化的source-annotation三元组也就是代谢物ID、来源类别、证据等级这三项。来源类别至少包括内源性、饮食来源、微生物来源、外源暴露来源等证据等级则用来区分这个来源是数据库直接标注的还是文献报道支持的避免把不同可信度的信息混在一起讲。分析的时候你拿做好的峰表里的代谢物注释去匹配这个来源映射表。匹配字段可以是代谢物名称、HMDB ID、KEGG ID甚至InChIKey也可以。匹配完成后每个代谢物会带上一组来源标签和证据等级。因为你全程是在tidy数据框架里操作溯源结果可以直接按样本分组做统计比如算一算某个差异模块里有多少比例代谢物是微生物来源还可以做卡方检验比较不同分组间的来源结构差异。我觉得这个逻辑设计得务实的地方在于它没有试图给每个代谢物强行指定唯一来源而是保留“一个代谢物可能有多个来源”的真实性把判断的主动权交给研究者。毕竟胆汁酸的例子摆在那里强行去重反而会丢掉生物学信息。2.3 功能模块分析网络、聚类、富集、关联的四步流水线功能模块分析部分TidyMass2的管线思路很像把WGCNA和Cytoscape的常用流程做成了R原生的四步流水线。第一步是构建共丰度网络。它计算所有代谢物两两之间的相关矩阵通常用Spearman相关系数对非正态分布更稳健然后设定相关系数阈值和p值阈值只保留显著相关的边。第二步是模块检测用Louvain这类图聚类算法把网络划分成若干个模块每个模块里的代谢物理论上共享相似的丰度变化模式。第三步是模块富集注释把每个模块里的代谢物映射到KEGG通路和HMDB化合物类别看这个模块在功能上是不是富集到某种特定通路。第四步是模块与表型的关联分析计算每个模块的特征值用这个整体值去和你的临床分组或连续变量做相关分析。这四步单拆开没有什么惊人的创新但TidyMass2把它们放在同一个对象生命周期里最大程度减少了转换和对接成本。而且因为溯源结果也是同一套对象的一部分你会很自然地做到一步操作先把差异代谢物做溯源再对同一个对象做模块分析最后直接看某个模块里到底富集了哪个来源的代谢物。溯源和功能模块这两个原本割裂的分析就这么串成了完整链条。2.4 为什么是把宝压在R而不是Python或网页工具上有人可能会问同样的事情用Python的sckit-learn加networkx也能做为什么要用R包我的理解是代谢组学质谱数据的上游处理在R里太成熟了。XCMS、mzR、BiocParallel这些工具已经把峰检测、峰对齐、注释的生态做得很完整TidyMass2长在R生态里可以直接和上游数据无缝衔接不需要在Python和R之间来回搬运文件。再有就是tidyverse在数据处理上的表达能力确实强表格操作顺手而且整条分析流程里的每一步都是显式、可审计的代码对论文方法学部分非常友好。R确实有内存管理的短板但现在通过data.table、disk.frame这类工具也能缓解我下面会专门聊性能优化的做法。整体来说对代谢组学这个领域的研究者R的学习曲线比Python要友好不少TidyMass2押注R生态是聪明的选择。3. 实操从峰表到溯源再到功能模块分析完整走一遍3.1 环境准备和数据组织决定后面顺不顺先说安装。TidyMass2基于R 4.2以上版本和tidyverse生态安装时用GitHub或CRAN渠道都行。因为依赖里包括一些Bioconductor的包我建议先装好BiocManager再装TidyMass2否则依赖解析容易出问题。install.packages(BiocManager) BiocManager::install(c(tidyverse, data.table)) # 以下命令以项目文档为准不同版本可能调整 remotes::install_github(TidyMass/TidyMass2)装好之后最重要的不是急着跑函数而是把数据整理成合适的结构。我自己常用的方式是准备三个表。第一是feature信息表每行是一个代谢物特征至少包含feature_id、mz、rt、代谢物名称、HMDB ID、KEGG ID。第二是样本信息表每行是一个样本包含sample_id、分组信息、批次、协变量。第三是定量表我推荐用长表结构每行是一个样本对某个feature的定量值列至少是feature_id、sample_id、intensity。如果你的数据目前是宽表也就是一行feature后面挂几百个样本列先用tidyr的pivot_longer转成长表再说。这一步看起来繁琐但整理完之后后面所有分析都会变得非常简单值得多花十分钟。主对象构建的代码大致长这样library(TidyMass2) obj - create_tidymass_object( peak_table peak_table, # 长表feature_id, sample_id, intensity feature_info feature_info, # feature级别的注释 sample_info sample_info # 样本级别的分组和临床信息 )3.2 代谢物溯源跑之前先把这些坑填平对象构建好之后溯源操作本身很快几行代码就出结果obj - trace_metabolite_origin( obj, source_levels c(endogenous, dietary, microbial, xenobiotic), evidence c(database, literature, prediction) ) origin_table - extract_origin_table(obj)这个结果表里通常每个代谢物一行或多行记录它命中了哪个来源类别、证据等级是多少、依据来自哪个数据库。我建议你拿到这个表的第一件事不是急着统计而是先随检一下你的代谢物名称质量。实际操作中常见的翻车原因就是feature表里的代谢物名字还挂在“XXX_ESI_POS”这种加合物后缀或者带了“(E)”这种立体异构标记匹配时很容易失败。碰到这种情况提前做一个名称清理函数把盐加合物、异构体后缀去掉会显著提高匹配率。统计来源比例时我推荐按“证据等级”分层报告。比如可以说“在差异代谢物中有75%可以溯源到至少一个已知来源其中数据库直接标注的比例约占40%”。如果某个来源类是“prediction”预测出来的比例很高写论文时要谨慎不要和实验验证混为一谈。溯源的核心价值是帮你说故事但前提是别把故事说得超过证据能支撑的程度。3.3 共丰度网络构建与模块发现参数要跑出“稳定感”模块分析的代码同样直接net_obj - build_coabundance_network( obj, method spearman, cor_cutoff 0.6, p_cutoff 0.05 ) net_obj - detect_modules(net_obj, algorithm louvain) module_table - get_module_table(net_obj)参数上最常见的争议是相关阈值选多少。我做过不少次之后倾向于这样的原则先看你的样本量样本量小于50时0.6的Spearman相关系数其实已经算比较强的信号了可以先用0.6跑一版看整体网络规模和模块数量的分布。如果模块数大于10个可能阈值偏低网络太碎试着把cor_cutoff往上提到0.7如果模块数只有两三个大模块里塞了上百个代谢物说明阈值太高或者说网络过度聚合可以往下调到0.5再试试。没有绝对正确的阈值能解释得通、模块结果稳定最重要。模块检测算法方面我通常先跑Louvain因为它在生物网络上的表现比较稳健速度快社区结构找得好。如果你发现Louvain结果不稳定可以再加一种算法比如walktrap做交叉验证两种算法都落在同一模块里的代谢物作为这个模块的“核心成员”做下游富集时优先看核心成员。这个办法能显著减少模块结果的抖动写论文时也更经得起审稿人追问。富集注释和表型关联的代码如下net_obj - enrich_modules( net_obj, pathway_db KEGG, organism hsa, pvalue_cutoff 0.05 ) enrich_result - get_enrichment_table(net_obj) mod_trait - correlate_module_phenotype( net_obj, phenotype disease_status, method spearman )3.4 可视化与联合解读溯源和模块是互补的分析做完可视化是最后一步好的图也直接决定你能不能说服别人。TidyMass2在这块基于ggplot2封装了一些默认出图函数但我的建议是你可以用自己的语法基于模块表去画更精细的图。无非是几张固定套路图模块-样本丰度热图、模块网络关系图、模块富集的点图以及模块特征值与指标的相关性图。真正出彩的是把溯源结果嵌到模块图里。我举个例子之前我做一个肠道菌群相关项目跑完模块分析发现模块3里有十几个代谢物单独看KEGG富集结果很平淡只是“bile acid biosynthesis”这类词。但当我把溯源表关联到模块成员上时才发现模块3里70%的代谢物都带“microbial”来源标签再结合文献一查这群代谢物正好是肠道菌群修饰胆汁酸的典型产物。有了“模块3是微生物来源的功能模块”这个结论整个故事一下子就立住了后面再做多组学验证也有了明确切入点。这就是我一直强调的观点溯源单独做只是给每个代谢物贴标签功能模块单独做只是把代谢物分成组只有把两者联合起来看才真正回答“这些代谢物从哪里来、为什么一起变、可能参与什么功能”的完整问题。4. 实际运行中容易踩的坑以及我常用的排查方法4.1 溯源匹配率太低先查这四件事匹配率低是溯源模块最高频的抱怨。如果你拿到结果发现一大半代谢物没有来源标签先别急着怀疑工具按顺序排查第一代谢物名字是否标准化去掉加合物、盐类、立体化学后缀第二注释ID是否填对了HMDB ID和KEGG ID不要混用先用主ID匹配再用名字兜底第三是否有大量同位素标记或碎片离子混在峰表里这类特征本身没有生物学来源匹配不上是正常的第四数据库中确实没有覆盖比如一些新型污染物或结构全新的天然产物这类代谢物溯源不上不是bug是领域现有知识的边界。我建议把匹配率预期设在40%到70%之间超过70%已经算很好。低于40%时再重点检查前面四个原因做过一轮清理后通常能提升十几个百分点。写论文时明确说“在可注释代谢物中有XX%可溯源”别拿所有features当分母不然数字很难看逻辑也不公平。4.2 模块划分不稳定用“多算法取交集”来兜底如果你连续跑两次模块划分结果差异明显多半是相关矩阵本身噪声大或者样本量不够支撑稳健的图聚类。有三个办法依次试。第一步对共丰度网络做bootstrap重采样只保留在大多数重采样中依然显著的边类似网络剪枝第二步同时跑Louvain和walktrap两种算法保留两种算法一致分组一致的核心成员第三步看核心成员的个数是否足够做下游富集如果核心成员少于五个说明这个模块本身信号弱建议提高相关阈值重跑。还有一个小技巧可以先不管p值直接观察相关性直方图。如果大多数代谢物对之间的相关性在0.3左右徘徊网络本身就缺乏清晰结构这时候任何阈值下的模块划分都不会有生物学意义。这种情况我会先回到数据清洗环节确认是否有批次效应或峰对齐质量不足的问题。4.3 大矩阵计算卡到爆优先过滤再算相关几千个feature和几百个样本做两两相关矩阵在R里内存消耗很大。我遇到过跑了好几分钟还没出结果、最后直接内存溢出的情况。现在我的处理顺序是先过滤低丰度特征比如在任何样本里都没有稳定检出的直接删掉再过滤低变异特征用组内变异系数或标准差排序只保留变异够大的代谢物完成过滤后features数量通常会从几千降到几百算相关矩阵的压力会小一个量级。如果过滤后矩阵还是大就把数据转成稀疏矩阵再算相关或者对feature分块计算边算边合并。TidyMass2的代码里一般会封装并行计算但你自己也可以提前用BiocParallel启动多个核能快不少。4.4 结果接不上下游统计就靠统一对象兜着最后也是我最看重的是TidyMass2输出的结果和下游其他统计方法怎么对接。我验证下来最稳妥的出口格式就是把模块特征值、代谢物来源标签和分组信息导出成一张干净的长表然后你想做logistic回归、生存分析、随机森林都行。因为模块特征值本身是每个样本一行、每个模块一列很容易merge回临床表。这里唯一要注意的是别把“溯源比例”“模块特征值”当成原始变量直接丢进模型而不做标准化。模块特征值默认是中心化和标准化的结果解释系数时要说清楚单位。溯源比例作为百分比数据进回归前建议做logit转换直接放百分比容易违反模型假设。5. 一些使用体会以及还能往哪个方向扩展5.1 这套工作流帮我节省了最多时间的地方我最直接的感受是TidyMass2省下的不是某个单步操作的时间而是“数据在不同工具之间搬运”的时间。过去做一个溯源加模块分析的完整流程我手动在HMDB、Cytoscape、MetaboAnalyst之间倒腾数据至少要两天现在在同一个R环境里跑一两个小时就能出全部结果表和初版图。最宝贵的是整条流程是对操作系统级的可复现不会因为某一步在网页工具里点了不同选项导致结果对不上。5.2 最值得关注的是“溯源模块”的联合分析姿势如果你只把TidyMass2当成一个替代Cytoscape的网络工具那有点亏。我更推荐的做法是在拿到溯源表后不要停马上把它作为模块特征的注释层融合进网络图里。比如把不同的来源类别赋予不同颜色叠加在模块网络上你几乎一眼就能看出哪些模块是“菌群来源的功能单元”哪些模块是“宿主内源代谢物聚集体”。这个图讲故事的效率比单看通路富集高很多尤其面对偏临床背景的读者时特别管用。5.3 未来如果和转录组做多组学联合怎么接还有一个值得提前布局的方向是如果你手上有同一个批次的转录组或宏基因组数据可以把TidyMass2得到的代谢物模块特征值和基因模块、菌群丰度做跨组学相关分析。原理也很朴素既然代谢物模块代表一组协同变化的代谢物那么和它显著相关的基因模块很可能参与同一生物学过程。溯源结果这时候又能提供另一层线索帮你判断代谢物模块里里外外的因果关系比如微生物来源的代谢物模块自然更容易和肠道菌群的丰度数据形成呼应。我在实际使用中最后养成了一个习惯每次拿到TidyMass2的溯源结果我都会随机抽20个代谢物人工去HMDB和文献里核对一遍来源标签看看是不是有过度预测的嫌疑。这不是对工具不信任而是代谢组学的注释错误率本来就不可忽略尤其碰到同分异构体和加合物干扰时任何工具给的结果都值得抽查一遍。你把这个习惯保持住再让TidyMass2帮你把脏活累活自动化整个分析流程就能又快又稳地把故事讲完整。