MelonnPan:Predictive metabolomic profiling of microbial communities using amplicon or metagenomic sequences
- 前言
- 引言
- 结果
- 1.The MelonnPan algorithm
- 2.MelonnPan accurately predicts metabolites from metagenomes
- 3.Estimating MelonnPan accuracy in new microbial communities
- 4.MelonnPan predictions outperform existing methods
- 5.MelonnPan uncovers meaningful biological relationships
- 6.Predicted metabolites reveal global structure in the IBD metabolome
- 7.Inference across the human body and environmental microbiomes
- 方法
- 1.Training and validation cohort descriptions
- 2.Taxonomic and functional profiling
- 3.Metabolite profiling
- 4.Filtering, transformation, and normalization
- 5.Elastic net regularization
- 6.Cross-validation and evaluation metric
- 7.Significance testing with shuffled data(不懂！）
- 8.RTSI score
- 9.Gene set enrichment analysis
- 10.Differential abundance analysis
- 11.Non-gut and non-human microbial profiles
前言
摘要
微生物群落代谢组，特别是在人体肠道中，开始提供一种新的途径来识别疾病中被破坏的功能和生态。然而，代谢数据是昂贵的且很难大规模获得，而数千人群体的扩增子或鸟枪式元基因组测序数据很容易获得。本文描述了一种计算方法，在给定新的微生物群落中来预测潜在的未观察到的代谢物，给出了一个模型，该模型训练来自感兴趣环境的成对代谢体和元基因组。聚焦于两个独立的人体肠道微生物组数据(PRISM和NLIBD)，我们展示了我们的框架成功地恢复了超过50%的相关代谢物的社区代谢趋势。利用珊瑚相关的、小鼠肠道和人类阴道微生物群的扩增图谱也保持了类似的准确性。我们还提供了一个预期性能得分(RTSI) 来指导模型在新样本中的应用。因此，我们的结果表明，这种“预测代谢体”方法可以帮助实验设计，并为目前仅有元基因组可用的数千个社区概况提供有用的（代谢物）解释。
结论
MelonnPan是一种推断与微生物群落相关的近似代谢物特征丰度的方法，它的验证和应用表明，微生物组分类和功能图谱中包含的信息与代谢组含量有足够的相关性，可以作为微生物群落生化环境的可操作预测。这不仅对人类的肠道有意义，在给定来自感兴趣环境的足够的训练测量（样本）情况下，也能推广到包括环境微生物在内的广泛栖息地，并且该模型提供了对新样本的预期性能(RTSI得分)的估计，以指导实验者。尽管不能取代代谢组分析，但它可以用代谢组成本的一小部分样本来近似预测和比较许多样本中可能的代谢情况，从而为更具成本效益的分层研究设计开辟道路，并在目前仅有元基因组数据的数千个现有样本中提供代谢解释和假设生成。
为了指导用户在将MelonnPan假设与下游实验验证相结合时，MelonnPan专门为每个新的微生物组提供了置信度分数(RTSI)，置信度分数较低表示与训练元基因组的高度不同。由于元基因组之间的训练差异(由RTSI捕获)会影响MelonnPan的准确性，因此RTSI值可以作为指导方针，指示在新的环境中可能需要多少额外的代谢数据来补充预先训练的MelonnPan模型 。这一信息尤其重要，因为MelonnPan以数据驱动的方式捕获元基因组-代谢体关联，即使在没有任何微生物生化注释的情况下也可以操作，这产生了比目前依赖于非常有限数量的良好表征的酶-代谢物关系的预测精度高得多的预测准确性。
与许多最近的研究一样，这项调查支持以尽可能高的分类分辨率表征微生物特征的重要性，因为与次生代谢物产生相关的主要微生物表型差异通常是物种或菌株特有的。具体地说，我们的分析证实了几个以前记录的与IBD相关的物种是肠道微生物-代谢物动态的重要驱动因素(图3A，补充数据7)。为了进一步询问物种丰度数据在人类肠道应用中是否会导致类似的代谢预测，我们使用物种丰度数据作为MelonnPan模型的输入预测因子 进行了额外的分析。虽然物种水平的预测因子导致了训练队列中类似的表现，但这些基于分类学的预测并不适用于独立的队列(补充图18)。验证队列中的这种显著较低的可预测性可能反映了由基因家族数据捕捉到的菌株水平效应，因为不同人群中的菌株差异会显著影响代谢物预测的概括性 。这突显了将基因水平图谱作为预测因子的重要性，因为特定菌株的特定代谢以及其他表型相关特征(例如抗生素耐药性)可能不能仅从物种丰度数据中捕捉到。
在使用MelonnPan预测时，必须考虑其局限性。它不直接预测代谢物流量或峰值(与基于约束的方法相反)；相反，它通过合成和组合微生物序列特征来估计每种代谢物在社区范围内的相对丰度。这里显示的最初应用主要应用于非定向MS代谢物测量，但我们已经证明MelonnPan在学习定向MS或NMR代谢物测量时仍保持相对准确。虽然MelonnPan经过充分研究可以用来在人类肠道以外的广泛环境中进行预测，但每个学习到的模型都是特定于环境的。因此，在人类肠道上学习的模型可以推广到其他人类肠道表型(补充图19)，但没有一个单一的模型对于跨环境预测任务是准确的。MelonnPan因此被用作假说生成工具，因为预测的代谢物相对丰度和测量的代谢物相对丰度之间的一般一致性通常足以为后续的实验验证研究提供信息，这些研究绝对应该进行，以确认预测并获得感兴趣的代谢物的直接测量值。
MelonnPan虽然是独立于机制和分子来源的预测，某些特定代谢物的强可预测性也可能对提出这样的机制有价值。例如，代谢物的可预测性可能表明它要么是由一组微生物产生的，要么是在特定微生物存在的情况下在宿主细胞中被刺激的。这种专门的代谢物将被用作一组已定义微生物的标志物，这是一种基于模型的关联发现途径 ，比通量平衡分析等方法更快，资源密集度更低。因此，从独立于培养的种群水平数据进行鉴定的未来工作有可能集中于菌株特有的基因组，甚至生物活性类群之间的单核苷酸多态性水平的差异。 进一步完善MelonnPan预测准确性的未来研究的其他方向包括(1)整合其他类型的微生物测量，如宏转录数据，(2)利用纵向剖面进行动态预测，以及(3)采用更复杂的机器学习策略，如多变量或贝叶斯框架，其可以明确地纳入量化特征，如群落范围的酶特异性反应信息和零膨胀（zero-inflation）等。
独立于培养的元基因组测序已经分析了数以万计的样本，这些样本包含数百万个分类群和微生物基因——其中数百万个没有特征化。例如，UniProtKB中只有大约1.0%的蛋白质被实验表征49。因此，MelonnPan的人体肠道模型的一个重要发现是，大多数(>60%)未注释的基因家族与代谢物相对丰度之间存在关联(补充数据7)。基因家族和代谢物之间的这种联系为基因本身的下游表征提供了有希望的靶点，特别是当应用于其他特征较差的环境时，因为它们可能在这些化合物的生成或代谢中起作用。因此，该计算方法得生物学意义(1)生成用于未来验证的生化和功能基因组假说，(2)有助于整个系统对微生物的理解50，(3)作为现有代谢重建模型的补充工具，(4)有助于为代谢组学在微生物群落中的翻译应用（translational applications）奠定实验设计基础。随着能与已知标准匹配的参考数据库继续饱和51,52以及训练数据集继续扩大53，MelonnPan的预测准确率将默认随着时间的推移而提高。综上所述，这一分析框架是迈向种群水平元组学数据集成的必要的第一步，最终使我们能够更好地理解微生物组的动态，超越分子目录，转向微生物组研究的健康应用。
引言
高通量代谢组学技术的进步，使微生物群落中的大量小分子代谢物得以全面覆盖。分析不同调控的生化途径之间的代谢差异，可以发现与疾病相关的潜在生物标志物，并为深入了解潜在的发病机制提供见解。这一点已被越来越多的研究突显出来，这些研究依赖于多组学图谱来同时表征人类微生物群或其他环境的群落生态、代谢特征和功能属性。例如，研究最多的影响人类健康的生物活性微生物的代谢物是短链脂肪酸(SCFAs)，包括丙酸、丁酸和醋酸盐（propionate, butyrate, and acetate），它们与炎症性肠病(IBD)和结肠直肠癌等多种疾病的发病机制有关13 - 15。其他的例子包括胆汁酸（bile acids）16、鞘脂（sphingolipids）17和色氨酸衍生物（tryptophan derivatives）18都有证据表明肠道中的微生物相互作用和生物活性。
因此，推断微生物群落产生分子的能力，并使用大规模数据集将新的特定基因与代谢物联系起来，是实现了解肠道微生物代谢如何以及为什么影响人类健康的关键的第一步。肠道微生物和代谢谱之间的关联强度表明，根据微生物群落的分类学或功能谱，从元基因组中近似预测微生物群落的代谢活动或特征是可能的。仅仅基于酶的作用很容易识别这种关联，这在很大程度上受到目前不饱和的基因-代谢物反应谱系的限制，以及由典型的测序和代谢组技术提供的相对(而不是绝对)丰度指标的限制。然而，尽管有这些限制，预测与肠道微生物谱相关联的代谢物特征的方法可以用作假说生成器，其可以促进新关联的群体规模的发现(例如，在大型元基因组数据集合中)，并导致新的可检验假说集合，作为实验验证研究的补充(例如，从扩增数据预测功能简档的情况)。
最近，已有研究采取了初步步骤，在事先知道微生物组和代谢组之间的联系机制的情况下(例如，从KEGG21等数据库得出的化学计量酶反应矩阵)进行这样的预测。一组方法，统称为预测的反应代谢周转率(PRMT)，计算基于群落的代谢物潜力(CMP)得分，它代表给定样本中群落产生或消耗每种代谢物的相对能力22-24。其他方法在基于约束或基于网络的建模框架25-27中重建社区新陈代谢的预测代谢模型。这两种方法的一个共同缺点是它们不能确定是由于参考数据库中缺少注释或准确的反应信息而导致的预测失败，还是由于替代生物学机制而导致的失败，使得它们难以以数据驱动的方式应用或验证。此外，这些方法依赖于对物种( 甚至菌株)-特异性代谢物间的准确描述和注释，而且它们不能很好地扩展到具有部分参考的类群或代谢物的复杂群落。因此，所有这些研究都将微生物功能潜力与代谢活性联系起来，但它们的范围有限，缺乏严格的外部(独立)验证，特别是在人体肠道等环境中，在这些环境中，代谢组训练测量是可行的，准确和新颖的生物活性发现可能会产生特别的影响。
我们定义了MelonnPan(Model-based Genomically Informed High-dimensional Predictor of Microbial Community Metabolic Profiles，基于模型的基因组信息高维微生物群落代谢谱预测器)，这是一个从微生物群落分布预测群落代谢组的计算框架。MelonnPan能推断出复合代谢组通过(1)对最优预测微生物特征的数据进行识别，(2)对预测良好的代谢产物的预测准确性进行稳健的量化。这使得研究人员可以从目前只有宏基因组可用的群落可再生产地推断代谢物。我们将MelonnPan应用于两个独立的肠道宏基因组数据集，包括>200例克罗恩病(CD)患者、溃疡性结肠炎(UC)患者和健康对照(HC)参与者。结果表明>50%代谢物在预测和观察到中的代谢趋势高度一致，这些代谢物的身份与实验室标准相一致，包括与胆汁酸、脂肪酸、类固醇、戊烯醇脂类和鞘脂类相关的代谢转移预测。当使用扩增子测序谱的分类特征时，珊瑚相关、小鼠肠道和人类阴道微生物群落也保持了类似的准确性。MelonnPan的实现、相关文档和示例数据集可以在MelonnPan软件包http://huttenhower.sph.harvard.edu/melonnpan中免费获得。
结果
1.The MelonnPan algorithm
我们开发了MelonnPan作为一种计算方法，通过结合分类学或功能谱形式的生物学知识，从扩增子或元基因组测序数据中预测代谢物的特征。与现有的基于化学计量学的方法依赖于有限数量的特征良好的分类群、酶和代谢物不同，它不需要功能注释，因为该工具的设计目的是使用机器学习捕获两者间关系，即使是对于未表征的微生物特征。本文我们专门讨论了它在人体肠道微生物群中的应用，但该方法可以推广到任何适当的微生物环境。简而言之，MelonnPan使用弹性网络规则化28来确定对给定的代谢物而言，哪些特征(分类或功能)是可预测的。给出一个新的分类图谱(来自扩增片段或元基因组)或元基因组功能图谱(即基因家族丰度)，然后它结合序列特征的子集来估计相关的复合代谢组。得到的预测代谢物每个都是预测特征(分类群或基因家族)的相对丰度的加权和，其中训练的弹性网络模型的回归系数被用作预测算法中的权重(图1)。
在拟合阶段，MelonnPan使用既有测序数据又有实验测量的代谢物丰度的样本进行训练(图1a)。这两个指标都是有效的相对丰度——分别是归一化读数和光谱计数。训练(以及后来的推断)元基因组可以用任何量化分类群或功能相关微生物基因家族的相对丰度的系统来描述；这里我们使用之前描述的扩增数据 和由HUMAnN2与UniRef90功能描述的元基因组作为参考目录，即来自UniProt的至少90%氨基酸同一性的序列簇集合 。
该拟合过程使用按代谢物弹性网络正则化来优化少量序列特征的系数 。特定环境(人体肠道或其他)的最终模型是基于严格的内部验证(交叉验证)选择的，该内部验证对应于每种代谢物的最大交叉验证似然。基于预测和训练测量之间的rank相关性(Spearman相关系数<0.3)，对任何通用模型都不能很好地预测的代谢物进行标记。最后，训练好的模型可以应用于来自类似环境的新微生物群落，使用简单的线性回归，将学习到的系数值乘以序列特征丰度(图1b)。在新数据集中的预测应用之前的模型评估期间，不使用来自测试集的信息来训练模型。在这里，对于我们的主要特定于肠道的模型，我们还将最终的内部验证模型应用于独立的外部验证队列。性能被总结为每个代谢物在所有样品与相应的被测代谢物(“方法”)之间的Spearman等级相关系数。
2.MelonnPan accurately predicts metabolites from metagenomes
我们使用两个独立的元基因组和代谢组学数据集——分别包括155名和65名IBD患者和对照，CD(n=68和20)，UC(n=53和23)和HC参与者(n=34和22，补充表1)，验证了最初的MelonnPan人体肠道模型。在每个横断面队列中[马萨诸塞州总医院(PRISM)和荷兰IBD队列(NLIBD)]，粪便样本通过猎枪元基因组测序和四种液相色谱串联质谱(LC-MS)方法(包括正负离子模式下的极性化合物、脂类以及游离脂肪酸和胆汁酸)进行分析。LC-MS图谱产生了约8000个聚类特征，以色谱保留时间和准确质量为特征。所有带有HUMAnN2的样本都生成了元基因组功能图谱，产生了大约100万个UniRef90基因家族。在建模之前，这两种数据类型都经过了质量控制和过滤：当特征的值在可用样本上没有变化时，它们将被移除。特别是，相对丰度和通用率极低的基因家族和代谢物(在≥10%的样本中<0.01%)都被去除了，只剩下2,818种代谢物和814个基因家族供最终建模(“方法”)。所有训练数据，包括内部交叉验证，仅使用PRISM的数据进行；绩效评估和外部验证使用荷兰的held-out样本(NLIBD)进行。
在这些数据中，我们模拟了一个由466种代谢物组成的小组，这些代谢物的身份已经根据实验室标准31进行了实验确认。经过初步过滤、模型拟合和内部交叉验证，MelonnPan在独立验证期间>50%能很好被预测(n=107，53.8%)(预测代谢物丰度与实测代谢物丰度之间的Spearman相关性≥0.3)(补充图1)。预测良好的代谢物(图2a)包括鞘磷脂[例如神经酰胺和植鞘糖苷，脂肪酸(例如二十二碳五烯酸和己酸)，B类维生素(例如泛酸，图2d)，以及胆固醇和胆汁酸的衍生物[例如胆甾酮(图2c)和胆酸]。这些化合物越来越被认为是调节系统宿主-微生物共代谢的重要信号分子。在补充数据1-4中提供了完全拟合的代谢物和受试者的模型参数、预测和性能总结，实施和模型可以在http://huttenhower.sph.harvard.上在线获得。
在交叉验证期间和独立的held-out元基因组中，都能很好地预测到大量的代谢物(补充数据5)；因此，MelonnPan可以成功地在整个元基因组样本中恢复适度大量的实验验证代谢物的代谢组谱(图2b，补充数据6，补充图2)。虽然这些验证主要测试使用实验标准验证的已识别代谢物的子集，但该模型可以训练并适用于所有独特的代谢物，即使是那些没有确定分配标签的代谢物(补充图3)。在具有唯一簇ID的2818种代谢物中，大约60%(n=1679，59.6%)在训练期间具有Spearman≥0.3的准确率(补充图1)。其中，大量(n=933，55.6%)与标记化合物一样，在交叉验证期间和独立验证数据中都得到了很好的预测。考虑到肠道微生物群落的复杂性，以及可能潜在影响代谢物丰度的众多外部和内部影响，这种高度的可预测性有点令人惊讶。这突出了MelonnPan的预测能力的健壮性，即使在复杂、不稳定和动态的社区(如人类的肠道)的环境中也是如此。
在每个预测良好的代谢物中，MelonnPan模型(补充数据3，补充图4)平均选择了<2%的基因家族(中位模型大小=12，中位正权重=9，中位负权重=11)。这在统计学上和生物学上都是有意义的，因为MelonnPan由强加正则化来估计稀疏模型并识别一小部分相关序列特征(其中弹性网络混合和稀疏参数都是基于内部十重交叉验证的“方法”来选择的)。从MelonnPan学习的权重表示给定样本中的特征与每个代谢物相关(正或负)的相对能力(假设样本之间有一些基线代谢物分布)。在生物学上，这允许我们通过逐步将相关性较小的特征的贡献设置为零并仅保留少量特征来建立可解释的模型——在这种情况下，具有潜在酶作用的基因或与目标代谢物有其他关联的基因。
作为额外的验证，为了进一步排除数据由人工制作的可能性，我们测试了MelonnPan在试图将随机的、无效的微生物图谱与代谢物图谱联系起来时的行为。我们独立地对样本中的代谢物和基因家族训练数据进行了排列，然后在排列之后对其进行重整化，以保留每个单独数据集的核心特征 (“方法”)。我们重复了1000次这个过程，每次都从训练好的MelonnPan模型中收集结果系数，并对迭代中预测良好的代谢物的数量进行平均，以得出最终的可预测性(在没有过度拟合的情况下，它将保持在接近零的水平)。我们发现，在对训练和独立验证数据集进行评估期间，与真实的未随机打乱的数据相比，随机化的空轮廓产生被认为预测良好的代谢物集合始终非常低(59.6%的真实化合物被良好预测，而经过置换后训练数据的真实化合物预测良好的比例为3.2%；55.6%的真实化合物在验证数据中得到很好的预测，而在置换之后则为4.4%，McNemar的精确检验P<0.0001，补充图5-6)。这既表明MelonnPan对过度拟合的健壮性 ，也证明了其任意阈值(即Spearman r>0.3)的合理性，对于总结足够多的“预测良好”的代谢物来说，这一点大体上是一致的。
3.Estimating MelonnPan accuracy in new microbial communities
MelonnPan的实用性还在于它从新的微生物群落样本中预测代谢物的准确性，以及相应地概括代谢组研究结果的能力。为了表征这种影响，我们开发了代表性训练样本指数(RTSI，Representative Training Sample Index)(按照PICRUSt的最近序列分类指数的精神)来量化新样本对于训练数据集的代表性(“方法”)。简而言之，MelonnPan首先标记训练元基因组中不存在的任何特征(分类单元或基因家族)，对于剩余的共同特征(训练样本和测试样本之间)，它基于主成分分析(PCA)计算平均相似性分数(每个微生物群落样本)。具体地说，RTSI分数是通过顺序寻找与最优主成分(PCs)的最高相关系数来计算的，这些PC解释了训练微生物群中的大部分差异(“方法”)。我们比较了所有预测良好的代谢物(补充图7)中所有NLIBD验证样本的RTSI分数和MelonnPan精确度，发现代谢物的代表性(更高的RTSI)和样品的可预测性(通过测量和预测的代谢物丰度之间的Spearman相关性来衡量)之间有很强的统计上的显著相关性(Spearman相关性=0.40，P=0.003)。这个相关值本身是保守的低值，这是由几个离群值样本造成的，如果没有这些离群值样本，MelonnPan性能的RTSI在新样本上的预测性甚至更大。
这些见解具有多个层面的潜在影响。首先，该方法为MelonnPan用户提供了一种适当谨慎的方式，能应用该模型来预测非人体肠道模型（默认的）或用户可能已经训练的其他模型的群落代谢物特征。MelonnPan的RTSI值使用户可以在运行分析之前确定他们的样本是否足够相似来进行可靠的MelonnPan预测。其次，NLIBD队列中的这些评估结果证实，在没有测量代谢组学数据的情况下，MelonnPan预测可以被认为是大量代表性元基因组群落代谢谱的替代品。这既可以作为一种假设生成工具，也可以作为一种有效、成本效益高的首次通过分析驱动实验设计，我们建议将其与后续实验配对，以证明推断的代谢物图谱 (就像使用PICRUSt进行扩增和元基因组数据分析的情况一样)。
4.MelonnPan predictions outperform existing methods
将MelonnPan与最近开发的代谢物预测方法MIMOSA所确定的代谢物24进行比较。选择与MIMOSA进行对比，因为它是目前唯一能够从种群水平的元基因组数据中预测整个群落代谢相对丰度的方法，并提供了软件实现。MIMOSA建立在以前提出的代谢网络模型(PRMT)22的基础上，根据分类组成和元基因组含量来估计微生物群落的代谢潜力。简而言之，MIMOSA首先生成一个化学计量矩阵，描述基因和代谢物之间的定量关系，以提供感兴趣群落的CMP得分的估计。接下来，它将所有样本之间CMP得分的差异与相应测量的代谢物的差异进行比较。为了识别具有统计意义的预测准确的代谢物，MIMOSA依靠基于Mantel测试32的错误发现率(FDR)校正的P值来确定两个距离矩阵之间的相关性。与MelonnPan相似，MIMOSA也依赖于用户提供的代谢物和微生物序列特征的配对表。然而，与MelonnPan不同的是，它没有显式地使用数据挖掘和模型构建来构建和验证预测模型。因此，MIMOSA不能对算法以前没有看到的新的元基因组样本进行预测。
为了评估MIMOSA在NLIBD队列中的预测性能，我们首先将化合物id映射到人类代谢组数据库(HMDB)，将对应的代谢物名称映射到KEGG标识符，包括对KEGG化合物标识符的交叉引用，得到303个KEGG化合物(149个唯一标识符)。接下来，我们通过将UniProt-KO注释赋给UniRef90中相应的蛋白家族，将UniRef90基因家族的丰度转化为与相应KO (KEGG Orthology)丰度近似的丰度。为了应用MIMOSA，我们首先使用MUSiCC34(MIMOSA软件中的默认和推荐选项)标准化KO丰度。只有少量的代谢物化合物被MIMOSA预测得很好(补充图8A;n=20 (23%)， Mantel检验Q< 0.05)。相比之下，MelonnPan能够准确预测这些代谢物中的绝大多数(n=130(84%)，Spearman r>0.3)。此外，与MelonnPan不同的是，少数代谢物被MIMOSA反向预测 ，MelonnPan在两种方法预测的共同代谢物之间通常产生更高的置信度(测量和预测的丰度之间更大的Spearman相关性)(补充图8B)，并且即使限于特征良好的代谢物的一小部分时，也能获得更多预测良好的代谢物(补充图8C)。这表明，带注释的微生物酶活性的比例相对较小，与使用机器学习新识别的微生物酶活性之间可能存在重大差距。
5.MelonnPan uncovers meaningful biological relationships
为了深入了解在这一背景下最具预测性的基因家族的分类和功能构成，我们接下来量化了每个基因家族对人类肠道的MelonnPan预测做出了多少贡献。我们执行了基因集富集分析(GSEA)35来确定代谢模型最频繁选择的基因类别，即在代谢物预测期间显著过高或过低的基因(根据基因家族特征在任何预测良好的代谢物的模型中使用的频率对基因家族特征进行排序)。我们使用HUMANN2对这些基因家族进行了分类，根据基因的分类分配(如果不是属水平上的独有基因，则是“未分类”)，选择大多数同源性基因的最低共同祖先来对这些基因家族进行分类。然后，我们用基于排列的Kolmogorov-Smirnov(KS)检验(“方法”)比较了每个基因集合(归类为属)中基因家族累积分布的总体差异与MelonnPan排名表中的基因家族分布的总体差异。
在MelonnPan基因列表中有8个属明显过多，其中在假黄酮菌、梭状芽孢杆菌、Coprorics、厌氧菌、Blautia、Collinsella、反刍球菌和厌氧菌中观察到的影响最强(GSEA结果有统计学意义，Q<0.25，图3A，补充图9)。除Collinsella属属于放线菌门外，大多数属都来自Firmicuts门。其中一些属，包括梭状芽孢杆菌和瘤胃球菌，编码了几种属于梭状芽孢杆菌簇IV或XIVa36的物种，它们优先定植于粘液层，从而提高了丁酸对结肠上皮细胞的生物利用度。这些物种在人类结肠中相对丰度的减少与几种疾病有关，包括IBD。此外，梭状芽孢杆菌第四簇的物种也被认为是人类结肠中单链脂肪酸的主要生产者，这种单链脂肪酸越来越被认为是肠道微生物区系和宿主之间的关键信号分子。对这些属的分解表明，它们进一步典型地是由几个有代表性的种或菌株贡献的(补充数据7)。因此，MelonnPan能够识别与物种甚至菌株特定的代谢相关的功能相关的基因家族，便于在更精细的分类学分辨率下进行生物相关的机制研究。
为了以更高的分辨率进一步破译这些基因家族，我们使用功能注释重复了富集分析，以识别在代谢物预测中显著过多或过少的生物过程(“方法”)。我们关注的是Pfam数据库，该数据库将这些代谢预测基因家族归类为蛋白质家族。令人惊讶的是，虽然在该测试中没有单个Pfam家族被富集，但在预测基因家族中存在显著的未表征蛋白结构域的过度表达(Fisher精确检验P=3.46e−52，图3B)，这在我们考虑的所有功能注释目录中都是一致的(补充图10)。这相当大数量的未注释基因可能包括尚未确定的代谢酶类，它们在群落新陈代谢中具有潜在的作用。这与许多未确定特征的微生物基因在解释人类肠道内绝大多数微生物多样性和功能上的拟议作用是一致的。这表明，将计算功能预测与微生物生理生化验证相结合的协同方法将是必要的，以将特定的微生物化学与来自代谢组源的新的个体生物、基因和酶联系起来。
6.Predicted metabolites reveal global structure in the IBD metabolome
最近的几项研究表明，IBD患者和健康人，以及IBD亚型(UC和CD)可以通过代谢图谱区分开来，这表明IBD代谢组将成为测试MelonnPan代谢物预测准确性的有意义的基准。为了证明MelonnPan可以在不直接测量代谢物的情况下捕捉代谢谱中的生物变异，我们比较了测量的代谢物中代谢变异的前两个成分，并将推断代谢物的预测变异叠加在同一个二维空间中 。具体地说，我们根据斯皮尔曼(Spearman)对前50个独特代谢物簇的预测和测量代谢物成分的不同，对NLIBD队列中65名受试者的主坐标进行了协调，这些代谢物簇的身份已根据实验室标准进行了确认。排序图显示了IBD微生物群中类似的全局结构，这反映在跨化合物的测量和预测剖面的接近上(图4)。
预测良好的与IBD相关的代谢物跨越了广泛的代谢类别，包括氨基酸、胆汁酸、脂肪酸和鞘脂等。特别是，排序显示在IBD代谢体结构中大约有三个簇得到了强有力的支持。共变代谢物的主要排序基团包括来自同一母体化合物或通过共同途径相互转化的化合物，包括：(1)几种胆汁酸和在IBD中枯竭的超长链脂肪酸基团(补充图11，右簇)，(2)几种富含IBD的胆固醇和四吡咯衍生物(左簇)，以及(3)一小部分非差异性丰富的代谢物，如氨基酸、肽、嘌呤及其衍生物(中心)。对预测和测量的代谢组学成分使用相同的差异丰度分析(“方法”)，在代谢物中产生高度相似的定量结果(基于测量和预测的曲线的效应大小估计之间的Spearman相关性分别为CD与HC比较= 0.70,UC与HC比较= 0.45;P < 2.2 e−06;补充图12)，表明即使没有全面的代谢组学分析，MelonnPan预测也可用于从宏基因组中推断疾病相关的代谢组学成分差异。
下面我们开始确定在MelonnPan预测IBD中，明显过量或不足的化合物的大类。我们关注的是在我们的数据集中至少有一个成员的类，并确定了在FDR校正后(Fisher’s exact test Q< 0.25，“Methods”)在MelonnPan预测中具有统计学显著性的强化。MelonnPan预测的化合物中有两种代谢类明显过量，其中胆汁酸和四吡咯的影响最强(补充图13)。在预测良好的代谢产物中，胆汁酸相关产物的富集凸显了群落生态学在胆汁酸微生物代谢中的重要作用。与先前的研究一致，IBD患者的初级胆汁酸显著升高，次级胆汁酸显著降低。胆汁酸的生物合成是由微生物酶活性直接介导的，在IBD中，微生物酶活性不能与初级胆汁酸解偶联，导致次级胆汁酸的减少及其对肠上皮细胞的抗炎作用。在其他丰富的代谢类别中，与对照组相比，IBD受试者中的四吡咯趋于持续消耗。总之，这些发现证实了MelonnPan能够提供广泛化合物的代谢相关预测，确定肠道微生物代谢轴的重要贡献者，进而促进微生物组的大规模集成多组分析。
7.Inference across the human body and environmental microbiomes
将MelonnPan应用于三个配对的16S和代谢组学数据集，这是MelonnPan从各种微生物环境和分析类型中产生生物学解释能力的最后一个例证。这些数据之前是在以下背景下产生的:(i)与生态临界造礁珊瑚相关的微生物群落的代谢组学和分类学特征42，(ii)阴道微生物群落的细菌群落和代谢组学特征43，以及(iii)配对的小鼠肠道分类和代谢组学样本44。通过16S rRNA基因扩增子(16S)测序，对每个数据集的样本进行分类组成分析。对于代谢产物，数据集1采用质子核磁共振(1H-NMR)技术，数据集2和3分别采用靶向LC-MS和非靶向LC-MS和GC -MS代谢组学技术(“Methods”)。
在这些数据集中，我们使用MelonnPan学习一个模型，从可用的微生物特征(分类剖面来自16S rRNA基因测序)中预测相关代谢物的相对丰度特征。由于这些数据集的样本量较小，我们使用留一交叉验证(LOOCV，leave-one-out cross-validation )对MelonnPan进行训练，然后对代谢物和操作分类单元(OTU)特征(“Methods”)进行独立过滤(即当特征值与可用样本之间没有变化时，将其删除)。在每个数据集中测定的代谢物中，>50%没有通过单独的预过滤(见“方法”，补充表2)，因此从下游分析中被丢弃。
我们发现，在每个数据集中，>60%的分析代谢物都得到了很好的预测(补充图14，补充表2)，这些通常与少量的OTU有关(数据集1-3的模型中值分别为29、14和32)，这表明，对于相当大一部分化合物，几个分类特征中包含的信息足以解释代谢物丰度的大部分变化(与人类肠道基因家族的应用一致)。再一次，MelonnPan学习的基于数据的模型比MIMOSA使用的机理模型(补充表2)要精确得多，后者仅限于这些更具挑战性的微生物环境中极少数已经表征的化合物(补充图15-17)，进一步强调了MelonnPan在各种生态环境中用于微生物代谢物假说生成的实用性。需要提醒的是，与人类肠道样本不同，我们无法在这些环境中获得独立的验证数据集。因此，我们认为这些应用仅仅是对MelonnPan预测在研究较少的环境中的潜在概括性(外部有效性)的初步评估。
方法
1.Training and validation cohort descriptions
Franzosa等人详细描述了培训和验证队列。简而言之，训练队列的受试者来自PRISM，这是一个以转介中心为基础的预期队列 。年龄大于18岁、根据标准内窥镜、放射学和组织学标准诊断为CD的患者有资格参加。共有155名成人患者入选，包括CD和UC患者和非IBD对照组(68名CD，53名UC，34名HC受试者)。PRISM研究方案由合作伙伴人类研究委员会(#2004-P001067)审查和批准，所有实验均遵守该审查委员会的规定。
验证队列的受试者来自荷兰的两个独立队列。队列1由22名非IBD(HC)受试者组成，他们参加了荷兰北部的普通人群研究LifeLines-DEEP(LLDeep)；队列2由43名IBD患者(UC=23，CD=20)组成，这些患者来自荷兰格罗宁根大学医学中心(UMCG)胃肠病和肝病科。共采集65份粪便样本。在这两个队列中，使用相同的方案来收集粪便样本。
2.Taxonomic and functional profiling
元基因组数据的生成和处理是在Broad研究所进行的。从粪便样品中提取DNA后，根据制造商的推荐方案来使用Nextera XT DNA文库制备试剂盒(Illumina)来制备元基因组文库，并在Illumina HiSeq 2500平台上进行测序，目标是每个101bp样本约读取2.5 GB的序列，成对端阅读（paired-end reads）。使用Trimomatic过滤长度小于60nt的低质量读取以及使用bowtie2过滤人类污染读取；这些步骤是使用KneadData管道(https：//bitbucket.org/biobakery/runaddata)执行的。
使用MetaPhlAn2(https://bitbucket.org/biobakery/Metlan2)推断所有样本的种级分类丰度，并使用默认参数运行。使用HUMAnN2进行功能分析。简而言之，HUMAnN2将元基因组读数映射到分类学描述步骤中上游识别的物种的泛基因组。这些泛基因组中的蛋白质编码序列已经被预先注释到它们各自的UniRef90家族(一个全面的、非冗余的蛋白质序列数据库)中。与已知的泛基因组不一致的读数，通过DIAMOND翻译搜索分别映射到整个UniRef90。所有的命中值(hits)都是基于比对质量和序列长度进行加权的，每个物种和未分类的命中值结合在一起，以RPK(每千个碱基的读数)为单位生成每个蛋白质家族的群落总数(除了物种分层的总数)。RPK单位进一步归一化为RPKM单位(每千碱基读取每百万样本读取)，以考虑样本间序列深度的变化。
3.Metabolite profiling
来自两个队列参与者的粪便样本(PRISM横断面样本155份，NLIBD横断面样本65份;重量范围50.5-167.8 mg)按照Franzosa等人的描述进行处理。四种测量互补代谢物类别的LC-MS方法被用于测量每个样本的极性代谢物和脂质。使用Genedata Expressionist v9.0对原始LC-MS数据进行处理，用于化学噪声去除、RT比对、峰检测和同位素聚类。四种LC-MS方法的结合产生了8869个聚类特征，其特征是色谱保留时间和准确质量小于5ppm的精确度。根据HMDB精确的m/z匹配，3829个代谢组特征与假定的识别器相关联。使用从内部化合物文库产生的参考数据，更精确地鉴定了一个466种代谢物的子集。更多LC-MS代谢组学实验的细节见Franzosa等。
4.Filtering, transformation, and normalization
将原始测量结果归一化为相对丰度后[0,1]，我们的分析仅限于那些在至少10%的样本中平均相对丰度大于0.01%的普遍和丰富的特征(物种、基因家族和代谢物)。由于元组数据的特定属性会显著影响模型的构建，如组合性、稀疏性、偏斜性(skewness)、均值-方差依赖性和极值，我们将输入特征(物种或基因家族丰度)分位数变换为标准正态分布的分位数，以提高弹性网络模型60-61的检测能力（这一方法已被广泛应用于遗传关联研究，并已被证明是建模非正常表型62 - 64的稳健方法），我们使用R(版本3.5.1)的GENABEL包中的rntransform函数进行基于分位数的逆正态变换65。为了识别预测特征的最佳子集，在应用线性模型时，代谢物相对丰度被反正弦平方根变换(Arcsin)66以近似同方差 。模型对转换后的数据进行拟合，并对得到的预测结果进行反向转换(Sqsin)，以保持预测的代谢物成分的覆盖范围 。
5.Elastic net regularization
我们设计并实现了用于代谢组学预测模型构建的弹性网正则化技术28。注意，每个代谢物弹性网模型适合于秩转换特征(物种或基因家族丰度)。关于弹性网方法及其变体的更多细节以前已经发表过67。利用R (version 3.5.1)中的glmnet软件包拟合弹性网模型，并基于交叉验证选择调优参数(即弹性网络混合参数α和稀疏性参数λ)。
6.Cross-validation and evaluation metric
十倍交叉验证(除非另有说明)用于确定弹性网模型的调优参数(tuning parameters)。用真实代谢物组成和预测代谢物组成之间的Spearman相关系数®来评价每个化合物的可预测性。根据Cohen68，我们将r>0.3的代谢物称为“well predicted”，将其余的标记为预测不良的代谢物。
7.Significance testing with shuffled data(不懂！）
为了量化我们的框架是否识别出了比偶然情况下预期的更多的可预测代谢物(即当基因和代谢物之间的所有共享信号都被破坏时)，我们反复打乱代谢物和基因家族表中的样本标签 ，利用随机化数据应用MelonnPan模型将基因与代谢物联系起来，并将这些随机化数据得到的预测良好的代谢物数量与原始数据得到的预测良好的代谢物数量进行比较。随机数据是根据R(版本3.5.1)包pecante中概述的方法生成的，使用了randomizeMatrix函数(该函数在样本内采用排列后重正化 ，以保持原始数据集的核心结构特征)。该过程重复1000次，以估计在训练和验证队列中预测性能的零分布(null distribution)。
8.RTSI score
为了计算新样本的RTSI，我们使用主成分分析(PCA)，通过基于分位数的逆正态变换，提取出反映训练宏基因组中种群结构的连续变异轴。具体地说，我们根据Tracy-Widom的统计数据69选择了排名靠前的PC。基于排名靠前的PC，我们根据NLIBD样本与提取的PC的最高相关性将其分为相似或不相似两类。我们将由此产生的相似性(相关性)称为RTSI得分。我们使用R(版本3.5.1)中的AssocTests包和0.05的显著性阈值来选择*PC的数量。
9.Gene set enrichment analysis
在构建合适的基因集合后，我们进行了两种类型的富集分析 ：过度表征分析和GSEA。
对于过度表征分析，我们创建了2乘2的表格，将属于该类别的候选基因的数量与不属于该类别的候选基因的数量进行比较，并使用单尾Fisher精确检验评估了过度表征的重要性。
对于GSEA，我们基于KS检验统计量计算了富集分数，该统计量反映了基因集在整个排序的基因列表的极端(顶部或底部)被过度代表的程度(根据基因家族特征在代谢物特征中的总体可预测性对基因家族特征进行排序)。为了评估重要性，我们使用R(版本3.5.1)包gsEasy的功能对每个基因集进行了100,000次排列。对于这两种分析，我们使用Benjamini–Hochberg错误发现率(FDR)方法70对多重假设检验进行了校正。
10.Differential abundance analysis
为了对NLIBD队列中测量和预测的代谢组学数据进行差异丰度分析，在调整了疾病状态(以HC为参考类别)、年龄和药物[免疫抑制剂(是/否)和抗炎(是/否)]后，我们分别对每个对数转换的代谢物相对丰度分布拟合了线性模型。
11.Non-gut and non-human microbial profiles
我们从可公开获得的非人类肠道环境数据库中获得了几个先前发表的数据集，每个数据集都是基于16S rRNA基因的分类学数据与代谢组谱配对的24、42。在每个数据集中，我们都使用了通过作者公开获得的经过处理的图谱，但是相关的序列数据也可以通过NCBI获得。
对于数据集1，我们使用了16S rRNA基因测序提供的分类图谱，并结合了基于1H-NMR的代谢组。
对于数据集2中的阴道样本，从阴道拭子中进行16S rRNA基因分析，并收集成对的宫颈阴道灌洗液进行代谢组学分析(针对180种化合物进行LC-MS)。
在数据集3中，再次使用16S rRNA基因扩增序列分析分类组成，并使用全球LC-MS和GC-MS测定代谢物。
在下游分析中，超过10%的样品中相对丰度<0.0001%的OTU和代谢物特征都被丢弃。此外，还采用了方差滤波 步骤来去除方差非常小的特征。然后，将MelonnPan的弹性网正则化与LOOCV相结合应用于质量控制图谱 。