宏基因组微生物网络图（宏基因组学技术与微生物群落多样性分析方法）

自然界中大量的微生物无法利用纯培养手段在实验室条件下进行培养，传统微生物学的技术手段限制了环境微生物的研究。高通量组学技术的迅速发展使得人类对各类生态系统中的复杂微生物群落有了前所未有的认知。

本文从扩增子测序与宏基因组测序出发，对宏基因组学在微生物群落检测中的基本分析流程进行了介绍，指出利用大数据分析技术与手段来克服宏基因组学数据解析，并将分析结果用更易理解的形式展现出来是未来研究的重点和难点。

组学技术的出现实现了从分子水平对环境微生物及其功能进行检测与分析，为了解完整的环境微生物全貌提供了有效的途径。

微生物组学通常是宏基因组、宏转录组、宏蛋白组、宏代谢组等各类系统生物学技术和方法的总称，注重研究生物系统组成及群落中物种之间的相互关系、系统结构和功能的关联、以及群落结构与生态系统的关联等整体上的科学问题。

其中，以高通量测序技术为基础的宏基因组学是目前最为关键和成熟的组学方法，也为其他组学的研究提供了研究基础。

宏基因组学测序分析流程

宏基因组指为某特定环境中所有微生物的基因组的总和，宏基因组学研究则通过直接分析环境中微生物的DNA来获知微生物群落的遗传、功能与生态特性。

目前的宏基因组研究紧密依赖高通量测序技术，包括扩增子测序与宏基因组测序。

扩增子测序主要针对核糖体RNA基因（rDNA）和功能基因，前者对细菌或古菌16S rDNA 及真菌18S rDNA与内部转录间隔区（ITS）序列等分子标记进行扩增，后者对于微生物某些特定功能基因进行扩增。

宏基因组测序则是对环境中所有DNA进行测序。基因组测序成本较高，且对于后续数据分析的计算资源要求也相对较高。

相比之下，扩增子技术凭借测序与分析成本较低的优点目前成为了环境微生物组学研究的主要手段。

• 扩增子测序分析流程

环境微生物群落研究的随机取样导致其扩增子分析可重复性低，但通过增加生物学重复、删除单一样品仅出现一次的序列等手段可以较好地加以弥补。

研究者定义了多项指数对生物多样性进行量化，同时也为不同生物多样性的比较提供了方法。

在微生物生态学兴起后，不少宏观生态学的研究方法和手段也被逐渐应用于微生物生态学的研究中，为其提供新的研究思路。

微生物扩增子测序的分析方法多样，分析流程也不尽相同。

以16S rDNA测序为例，扩增子测序分析的主要方法和流程

随着测序技术的更迭，Illumina测序平台逐渐占据着微生物扩增子测序的大部分市场，目前以双端250 bp测序策略居多。

通过该方案拼接完成的序列可以作为挑选代表序列的起始文件。

目前OTU（可操作分类单元）与ASV（扩增子序列变异）是2种主要的代表序列形式。

• 宏基因组测序分析流程

近年来，随着测序技术在通量和读长方面的持续提高，其成本也不断降低，针对微生物群落的全部基因组DNA的鸟枪测序也在不断增加。

宏基因组测序数据的数据量大，需要更加专门的算法与软件来处理与分析。宏基因组学的分析常有一套通用的流程。

宏基因组测序分析的常用分析流程

虽然宏基因组测序分析的流程类似，但由于其测序数据量大，目前缺乏标准的分析工具进行统一化处理，不同的分析工具和方法在性能和速度方面差异较大，尤其是不同类型的微生物组数据往往还需要个性化的调整。

随着三代测序技术的普及，国内外现有针对宏基因组数据分析各个步骤设计的软件正处于飞速发展中。

微生物群落多样性分析方法

国际生物多样性公约对生物多样性的定义为“来源于包含陆地、海洋与其他水生生态系统以及它们组成的复合生态系统中的生物的可变性，其中包括物种内、物种间以及生态系统的多样性”。

微生物生态学中的微生物多样性按照描述物种的尺度进行层级划分，通常主要有分类多样性、谱系多样性、遗传多样性和功能多样性。

其中分类多样性与功能多样性常通过分析分类单元、功能基因或通路在不同环境下的分布情况进行衡量，谱系多样性通过计算不同分类单元在系统发育水平上的接近程度衡量，而遗传多样性需要通过更精细水平的组学研究技术进行相应的描述。

• 微生物多样性的数量描述

按宏观生态学描述习惯，多样性常根据空间尺度分为3个种类：α-多样性主要描述局部群落或斑块中的多样性，β-多样性主要描述不同群落间（或整个景观的）物种差异，γ-多样性则关注更大区域性尺度的多样性。

对于微生物生态学的研究，多样性分析常常聚焦于α-多样性与β-多样性。

由于取样和测序的随机性，分析结果并不能完全反应群落的真实状态。对于这类数据的物种累积曲线，随着样本大小的增加，序列数量以恒定速度线性增加，同时观测到的物种数目以递减的对数速率累积。

稀释化方法使得不同样本大小的物种累积曲线可以进行比较。使用该方法绘制的累积曲线称为稀释曲线，其绘制方式为保持样本中OTU百分比组成不变，构建具有相同OTU组成但具不同样本大小的样本物种累积曲线。

一般认为，当某一样本的稀释曲线末端趋于平缓时，即认为该样本的采样和测序已经近似完全。

稀释化方法的缺点在于对稀有种等信息会造成失真，故一般认为样本内物种符合随机分布或均匀分布时，稀释曲线才能有效工作。

在获得扩增子数据并根据这些数据计算α-多样性后，一般需要构建稀释曲线，并对OTU表进行重新抽取，以降低样本大小对于多样性指标间比较的影响。

1-α-多样性的数量描述

α-多样性的数量描述对象主要是物种丰富度与物种数量分布。以下以OTU为例进行说明，OTU表中OTU观察值（Sobs）可作为物种丰富度的观测指标。

除上表外，希尔数也是用于描述群落α-多样性的重要指数，是一类多样性指数组成的指数家族，整合了相对丰度、物种丰富度并消除了一些缺陷。

希尔数符合复制原则，即两完全相异群落的希尔数之和等于两群落混合后的希尔数。

2-β-多样性的数量描述

β-多样性所关注的是多个微生物群落或样本间的相似性或不相似性。

在对β-多样性的数量描述中，互补性是一个重要的描述角度，指2个样本之间包含对方所不包含物种的数量。2个样本的互补性越强，可以认为它们的β-多样性越高。

对于互补性的计算、描述以及延伸，多仿用了集合中的相应规则，互补性也可以用维恩图进行可视化表示，同时通过样本间的共享物种与特有物种计算样本间的相似性或不相似性。

样本间不相似性可以用距离指数来衡量，对于OTU表来说，全体样本的成对距离形成的矩阵称为距离矩阵或不相似性矩阵。

常见的OTU表中的数据代表了各样本中各OTU下的序列数，即每个样本不仅有OTU种类信息，同时还含有每个OTU的丰度信息，这类数据常被称为定量数据。应用中的另外一类数据，不包含每个OTU的丰度信息，常称为存在-缺失数据，也被称为1-0数据。

常用的相似性-不相似性指数有很多，每种形式的指数对于定量数据与存在-缺失数据的计算方法也不同。

其中，Jaccard距离是典型的存在-缺失数据距离指数，是以相似系数形式提出的。相比于Jaccard指数，Sørensen指数加大了2个样本共有OTU的权重。

针对生态数据分析中的“双零问题”（某些OTU同时在2个样本中均表现为缺失），Jaccard指数与Sørensen指数计算中，双零数据并不参与到数据计算中，因此称为非对称指数。

Bray-Curtis不相似度则为定量数据距离指数，在计算中考虑了丰度信息。由于考虑了双零数据的指数，Bray-Curtis不相似度也属于非对称指数。

针对类别多样性，Faith提出了谱系多样性（PD）的概念，并定义其为待观察物种在分支树上的最短进化分支长度之和。PD考虑了物种之间在进化水平上的差异，包含了物种表型性状与生态位等信息。

谱系α-多样性的计算方式基于Faith对于PD的基本定义。计算谱系β-多样性时，UniFrac指数是常用的计算指标，该指数根据不同群落包含的共有与特有谱系结构计算群落间的不相似性。

• 群落结构的解析方法

多数微生物生态学研究主要关注于不同生境或不同环境梯度下微生物群落的变化，一组高通量扩增子数据是对于微生物在多组时空样本下的观测。

这类研究所产生的数据集，现多用多元统计方法进行分析。

1-探索性方法

探索性方法提供了样本变化的主要梯度以及样本的相似程度，但即使样本经过分析后体现出了某些规律，仍然需要验证。

主成分分析（PCA）是最常见、应用最广泛的多元统计方法之一。在数学上，PCA实际是降维的过程。PCA采用欧几里得距离来度量样本之间的差异性，但样本覆盖梯度太长时（即多样本中有很多相同OTU）会出现马蹄效应等问题。

对应分析（CA）通常用于衡量由样本OTU数据反映的样本群落间的差异。CA规避了马蹄效应，然而CA排序常会伴有弓形效应，去趋势对应分析（DCA）可用于尽量减小弓形效应。

主坐标分析（PCoA）在概念上由PCA衍生而来，同样遵循了降维的基本思路，将样本空间压缩并投射到低维空间。由于使用样本间两两成对不相似性矩阵，PCoA的排序轴与原始变量间并不存在直接关系，但其方差解释度还是能够通过校正后的不相似性矩阵特征值给出。

非度量多维尺度分析（NMDS）是一种特殊的排序方法，进行分析时往往会进行多次迭代排序，以求取得尽量小的胁迫值（原的样品间不相似性被改的程度的量化指标），一般认为胁迫值小于0.15是可以接受的。NMDS分析中排序距离和样本之间的原始不相似性无关，排序轴不具有解释样本不相似性方差的作用，故NMDS排序图的排序轴上无法给出合理的解释度。

2-解释性方法

在分析不同样本组间的微生物群落差异时，往往还关注造成这种差异的环境因子，即将微生物群落差异看作响应变量（因变量），而将环境因子看作解释变量（自变量）。

由此，在探索性方法的基础上，解释性方法增加了一组解释变量。解释变量在每个排序轴上的分量表示该变量对于样品沿该轴分布的贡献。

冗余分析（RDA）与典范对应分析（CCA）是2种典型的解释性排序方法。

RDA可看作是PCA排序的一种拓展，加入解释变量后使得排序轴（主成分）被约束为解释变量的线性组合。同PCA类似，RDA也不适合于处理样本覆盖梯度长的数据集。

CCA是RDA不适用的情况下更好的选择，利用解释变量约束响应变量后进行对应分析的典范形式。

RDA与CCA的可视化展示则是在PCA排序图的基础上，增加了代表解释变量的向量（数量变量）或点（类别变量）。

3-统计检验方法

常见的对样本间差异进行统计检验的方法有：ANOSIM（分组相似性分析）、PERMANOVA（置换方差分析）与MRPP（多响应置换过程）等。

但在进行环境因子与样本群落间差异的相关性统计性检验时，传统的相关系数检验往往不能很好地实施。

Mantel检验方法可以对两个矩阵之间的相关性显著性进行检验。对多变量矩阵或控制变量矩阵进行相关性统计检验时，还可以使用偏Mantel检验，即选取一个解释变量，其余解释变量作为协变量。

变量分解分析（VPA）利用了偏分析的思想，将响应变量数据集中的总方差划分为单个解释变量的独立解释贡献以及联合解释贡献，常用于在确定了对微生物群落间差异有显著影响的环境因子之后，进一步说明不同环境因子对于不同群落的差异的贡献度。

结论

多组学技术的联合应用已经逐渐成为认识环境微生物群落及其功能的重要手段，而通过组学技术的应用，研究者逐渐意识到生活在土壤、淡水、海水、空气，甚至人体等环境中的微生物，其系统发育的多样性和功能的复杂度远远超过以往的认识。

目前，如何利用正处于蓬勃发展时期的大数据分析技术与手段来克服宏基因组学数据解析这一难关，并将分析结果用更易理解与操作的形式展现出来，这是从事环境微生物学、生物信息和统计学研究人员共同的挑战。

本文作者：彭玺，冯凯，厉舒祯，邓晔

作者简介：彭玺，中国科学院生态环境研究中心、中国科学院环境生物技术重点实验室、中国科学院大学资源与环境学院，博士研究生，研究方向为生物信息学、微生物生态学；邓晔（通信作者），中国科学院生态环境研究中心、中国科学院环境生物技术重点实验室、中国科学院大学资源与环境学院，研究员，研究方向为环境微生物生态。

论文全文发表于《科技导报》2022年第3期

本文有删减，欢迎订阅查看

,