CN104039982A

CN104039982A - 一种分析微生物群落组成的方法和装置

Info

Publication number: CN104039982A
Application number: CN201280064063.2A
Authority: CN
Inventors: 李胜辉; 冯强; 覃俊杰; 朱剑锋; 官远林; 王俊; 汪建; 杨焕明
Original assignee: BGI Shenzhen Co Ltd
Current assignee: BGI Shenzhen Co Ltd
Priority date: 2012-08-01
Filing date: 2012-08-01
Publication date: 2014-09-10
Anticipated expiration: 2032-08-01
Also published as: CN104039982B; WO2014019164A1; HK1196642A1; US20150242565A1

Abstract

本发明提供了一种用于分析环境样品中的微生物群落组成的方法和装置。该方法包括测序、初级组装、分栈、基于栈的高级组装和鉴定等步骤。

Description

一种分析微生物群落组成的方法和装置技术领域

本发明涉及宏基因组学（metagenomics)和生物信息学领域。特别地，本发明涉及用于分析环境样品中的微生物群落组成的方法和装置。背景技术

宏基因组学又称为环境基因组学，元基因组学，生态基因组学，或者群落基因组学，其是一门直接研究各种环境（例如自然环境）中的微生物群落（包含了可培养的和不可培养的细菌、真菌和病毒等的总和）的学科。研究各种环境中的微生物群落和物种多样性具有特别的益处。例如，对人脉道环境中的微生物群落和物种多样性的研究，对于菌群的临床药物开发以及人菌代谢途径的了解是非常有用的。然而，由于传统研究方法的限制，我们对环境（例如肠道环境）中的微生物组成知之甚少。特别地，由于环境中可能包含了不可培养的细菌、真菌或病毒，许多物种无法通过涉及培养的传统研究方法进行鉴定。

如今，全基因组鸟枪法（whole genome shotgun, WGS)已经在宏基因组学研究中逐渐兴起。这种方法一般通过高通量测序得到大量的测序片段（ reads )，然后通过组装得到较大的连接片段（ contigs )、拼接片段（scaffolds ) 、或甚至是全基因组。与此同时，新一代高通量测序技术得到了长足的发展，这为利用 WGS策略来认识群落结构、研究群落差异和功能提供了良好契机。例如，最近的宏基因组学研究已经在多种环境中，在发现新物种、解析微生物群落多样性与交互关系方面取得了初步的成果：关于海洋环境，参见例如 Venter et al. 2004; 关于石油环境，参见例如 Daniel, 2005; 关于人体环境，参见例如 Gi l l et al. 2006。

然而，当利用宏基因组学研究（例如， WGS策略)来分析环境样品中的微生物群落组成时，仍然存在两个巨大的挑战，即，大量的短基因片段（例如，测序片段）的组装以及不同物种的识别。由于宏基因组学研究所收集到的是，一个特定环境中的所有物种的基因信息，因此，如何将这些大量的、混合的短基因片段组装成连接片段或拼接片段，是一个巨大的难题与挑战。同时，在组装得到较长的连接片段或拼接片段后，如何判别这些长片段的物种来源，又是一个巨大的难题与挑战。

目前，已经开发出一些程序，用于组装混合的短基因片段，例如

Velvet (Zerbino and Birney 2008) , EULER-SR (Chaisson and Pevzner 2008) , Newbler (Mergul ies et al. 2006)和 Soapdenovo (Li et al. 2009)。此外，分栈（binning)方法已被广泛用于判别连接片段或拼接片段的归属物种，其包括但不限于，基于相似度 (s imi larity-based)的 MEGAN (Husonet al. 2007)和 CARMA (Tzahoret al. 2009) , 这类方法是通过与参考基因组进行序列比对来进行片段分类的；基于组成特征（compos ition- based)的分栈方法，例如基于 GC 含量、 k~mer频率（Schbath et al. 1995)或者四核苷酸频率（Teel ing et al. 2004)等的分栈方法，这类方法在很大程度上受限于片段长度和序列特征的辨别能力；以及，基于片段丰度（abundance- based)的 AbundanceBin (Wu and Ye 2011) , 这类方法才艮据环境中不同物种的丰度进行片段分类，只适合用于短的测序片段。

然而，宏基因组学的研究目的是重建环境样品中的各种微生物的基因组，以分析环境样品中的微生物群落组成。上述方法把组装和分栈分开，而各自只关注于一个方面。因此，上述方法并不能充分达到宏基因组学的研究目的。另外，即使将上述的组装方法和分栈方法简单组合在一起，由于不同方法所采用的算法、步骤、兼容性并不一定匹配，其最终结果是否能达到宏基因组学的研究目的，以及最终结果的精确度和有效性也都难以预料。

因此，本领域仍然需要一种高效率、高精度的分析环境样品中的微生物群落组成的方法。发明内容

在本发明中，除非另有说明，否则本文中使用的科学和技术名词具有本领域技术人员所通常理解的含义。并且，本文中所使用的各种实验室操作步骤均为相应领域内广泛使用的常规步骤。同时，为了更好地理解本发明，下面提供相关术语的定义和解释。

如本文中所使用的，术语 "环境" 是指广义上的各种环境，其包括但不限于，自然环境（例如土壤环境，海洋环境，河流环境）和体内环境（例如口腔环境，肠道环境）。更确切而言，术语 "环境" 是指，可能存在微生物 /微生物群落的任何区域。

如本文中所使用的，术语 "环境样品" 是指，来自各种环境的可能含有微生物 /微生物群落的样品。

如本文中所使用的，术语 "微生物" 具有本领域技术人员所通常理解的含义，其包括但不限于，细菌、真菌和病毒。

如本文中所使用的，术语 "微生物群落" 是指，在特定环境中生活在一起的各类微生物的联合。通常，同一微生物群落中的各种微生物不仅相互之间具有直接或间接的交互关系，而且与它们所生存的环境之间也具有交互关系：环境的改变会导致微生物群落的组成（包括，微生物的种类和 /或丰度）发生改变；反过来，微生物群落组成的改变也影响着环境。

如本文中所使用的，术语 "宏基因组" 是指，群落中的各种生物的基因组的总和。特别地，在本发明的方法和装置的背景下，术语 "宏基因组"是指，微生物群落中的各种微生物的基因组的总和。相应地，术语 "宏基因组测序数据"是指，对整个宏基因组进行测序所获得的数据。由于宏基因组所包含的 DNA信息十分巨大，因此，通常使用高通量测序技术 (例如，第二代测序技术或第三代测序技术）进行测序。然而，也可以通过其他方法或其他来源，获得所需的宏基因组测序数据。测序数据通常由大量的测序片段（read)构成。

第二代测序技术是本领域技术人员公知的，其包括例如， 454 测序法（Roche ) ， Solexa测序法 ( I l lumina ) ， S0UD测序法 ( ABI ) 和单分子测序法。关于第二代测序技术的详细综述，参见例如，

Michael Metzker (2010) , Sequencing technologies— the next generat ion, Nature Genet ics。关于第三代测序技术，可参见例如， Er ic E. Schadt等人， A window into third-generat ion sequencing, Human Molecular Genetics, 2010, Vol. 19， Review I ssue 2, R227- 240。

表述 "测序质量低的序列" 的含义是本领域技术人员已知的，其例如可在测序过程中由测序平台和测序软件确定（参见， Qual i ty Scores for Next-Generat ion Sequencing, Technical Note: Sequencing, I l lumina ) 。

如本文中所使用的，表述 "去冗余" 是指，对于彼此相似度达到 95%或以上的序列，只保留一个，例如，将重复的连接片段和拼接片段去除。

如本文中所使用的，表述 "参考集" 是广义上的组装片段集或基因集，其中，组装片段是指由测序片段组装得到的长片段，例如连接片段（cont igs ) 、拼接片段（scaffolds ) ；基因集是指在组装片段上预测出来的基因的集合。所述组装片段或基因构成且被称为参考集的 "元素" 。

如本文中所使用的，术语 "分栈（binning) " 和 "聚类" 具有相同的含义， "栈（bin) " 和 "类" 具有相同的含义。它们可互换使用。

如本文中所使用的，术语 "多元正态分布模型" 和 "最大似然函数法" 具有本领域技术人员所通常理解的含义。关于这 2个术语的详细描述，可参见例如 Fraley and Raf tery, 1998。

如本文中所使用的，术语 "基于相似度的聚类方法" 是指，通过比较两两序列之间的序列同一性来度量序列之间的相似度（或距离），并基于这个相似度（或距离）进行聚类；术语 "基于组成特征的聚类方法" 是指，通过比较两个序列自身组成特征的相似性，如寡核苷酸频率， GC含量等，来度量序列之间的相似度（或距离），并基于这个相似度（或距离 )进行聚类。基于相似度的聚类方法例如但不限于，基于相似度（simi lari ty- based)的 MEGAN (Husonet al. 2007)和 CARMA (Tzahoret al. 2009)。基于组成特征的聚类方法例如但不限于，基于 GC含量、 k- mer频率（Schbath et al. 1995)或者四核苷酸频率 (Tee l ing et al. 2004)的聚类方法。本发明所要解决的一个技术问题是，提供一种能有效分析环境样品中的微生物群落组成的方法和装置。基于此，发明人创造性地将组装方法和分栈方法结合在一起，开发了能够高效率且高精度地分析从环境样品获得的宏基因组数据，并进而确定环境样品的微生物群落组成的方法和装置。特别地，本发明的方法也被命名为 Soap series of Met a genome analysis (在下文中简称为 SoapMeta ) 。因此，在一个方面，本发明提供了一种用于分析环境样品中的微生物群落组成的方法，其包括以下步骤：

1 )测序：

对来自环境样品的基因组 DNA进行构建文库和测序，从而获得由测序片段池（reads pool )构成的宏基因组测序数据；

2 )初级组装：

2a )构建或完善参考集：对测序片段进行组装以得到组装片段，然后去冗余，从而构建非冗余的参考集（即，组装片段集）；任选地，可在所获得的组装片段上预测基因，并将预测出来的基因的集合作为参考集（即，基因集）；或者，如果针对所述环境样品，存在已知的参考集，那么直接将它作为参考集，或者将该已知的参考集与如上所述构建的参考集组合并去冗余，从而获得最终的参考集

2b )构建元素相对丰度讲矩阵：将所述测序片段与参考集进行比对，并计算参考集中的各个元素在样品中的相对丰度；

3 )分栈，即，通过下述步骤确定参考集中的每一个元素所归属的栈，得到聚类的栈：

3a )基于丰度的分栈：基于元素在样品中的相对丰度，使用聚类算法，如自底向上的层次聚类方法（HIERARCHICAL CLUSTERING SCHEMES, STEPHEN C. JOHNSON, 1967 ) ，确定各个元素的初始栈；和

3b )基于模型的分栈：

(i) 将每一个初始栈作为一个独立的多元正态分布模型，并基于丰度矩阵，利用最大似然函数法计算所述模型的参数；

(i i) 构建一个软矩阵（fuzzy matrix) , 用于存储每一个元素归属某一个栈的概率；和

(i i i) 迭代运算 E步和 M步，直至似然函数达到最大化： E 步，根据每一个栈的模型参数，分别计算每一个元素属于某一个栈的后验概率，并且修改软矩阵中所述元素属于所述栈的概率；

M 步：根据软矩阵，用最大似然函数法计算每一个栈的模型参数；

4 )基于栈的高级组装：

4a )通过将测序片段与已分栈的元素进行比对，从宏基因组测序数据中寻找对应到之前确定的各个栈的测序片段；

4b )使用 SOAPdenovo 或者使用其他针对微生物测序数据的组装软件，将对应到各个栈的测序片段分别进行组装；

4c )使用基于相似度的聚类方法和 /或基于组成特征的聚类方法，对每一个栈所包含的元素的分栈进行校正；任选地，还在已获得的栈内部进行再次聚类，然后根据聚类的结果，对已获得的栈进行拆分或保持不变，从而使结果更加准确可信；

4d )重复步骤 4a ) - 4c )，直到各个栈的基因组序列的大小无明显变化为止（总长度增长率小于 5% ) ；

5 )鉴定：

利用各个栈的基因组序列，确定各个栈所对应的微生物的类别，从而确定所述环境样品中的微生物群落组成。关于测序在一个优选的实施方案中，环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境。在另一个优选的实施方案中，环境样品来源于体内环境，例如口腔环境和肠道环境。

在一个优选的实施方案中，使用第二代测序技术（例如， 454 测序法， Solexa测序法， SOLiD测序法或单分子测序法）或第三代测序技术对环境样品所包含的微生物群落的宏基因组进行测序，从而提供来自环境样品的宏基因组测序数据。

在一个优选的实施方案中，通过下列步骤来获得宏基因组测序数据：

la )提供环境样品；

lb )从所述环境样品中提取宏基因组 DNA;

lc )利用所述宏基因组 DNA构基因组文库；

Id )对所述宏基因组文库进行测序，优选使用 Solexa 测序法进行测序，从而提供所述环境样品的宏基因组测序数据。

在一个优选的实施方案中，宏基因组测序数据是由测序片段构成的测序片段池（reads pool ) 。此类测序片段通常通过第二代测序技术（例如 Solexa测序法）或第三代测序技术获得。

在一个优选的实施方案中，测序片段是末端配对的测序片段 ( aired end reads ) 。

测序片段中可能包含测序过程中所使用的接头（adapter)的序列，测序质量低的序列和 /或在分析来自体内环境的样品的情况下，来自宿主基因组的序列。此类序列可能会影响后续的处理和分析，因此，此类序列的去除可能是有利的。

因此，在一个优选的实施方案中，在进行步骤 2 )之前，对测序数据进行预处理，即，去除接头序列、测序质量低的序列和 /或宿主基因组序列。

在一个优选的实施方案中，对来自相同或相似环境的多个样品进行测序，并将所有样品的测序数据组合在一起，构成宏基因组测序数据。在一个优选的实施方案中，宏基因组的测序深度为至少 10 χ，优选至少 20 X，优选至少 30 X，优选至少 40 X，更优选至少 50 χ 。关于初级组装

在一个优选的实施方案中，使用 Soapdenovo将所述测序片段组装成组装片段（例如，连接片段和 /或拼接片段）。此类组装方法是本领域技术人员已知的，参见例如， Li et al. 2009。

在一个优选的实施方案中，使用多个环境样品来进行本发明的方法，并且针对每个样品分别获得了各自的参考集。在此情况下，将所有样品的参考集组合在一起，并去冗余，从而构建最终的非冗余的参考集。也即，将来自多个样品的参考集组合在一起，并去冗余，从而构建最终的非冗余的参考集。

在一个优选的实施方案中，如果针对所述环境样品，存在已知的参考集，那么可以直接将它作为参考集，也可以将该已知的参考集与步骤 2a )中利用测序片段构建的参考集组合并去冗余，从而提供最终的参考集。

例如，在人脉道微生物群落的 MWAS研究中， Junj ie Qin et al. (2010) A human gut microbial gene catalogue establ ished by metagenomic sequencing. Nature, 464: 59-65已构建并公开了 3. 3M 欧洲人脉道微生物群落的非冗余基因集（即，参考集）。因此，在一个优选的实施方案中，所述环境样品是人脉道样品，并且将所述 3. 3M 欧洲人脉道微生物群落的非冗余基因集与步骤 2a )所构建的参考集组合并去冗余，从而提供最终的参考基因集。

在一个优选的实施方案中，通过使用 S0AP2或 MAQ比对软件，将所述测序片段与参考集进行比对。 S0AP2和 MAQ是本领域技术人员是已知的，参见例如， R Li et al. 2009和 Li et al. 2008。

在一个优选的实施方案中，使用 S0AP2将测序片段与参考集进行比对，并按照下列公式计算出参考集中的各元素的相对丰度： xJ L

a- = —— '■——

∑ ( A)，其中

A : 元素 i在样品中的相对丰度；

^L". 元素 i的长度；

^Xi：元素 i在样品中被检测到的次数。关于分栈

在一个优选的实施方案中，通过下列步骤来确定元素的初始栈：首先，基于元素在样品中的相对丰度，计算两两元素之间的相关性，例如 pearson相关系数， spearman相关系数， kendal l相关系数，欧几里得距离，曼哈顿距离等；然后，根据两两元素之间的相关性，通过聚类算法，如自底向上层次聚类等，将相关性密切的元素聚到一个类中，从而确定各个元素的初始栈。

在步骤 3 ) 的分栈之后，同一个栈里面的各个元素在所有样品中的丰度符合一定的分布模型，如正态分布。因此，聚到同一个栈里面的多个元素具有以下几种可能：（1) 这些元素属于同一个物种；（2) 这些元素来自共生的物种，因为共生物种的丰度分布相似；（3) 这些元素是几个物种共有的，因为几个物种共有的元素的丰度不同于每一个物种各自的丰度。关于基于栈的高级组装

在一个优选的实施方案中，使用 S0AP2来将测序片段与已分栈的元素进行比对。

在一个优选的实施方案中，使用 GC- depth spectra class if ier 和 /或 tetranucleotide frequencies (TNFs) class if ier (Teel ing et al. 2004)进行校正。关于鉴定在一个优选的实施方案中，通过将各个栈的基因组序列与已知的基因组数据库进行比对，从而确定各个栈所对应的微生物的类别。

在一个优选的实施方案中，所述基因组数据库包括但不限于，

NCBI/IMG已测序细菌库， CBI的 NR库等。

在一个优选的实施方案中，所述比对是核酸水平和 /或蛋白水平的比对。在另一个方面，本发明提供了一种用于分析环境样品中的微生物群落组成的装置，其包括以下模块：

1 )测序模块，其用于对来自环境样品的宏基因组 DNA进行测序，提供由测序片段池构成的宏基因组测序数据；

2 )初级组装模块，其与测序模块相连，且包括彼此相连的下列模块：

2a )组装构建模块，其用于对测序片段进行组装以得到组装片段，然后去冗余，从而构建非冗余的参考集（即，组装片段集）；任选地，所述组装构建模块还可在所获得的组装片段上预测基因，并将预测出来的基因的集合作为参考集（即，基因集）；和

2b ) 比对计算模块，其用于将测序片段与参考集进行比对，并计算参考集中的各个元素在样品中的相对丰度；

3 )分栈模块，其与初级组装模块相连，用于确定参考集中的每一个元素所归属的栈，得到聚类的栈，且包括彼此相连的下列模块：

3a )丰度分栈模块，其基于丰度确定各个元素的初始栈；和 3b )模型分栈模块，其基于模型确定各个元素所归属的栈；

4 ) 高级组装模块，其与测序模块和分栈模块相连，其用于从宏基因组测序数据中寻找对应到各个栈的测序片段，并将对应到各个栈的测序片段分别进行组装，且对组装结果进行^ £及调整；和

5 )鉴定模块，其与高级组装模块相连，用于通过各个栈的基因组序列，确定各个栈所对应的微生物的类别，从而确定所述环境样品中的微生物群落组成。在一个优选的实施方案中，环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境。在另一个优选的实施方案中，环境样品来源于体内环境，例如口腔环境和肠道环境。

在一个优选的实施方案中，所述测序模块使用第二代测序技术

(例如， 454测序法， Solexa测序法， SOLiD测序法或单分子测序法）或第三代测序技术对环境样品所包含的微生物群落的宏基因组进行测序，从而提供来自环境样品的宏基因组测序数据。

在一个优选的实施方案中，所述装置还包括彼此相连的 DNA提取模块和文库构建模块，其中，所述 DNA提取模块用于从所述环境样品中提取宏基因组 DNA, 并且，所述文库构建模块与测序模块相连，且利用所述宏基因组 DNA构建基因组文库。

在一个优选的实施方案中，所述测序模块所获得的测序片段是末端配对的测序片段 ( paired end reads ) 。

在一个优选的实施方案中，所述装置还包含过滤模块，其与测序模块和初级组装模块相连，用于在进行初级组装之前，去除测序片段中的接头序列、测序质量低的序列和 /或宿主基因组序列。

在一个优选的实施方案中，所述测序模块对宏基因组的测序深度为至少 10 ，优选至少 20 X，优选至少 30 ，优选至少 40 χ，更优选至少 50 X 。

在一个优选的实施方案中，所述组装构建模块使用 Soapdenovo 将测序片段组装成连接片段和 /或拼接片段。

在一个优选的实施方案中，所述组装构建模块还包含接收亚模块，其用于接收已知的参考集。在一个优选的实施方案中，所述组装构建模块将所接收的已知参考集作为最终的参考集。在另一个优选的实施方案中，所述组装构建模块将所接收的已知参考集与利用测序片段构建的参考集组合并去冗余，从而提供最终的参考集。

在一个优选的实施方案中，所述组装构建模块能够将来自多个样品的参考集组合在一起，并去冗余，从而构建最终的非冗余的参考集。

在一个优选的实施方案中，所述比对计算模块通过使用 S0AP2或 MAQ, 将测序片段与参考集进行比对。

在一个优选的实施方案中，所述比对计算模块使用 S0AP2将测序片段与参考集进行比对，并按照下列公式计算出参考集中各元素的相对丰度：

其中

元素 i在样品中的相对丰度；

^∑'：元素 i的长度；

^Xi：元素 i在样品中被检测到的次数。

在一个优选的实施方案中，所述丰度分栈模块基于元素在样品中的相对丰度，计算两两元素之间的相关性，然后通过聚类算法，确定各个元素的初始栈。

在一个优选的实施方案中，所述模型分栈模块通过下列来确定元素所归属的栈：

(i i i) 迭代运算 E步和 M步，直至似然函数达到最大化：

E 步，根据每一个栈的模型参数，分别计算每一个元素属于某一个栈的后验概率，并且修改软矩阵中所述元素属于所述栈的概率； M步：根据软矩阵，用最大似然函数法计算每一个栈的模型参数。在一个优选的实施方案中，所述高级组装模块通过下列来实现其功能：

( a )通过将测序片段与已分栈的元素进行比对，从宏基因组测序数据中寻找对应到所述分栈模块所确定的各个栈的测序片段；

( b M吏用 SOAPdenovo或者使用其他针对微生物数据的组装软件，将对应到各个栈的测序片段分别进行组装； ( c )使用基于相似度的聚类方法和 /或基于组成特征的聚类方法，对每一个栈所包含的元素的分栈进行校正；任选地，还在已获得的栈内部进行再次聚类，然后根据聚类的结果，对已获得的栈进行拆分或保持不变，从而使结果更加准确可信；

( d )重复步骤（a ) - (c ) ，直到各个栈的基因组序列的大小没有明显变化为止（总长度增长率小于 5% ) 。

在一个优选的实施方案中，所述高级组装模块使用 S0AP2来将测序片段与已分栈的元素进行比对。

在一个优选的实施方案中，所述高级组装模块使用 GC- depth spectra class if ier 和 /或 tetranucleotide frequencies (TNFs) classif ier进行校正。

在一个优选的实施方案中，所述鉴定模块通过将各个栈的基因组序列与已知的基因组数据库进行比对，从而确定各个栈所对应的微生物的类别。

在一个优选的实施方案中，所述基因组数据库包括，但不限于， NCBI/IMG已测序细菌库， CBI的 NR库等。

在一个优选的实施方案中，所述鉴定模块在核酸水平和 /或蛋白水平上进行比对。在另一个方面，还提供了本发明的装置用于分析环境样品中的微生物群落组成的用途。在一个优选的实施方案中，所述环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境。在另一个优选的实施方案中，环境样品来源于体内环境，例如口腔环境和肠道环境。发明的有益效果

本发明的方法和装置基于高通量测序技术，利用相同或相似环境下多个样品的测序数据进行组装，聚类和再组装，从而得到微生物群落的物种组成信息和物种的基因组信息，有着非常广泛的应用前景。与现有技术中的传统组装方法相比较，本发明的方法和装置有如下优点：

1、系统地将各种测序序列的属性结合起来，用于构建微生物群落的宏基因组的参考集，这特别适合于微生物物种分类，以及从来自同一环境的多个样品的测序数据重因组；

2、创造性地将分栈和组装有效地结合在一起，使物种基因组的组装结果更加精确，从而能够实现高效率、高精度地确定微生物群落的组成；

3、首次基于多个样品进行聚类分析，并进行了迭代高级组装。利用多个样品进行聚类分析具有有两个显著的优点： a )可以覆盖更多的低丰度物种，从而更全面地研究微生物群落； b ) 由于环境因素，不同的样品可能具有不同的物种组成和丰度，从而可以有利地进行比较研究。相比之下，利用单一样本进行的宏基因组学分析通常只能获得精确的优势物种，而无法全面地分析微生物群落，特别是低丰度物种（参见例如， Hess et al. 2011 ) 。

下面将结合附图和实施例对本发明的实施方案进行详细描述，但是本领域技术人员将理解，下列附图和实施例仅用于说明本发明，而不是对本发明的范围的限定。根据附图和优选实施方案的下列详细描述，本发明的各种目的和有利方面对于本领域技术人员来说将变得显然。附图说明

图 1示意性地描述了本发明的 SoapMeta方法的流程图，其中，虚线空心框、实线空心框和实心框示意性表示源自三个不同的物种。

图 2示意性地描述了本发明的 SoapMeta方法的初级组装的流程图。图 3是示意性地描述了本发明的 SoapMeta方法的分栈的流程图。图 4是示意性地描述了本发明的 SoapMeta方法的高级组装的¾½ 图。

图 5是描述了用于实施本发明的 SoapMeta方法的装置的结构示意图。图 6-8展示了实施例 2中利用第一种策略获得的 3个样品（样品 A-C ) 的 GC含量-测序深度讲图。图 6: 样品 A; 图 7: 样品 B; 图 8: 样品 C；。结果显示，样品 B和样品 C中的一些细菌很难区分，因为他们的 GC含量和测序深度非常接近。

图 9展示了本申请实施例 3中通过 16S rRNA测序获得的物种分类的信息图。

图 10展示了利用 16S rRNA测序法获得的 Akke ansia 16S rRNA标签的数量与利用本发明的 Soapmeta方法组装出来的相应基因组的测序深度的相关性。

图 11展示了利用 16S rRNA测序法获得的 Lactobacillus 16S rRNA标签的数量与利用本发明的 Soapmeta方法组装出来的相应基因组的测序深度的相关性。

图 10-11的结^ ^示，利用 16S rRNA测序法获得的 rRNA标签的数量与利用本发明的 Soapmeta方法组装出来的相应基因组的测序深度之间具有很强的相关性。这些结果表明，本发明的 Soapmeta方法的结果与 16S rRNA 测序法的结果是基本上一致的，再次证实了本发明的 SoapMeta方法的可靠性、准确性和高效性。具体实施方式

现参照下列意在举例说明本发明（而非限定本发明）的实施例来描述本发明。

除非特别指明，本发明中所使用的分子生物学实验方法，基本上参照 J. Sambrook等人，分子克隆：实验室手册，第 2版，冷泉港实验室出版社， 1989 , 以及 F. M. Ausubel等人，精编分子生物学实验指南，第 3版， John Wi ley & Sons, Inc. , 1995中所述的方法进行；并且各种酶的使用依照产品制造商推荐的条件。那些在实施例中未详细描述的过程和方法是本领域中公知的常规方法。本领域技术人员知晓，实施例以举例方式描述本发明，且不意欲限制本发明所要求保护的范围。实施例 1. 模拟环境样品的分析

1、数据模拟

为了模拟环境样品，我们从 NCBI基因组数据库（Wheeler et al. 2007)中选取了 100个不同的物种，这些物种的基因组从变形菌门中随机选择。另外，为了简化模型，不选择同一物种的不同品系。

我们一共模拟了 10例样品，每个样品的测序量均为 720 M。模拟的末端配对的测序片段的长度为 90bp，插入片段的大小为 500 ± 20bp (均值 ±标准差），测序错误率为 0. 1%。通过 Broken- Stick 模型 (MacArthur 1957)的相对物种丰度 (relative species abundance, RSA) , 来确定每一个样品的物种丰度组成比例。每一个样品所包含的大多数细菌的测序量是比较低的（64%的细菌的 RSA < 0. 01 ) 。将 10 个样品的测序数据合并后，这些低丰度细菌的测序量达 13. 6-182. 0 Mbp, 且测序深度为 2. 7 - 160. 4X。

2、初级组装

我们将所有样品的测序数据（测序片段）合并在一起，并使用组装软件 Soapdenovo (Li et al. 2009)进行初步的组装（即，不单独对每一个样品的测序数据进行分别的组装）。在组装后，对组装结果进行去冗余，从而得到非冗余的参考集。

特别地，在本实验中，混合样品的初级组装结果（即，参考集）共包含 41754条连接片段（contigs ) ，且连接片段的长度范围为 200- 2, 001, 157 bp (N50=93, 353bp) ( N50是衡量基因组图讲盾量的一个判断标准，其是指，当将所有的组装得到的序列按照长度从大到小排列，并从大到小将序列的长度相加，直至相加得到的总长度为所有组装得到的序列的总长度的百分之五十时，那条组装序列的长度，参见例如， Mi l ler et al. 2010. Assembly algorithms for next generation sequencing data. Genomics. 95 (6) : 315-327 ) 。将这些连接片段与原始细菌基因组进行 BLASTN比对。结果显示，组装后的连接片段对原始细菌基因组的平均覆盖度为 88. 7%, 并且每个细菌的覆盖度与测序深度呈现正相关，但是，当测序深度高于 20x时，参考集的覆盖度不再发生显著的变化。

使用 S0AP2 , 将测序片段与非冗余参考集进行比对，并通过下式计算出参考集中的各连接片段的相对丰度：

其中，

连接片段 i在样品中的相对丰度；

^L". 连接片段 i的长度；

χ'：连接片段 i在样品中被检测到的次数。

3、分栈（bin)

3. 1 基于丰度的分栈（初始分栈）

首先计算丰度矩阵中各连接片段的两两 Kendal l' s tau秩相关系数；然后根据连接片段两两之间的相关性，采用自底向上层次聚类算法，将相关性比较密切的片段聚到一个类中，从而获得初始的栈。

在本实验中，我们还使用默认的聚类参数，过滤掉了包含小于 10 个连接片段的初始栈，最终得到 343个初始栈。这些栈覆盖了 96. 8% 的连接片段（40， 438/41, 754)。

对于每一个初始栈，我们还给它定义一个属性， "最优的比对细菌"。也即，如果栈里面大部分的连接片段来自于某一个特定的细菌，那么这个细菌就是这个初始栈的最优的比对细菌。另外，还将栈的精度定为，来自最优的比对细菌的连接片段的总长度占栈里面的连接片段的总长度的百分比。在本实验中，初始栈的精度为 50. 3% - 100. 0% (平均值为 95. 1%)。

3. 2 基于模型的分栈

我们接着用基于模型的分栈方法来最优化初始分栈的结果。简言之， 1 )将每一个初始栈作为一个独立的多元正态分布模型，并基于丰度矩阵，利用最大似然函数法计算所述模型的参数；

2 )构建一个软矩阵（fuzzy matrix) , 用于存储每一个连接片段归属某一个栈的概率； 3 )迭代运算 E步和 M步，直至似然函数达到最大化：

E步，根据每一个栈的模型参数，分别计算每一个连接片段属于某一个栈的后验概率，并且修改软矩阵中所述连接片段属于所述栈的概率；

M步：根据软矩阵，用最大似然函数法计算每一个栈的模型参数。

在该步骤后，所获得的栈减少到 135个。与初始分栈相比，这些栈的覆盖度下降到 91. 9% ( 38, 364/41, 754个连接片段），且精度下降到 33. 2% - 100. 0% (平均值 92. 3%)。在这 135个栈中，每一个栈代表一个物种。基于各个栈中的连接片段的序列，我们鉴定到了 86个物种（86% ) ，且每一个物种的基因组覆盖度超过 50%。

4、高级组装

高级组装分成以下 3步：

1 )使用 S0AP2, 通过序列比对，在模拟的测序数据中寻找对应到之前确定的各个栈的测序片段；

2 )使用 SOAPdenovo分别将对应到各个栈的测序片段进行深度组装；

3 )使用基于相似度和组成特征的聚类方法，对每一个栈所包含的连接片段的分栈进行校正，并且在已有的栈内部进行再次聚类，然后根据再次聚类的结果，对已有的栈进行拆分或保持不变，从而使结果更加准确可信；

4 )重复步骤 1 ) - 3 )，直到各个栈的基因组序列的大小没有明显变化为止（总长度增长率小于 5% ) 。

对之前获得的 135个栈进行高级组装之后，得到 148个经组装的栈。栈的数目的增加是因为，我们使用了基于组成特征的聚类方法，根据 GC含量，测序深度等特征，将一个栈里面的一些可以明显再细分的栈拆开了。

在高级组装后，栈的平均精度达到 94. 2% (参见，表 1)，略微高于前一步的结果。另外，当用原始细菌基因组覆盖组装的栈的基因组时，结果显示，覆盖度为 95. 5%; 反之，当用组装的栈的基因组覆盖原始细菌基因组时，覆盖度为 57. 4%。

在这 148个栈中，基于各个栈的组装的基因组序列，我们鉴定到了 100个初始细菌物种中的 95个（95% ) ，且如上所述，每一个物种的基因组覆盖度超过 50%。

上述结果表明，本发明的 SoapMeta方法的特异度较好，且能够有效地鉴别出模拟样品中所包含的绝大部分物种（95% ) 。

表 1. 每一步骤得到的栈的比较

覆盖的初始连连接片段的平覆盖的初始细栈的数量

接片段 (%) 均精度 (%) 菌种类 (%) 初始栈 343 96. 8 95. 1 90

EM迭代后的栈 135 91. 9 92. 3 86 高级組装后的栈 148 - 94. 2 95 实施例 2. 简单环境样品（纤维素降解菌群）的分析

本实施例以一个真实的简单环境为例，对本发明的 SoapMeta方法进行了进一步的解释说明，并且通过与传统的分析方法相比较，证实了本发明的 SoapMeta方法的优势。在本实施例中，我们收集了三个样品（样品 A、 B、 C ) ，它们分别来自不同培养^^下的纤维素降解菌群：从同一沼泽的土壤采集 3 个样品，并且分别用三种包含不同碳源（滤纸、纤维二糖、葡萄糖）的培养基在 37 Ό下培养 52小时，然后分别收获菌体，从而获得样品 A、 B、 Co针对每一个样品，我们分别构建了一个测序文库（参数设置：末端配对的测序片段的长度为 90bp, 插入片段的大小为 500 ± 20bp) : 首先用 HiSeq2000 对样品进行测序，从而得到原始测序片段（raw reads)；然后，过滤掉其中的低质量序列和接头序列，从而提供 3. 88Gb 的用于分析的宏基因组测序数据 ( 3个样品的测序数据的总和）。在在本本实实施施例例中中，，我我们们应应用用了了两两种种策策略略来来构构建建微微生生物物的的基基因因组组。。第第一一种种策策略略是是，，使使用用传传统统的的分分析析方方法法，，对对每每个个样样品品分分别别进进行行测测序序数数据据的的组组装装，，从从而而构构建建微微生生物物的的基基因因组组（（参参见见，， MMEEGGAANN ((HHuussoonneett aall.. 22000077)) ))；；第第二二种种策策略略是是，，使使用用本本发发明明的的 SSooaappMMeettaa方方法法，，将将所所有有样样品品的的测测序序数数据据混混合合在在一一起起，，然然后后进进行行初初级级组组装装，，分分栈栈和和高高级级组组装装，，从从而而构构建建微微生生物物的的基基因因组组。。将将第第一一种种策策略略用用作作对对照照，，以以证证实实本本发发明明的的 SSooaappMMeettaa方方法法在在多多个个样样品品的的混混合合组组装装方方面面的的优优势势。。

在在第第一一种种策策略略下下，，用用基基于于组组成成特特征征的的聚聚类类方方法法对对来来自自单单个个样样品品的的测测序序片片段段进进行行聚聚类类，，以以判判别别样样品品中中潜潜在在的的微微生生物物。。对对于于所所使使用用的的 33个个样样品品，，我我们们分分别别得得到到了了 66个个类类 ((样样品品 AA))，， 22个个类类 ((样样品品 BB))，，和和 33个个类类 ((样样品品 CC))。。这这 33个个样样品品各各自自的的 GGCC图图（（参参见见图图 66--88 ))显显示示，，样样品品 BB和和样样品品 CC 中中的的一一些些细细菌菌很很难难区区分分，，因因为为他他们们的的 GGCC含含量量和和测测序序深深度度非非常常接接近近。。

在在第第二二种种策策略略下下，，我我们们首首先先在在初初级级组组装装中中得得到到了了连连接接片片段段的的相相对对丰丰度度。。进进一一步步，，通通过过使使用用本本发发明明的的 SSooaappMMeettaa方方法法，，我我们们从从 33个个样样品品的的混混合合测测序序数数据据中中鉴鉴定定到到了了 1100个个栈栈，，其其中中有有 99个个栈栈的的组组装装的的基基因因组组序序列列大大于于 llMMbbpp,, 并并且且这这 1100个个栈栈的的基基因因组组序序列列总总长长覆覆盖盖了了所所有有样样品品测测序序数数据据的的 8899.. 55%%。。在在这这 1100个个栈栈中中，，每每一一个个栈栈对对应应一一个个潜潜在在的的物物种种。。随随后后，，我我们们对对每每个个栈栈的的组组装装的的基基因因组组序序列列进进行行了了 TTBBLLAASSTTXX比比对对，，以以确确定定各各个个栈栈所所对对应应的的潜潜在在的的物物种种，，结结果果见见表表 22。。

表表 22的的结结果果显显示示，，在在这这 1100个个栈栈中中，，有有 66个个栈栈的的组组装装的的基基因因组组序序列列很很纯纯（（即即，，基基本本上上对对应应至至同同一一个个微微生生物物物物种种的的基基因因组组））：：短短短短芽芽孢孢杆杆菌菌 NNBBRRCC 110000559999 ((BBrreevvii bbaacciilllluuss bbrreevviiss WWBBRRCC 1100005599 ))、、救救、、结结簧簧孢孢杆杆菌菌 22——66 ((BBaacciilllluuss ccooaagguullaannss 22--66))、、耐耐盐盐芽芽抱抱杆杆菌菌 CC—— 112255 {{BBaacciilllluuss hhaalloodduurraannss CC--1122SS)) ,,肉肉毒毒梭梭菌菌 ΑΑ22 KKyyoottoo CClloossttrriiddiiuumm bboottuulliinnuumm AA22 KKyyoottoo)) ,,热热解解纤纤维维梭梭菌菌 AATTCCCC 2277440055 {{CClloossttrriiddiiuumm tthheemmoocceelllluumm AATTCCCC 22774400$$)) ,,热热解解纤纤维维梭梭菌菌 AATTCCCC 2277440055 {{CClloossttrriiddiiuumm tthheemmoocceelllluumm AATTCCCC 227744 OOSS))，，、、CClloossttrriiddiiwwnn tthheerrmmoocceelllluuiinn))

菌 (Weimer and Zeikus 1977; Bayer et al. 1983; 和 Schwarz 2001)。此外，其中的短芽抱杆菌（5reW6a '//i/5 和芽抱杆菌（ c///i/5 也已知具有纤维降解能力（Liang et al. 2009; Li et al. 2006; 和 Rastogi et al. 2009)。

从上面的结果可知，本发明的 SoapMeta策略不仅在精度和覆盖度上显著优于第一种策略（即，基因组覆盖度更全，分类准确度更高），而且能够更有效、更精确地鉴定环境样品的微生物组成。表 2、纤维素降解菌群的组装基因组总表

^品^序度相当于第一神策略 '吳

我连接片段总长度覆盖.

最^近细 §基¾ i

數量手均相似度 )

未知种

未知物种

未知^种

^知种注：图中的 *表示，该栈包含有多个物种的序列，并且无法进一步明确区分。例如， B1*表示，栈 B1中含有无法进一步区分的多个物种的序列（在使用第二种策略的方法中，栈 B1 中的这些物种被进一步区分为短短芽孢杆菌 NBRC 100599和热解纤维梭菌 ATCC 27405 ) 。实施例 3. 复杂环境样品（小鼠肠道菌群）的分析

本实施例以一个真实的复杂环境为例，示例性地展示了本发明的 SoapMeta方法在小鼠肠道菌群的探测中的应用。本实验采用了两种常见的小鼠， SV- 129和 C57Black/6 (Fuj i i et al. 1997)。在现实中，小鼠肠道的菌群的相对丰度会随着年龄，性别，饮食等等因素的变化而变化，但是如果小鼠的饮食固定，且环境固定的话，这些菌群的微生物组成一般不会有太大的变动。因此，可以利用本发明的 SoapMeta方法来研究特定环境、特定饮食下小鼠的肠道菌群的微生物组成，并构建菌群物种的基因组。

收集了 13个粪便样品（其中 6个样品来自 SV- 129小鼠， 7个样品来自 C57Black/6小鼠），并构建了测序文库（参数设置：末端配对的测序片段的长度为 90bp，插入片段的大小为 350 ± 15bp) : 首先用 HiSeq2000对样品进行测序，从而得到原始测序片段（raw reads)；然后，过滤掉其中的低盾量序列、接头序列以及小鼠基因组序列，从而获得 3. 96 ± 0. 55Gbp (每个样品的平均测序数据）的用于分析的宏基因组测序数据。

根据本发明的 SoapMeta方法：

首先，对样品的宏基因组测序数据进行了初级组装，得到 246. IMb 的连接片段集（n=180， 056个， 50=2, 613bp)；

然后，进行了分栈，得到 325个栈（将序列含量低于 lOOKbp的栈过滤掉），这些栈的总序列含量为 213. 6Mbp (86. 8%) , 并且其中有 56个栈的序列含量大于 IMbp;

最后，对上述序列含量大于 IMbp的 56个栈进行了高级组装，最终得到 ⁵7个基因组（栈），其总序列含量达 14L 6Mbp (每个基因组的平均序列含量为 2. 48Mbp ) ，并且覆盖了 49. 5%的测序片段。结果概述于表 3中。

使用 BLASTN (核酸水平）和 TBLASTX (蛋白水平），将高级组装得到的栈与已知的基因组数据库进行比对。结果显示，有 8个栈在核酸水平上与已知的物种十分接近：它们均具有高于 90%的精度和高于 95%的序列相似度。此外，还有 48个栈在蛋白水平上与已知的物种高度同源：它们均具有高于 70%的精度和高于 50%的序列相似度。另外，还有 1个栈比对到未知的物种。表 3: ⁵7个基因组（栈）的详细信息

为了验证上述结果，我们通过 Solexa 测序法对这些样品的 16S rRNA ( V6高变区）进行了测序，得到高质量的 3· 63 ± 0· 68Μ (均值士标准差）的 16S rRNA标签（tags) (已过滤掉接头序列，低盾量序列，重叠序列和引物序列）。利用 BLASTN, 将这些 16S rRNA标签与 RefSSU 数据库（Huse et al. 2010)进行比对。结果示于图 9中。结果显示，小鼠肠道菌群中，丰度较高的微生物是：職 ^^Lachnospiracea )、 #L^f S>¾ (Lactobacillus) . j5'J (Allobaculum) . 阿克曼氏菌属 (Akker ansia)、 Ruminococca ea L头杆菌属 {Papillibacter)、拟杆菌属 {Bacteroides)和脱疏 ¾ 菌科 {Desulfo vi briona cea e ) 。这些细菌大部分能够被本发明的 SoapMeta方法组装出来的基因组覆盖，这充分表明，本发明的 SoapMeta方法能够高效、精确地鉴定环境样品中的微生物组成。

另夕卜，我们还将 Akkermansia属和 Lactobacillus 16S rR A 标签的数量与用 Soapmeta方法组装出来的基因组的测序深度做比较。结果显示，它们之间具有很强的相关性（参见图 10-11)。这再次表明了本发明的 SoapMeta方法的准确性和高效性。尽管本发明的具体实施方式已经得到详细的描述，但本领域技术人员将理解：根据已经公开的所有教导，可以对细节进行各种修改和变动，并且这些改变均在本发明的保护范围之内。本发明的全部范围由所附权利要求及其任何等同物给出。

Claims

利要求

1. 一种用于分析环境样品中的微生物群落组成的方法，其包括以下步骤：

1 )测序：

对来自环境样品的基因组 DNA进行构建文库和测序，从而获得由测序片段池构成的宏基因组测序数据；

2 )初级组装：

2a )构建或完善参考集：对测序片段进行组装以得到组装片段，然后去冗余，从而构建非冗余的参考集；或者，可在所获得的组装片段上预测基因，并将预测出来的基因的集合作为参考集；或者，如果针对所述环境样品，存在已知的参考集，那么直接将它作为参考集，或者将该已知的参考集与如上所述构建的参考集组合并去冗余，从而获得最终的参考集；

2b )构建元素相对丰度讲矩阵：将所述测序片段与参考集进行比对，并计算参考集中的各个元素在样品中的相对丰度；

3 )分栈，即，通过下述步骤确定参考集中的每一个元素所归属的栈，得到聚类的栈：

3a )基于丰度的分栈：基于元素在样品中的相对丰度，使用聚类算法，例如自底向上的层次聚类方法，确定各个元素的初始栈；和

3b )基于模型的分栈：

(i) 将每一个初始栈作为一个独立的多元正态分布模型，并基于丰度矩阵，利用最大似然函数法计算所述模型的参数；

(ϋ) 构建一个软矩阵，用于存储每一个元素归属某一个栈的概率；和

(i i i) 迭代运算 E步和 M步，直至似然函数达到最大化： E步，根据每一个栈的模型参数，分别计算每一个元素属于某一个栈的后验概率，并且修改软矩阵中所述元素属于所述栈的概率；

M步：根据软矩阵，用最大似然函数法计算每一个栈的模型参数；

4 )基于栈的高级组装： 4a )通过将测序片段与已分栈的元素进行比对，从宏基因组测序数据中寻找对应到之前确定的各个栈的测序片段；

4b )使用 SOAPdenovo或者使用其他针对微生物测序数据的组装软件，将对应到各个栈的测序片段分别进行组装；

4c )使用基于相似度的聚类方法和 /或基于组成特征的聚类方法，对每一个栈所包含的元素的分栈进行校正；任选地，还在已获得的栈内部进行再次聚类，然后根据聚类的结果，对已获得的栈进行拆分或保持不变；

4d )重复步骤 4a ) - 4c ) , 直到各个栈的基因组序列的大小无明显变化为止（总长度增长率小于 5% ) ；

5 )鉴定：

利用各个栈的基因组序列，确定各个栈所对应的微生物的类别，从而确定所述环境样品中的微生物群落组成。
2. 权利要求 1的方法，其中，

例如，所述环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境；或者所述环境样品来源于体内环境，例如口腔环境和肠道环境；

例如，在步骤 1 )中使用第二代测序技术（例如， 454测序法， Solexa 测序法， SOLiD测序法或单分子测序法）或第三代测序技术对环境样品所包含的微生物群落的宏基因组进行测序，从而提供来自环境样品的宏基因组测序数据；

例如，在步骤 1 ) 中通过下列步骤来获得宏基因组测序数据：

la )提供环境样品；

lb )从所述环境样品中提取宏基因组 DNA;

lc )利用所述宏基因组 DNA构基因组文库；

Id )对所述宏基因组文库进行测序，优选使用 Solexa测序法进行测序，从而提供所述环境样品的宏基因组测序数据；

例如，所述测序片段是末端配对的测序片段；例如，在进行步骤 2 )之前，对测序数据进行预处理，即，去除接头序列、测序质量低的序列和 /或宿主基因组序列；

例如，对来自相同或相似环境的多个样品进行测序，并将所有样品的测序数据组合在一起，构成宏基因组测序数据；

例如，宏基因组的测序深度为至少 10 x，优选至少 20 χ，优选至少 30 X，优选至少 40 ，更优选至少 50 。
3. 权利要求 1或 2的方法，其中，在步骤 2 ) 中

例如，使用 Soapdenovo将所述测序片段组装成组装片段（例如，连接片段和 /或拼接片段）

例如，将来自多个样品的参考集组合在一起，并去冗余，从而构建最终的非冗余的参考集；

例如，使用 S0AP2或 MAQ比对软件，将所述测序片段与参考集进行比对；

例如，使用 S0AP2将测序片段与参考集进行比对，并按照下列公式计算出参考集中的各元素的相对丰度：

其中

元素 i在样品中的相对丰度；

^∑'：元素 i的长度；

^Xi：元素 i在样品中被检测到的次数。
4. 权利要求 1-3任一项的方法，其中，

例如，在步骤 3 ) 中，通过下列步骤来确定元素的初始栈：首先，基于元素在样品中的相对丰度，计算两两元素之间的相关性，例如 pearson相关系数， spearman相关系数， kendal l相关系数，欧几里得距离，曼哈顿距离等；然后，根据两两元素之间的相关性，通过聚类算法，如自底向上层次聚类等，将相关性密切的元素聚到一个类中，从而确定各个元素的初始栈；

例如，在步骤 4 ) 中，使用 S0AP2来将测序片段与已分栈的元素进行比对；

例如，在步骤 4 ) 中，使用 GC- depth spectra class if ier和 /或 tetranucleotide frequencies (TNFs) classif ier进行校正。
5. 权利要求 1-4任一项的方法，其中，在步骤 5 ) 中，

例如，通过将各个栈的基因组序列与已知的基因组数据库进行比对，从而确定各个栈所对应的微生物的类别；

例如，所述基因组数据库选自下列： NCBI/IMG已测序细菌库， NCBI 的 NR库;

例如，所述比对是核酸水平和 /或蛋白水平的比对。
6. 一种用于分析环境样品中的微生物群落组成的装置，其包括以下模块：

1 )测序模块，其用于对来自环境样品的宏基因组 DNA进行测序，提供由测序片段池构成的宏基因组测序数据；

2 )初级组装模块，其与测序模块相连，且包括彼此相连的下列模块：

2a )组装构建模块，其用于对测序片段进行组装以得到组装片段，然后去冗余，从而构建非冗余的参考集；任选地，所述组装构建模块还可在所获得的组装片段上预测基因，并将预测出来的基因的集合作为参考集；和

2b )比对计算模块，其用于将测序片段与参考集进行比对，并计算参考集中的各个元素在样品中的相对丰度；

3 )分栈模块，其与初级组装模块相连，用于确定参考集中的每一个元素所归属的栈，得到聚类的栈，且包括彼此相连的下列模块：

3a )丰度分栈模块，其基于丰度确定各个元素的初始栈；和 3b )模型分栈模块，其基于模型确定各个元素所归属的栈； 4 ) 高级组装模块，其与测序模块和分栈模块相连，其用于从宏基因组测序数据中寻找对应到各个栈的测序片段，并将对应到各个栈的测序片段分别进行组装，且对组装结果进行^ £及调整；和

5 )鉴定模块，其与高级组装模块相连，用于通过各个栈的基因组序列，确定各个栈所对应的微生物的类别，从而确定所述环境样品中的微生物群落组成。
7. 权利要求 6的装置，其中，

例如，所述环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境；或者，所述环境样品来源于体内环境，例如口腔环境和肠道环境；

例如，所述测序模块使用第二代测序技术（例如 454测序法， Solexa 测序法， SOLiD测序法或单分子测序法）或第三代测序技术对环境样品所包含的微生物群落的宏基因组进行测序，从而提供来自环境样品的宏基因组测序数据；

例如，所述装置还包括彼此相连的 DNA提取模块和文库构建模块，其中，所述 DNA提取模块用于从所述环境样品中提取宏基因组 DNA, 并且，所述文库构建模块与测序模块相连，且利用所述宏基因组 DNA构建基因组文库；

例如，所述测序模块所获得的测序片段是末端配对的测序片段；例如，所述装置还包含过滤模块，其与测序模块和初级组装模块相连，用于在进行初级组装之前，去除测序片段中的接头序列、测序质量低的序列和 /或宿主基因组序列；

例如，所述测序模块对宏基因组的测序深度为至少 10 X，优选至少

20 X，优选至少 30 X，优选至少 40 X，更优选至少 50 χ。
8. 权利要求 6或 7的装置，其中，

例如，所述组装构建模块使用 Soapdenovo将测序片段组装成连接片段和 /或拼接片段；例如，所述组装构建模块还包含接收亚模块，其用于接收已知的参考集；优选地，所述组装构建模块将所接收的已知参考集作为最终的参考集，或者将所接收的已知参考集与利用测序片段构建的参考集组合并去冗余，从而提供最终的参考集；

例如，所述组装构建模块能够将来自多个样品的参考集组合在一起，并去冗余，从而构建最终的非冗余的参考集；

例如，所述比对计算模块通过使用 S0AP2或 MAQ, 将测序片段与参考集进行比对；

例如，所述比对计算模块使用 S0AP2将测序片段与参考集进行比对，并按照下列公式计算出参考集中各元素的相对丰度：

其中

元素 i在样品中的相对丰度；

^∑'：元素 i的长度；

^Xi：元素 i在样品中被检测到的次数。
9. 权利要求 6-8任一项的装置，其中，

例如，所述丰度分栈模块基于元素在样品中的相对丰度，计算两两元素之间的相关性，然后通过聚类算法，确定各个元素的初始栈；

例如，所述模型分栈模块通过下列来确定元素所归属的栈：

(i) 将每一个初始栈作为一个独立的多元正态分布模型，并基于丰度矩阵，利用最大似然函数法计算所述模型的参数；

(ϋ) 构建一个软矩阵，用于存储每一个元素归属某一个栈的概率；和

(i i i) 迭代运算 E步和 M步，直至似然函数达到最大化：

E步，根据每一个栈的模型参数，分别计算每一个元素属于某一个栈的后验概率，并且修改软矩阵中所述元素属于所述栈的概率；

M步：根据软矩阵，用最大似然函数法计算每一个栈的模型参数；例如，所述高级组装模块通过下列来实现其功能：

( a )通过将测序片段与已分栈的元素进行比对，从宏基因组测序数据中寻找对应到分栈模块所确定的各个栈的测序片段；

( b )使用 SOAPdenovo或者使用其他针对微生物数据的组装软件，将对应到各个栈的测序片段分别进行组装；

( c )使用基于相似度的聚类方法和 /或基于组成特征的聚类方法，对每一个栈所包含的元素的分栈进行校正；任选地，还在已获得的栈内部进行再次聚类，然后才艮据聚类的结果，对已获得的栈进行拆分或保持不变；

( d )重复步骤（a ) - (c ) , 直到各个栈的基因组序列的大小没有明显变化为止（总长度增长率小于 5% ) ；

例如，所述高级组装模块使用 S0AP2来将测序片段与已分栈的元素进行比对；

例如，所述高级组装模块使用 GC- depth spectra class if ier和 / tetranucleotide frequencies (TNFs) class if ier进行校正；例如，所述鉴定模块通过将各个栈的基因组序列与已知的基因组数据库进行比对，从而确定各个栈所对应的微生物的类别；优选地，所述基因组数据库选自 NCBI/IMG已测序细菌库和 /或 NCBI的 NR库；优选地，所述鉴定模块在核酸水平和 /或蛋白水平上进行比对。
10. 权利要求 6-9任一项的装置用于分析环境样品中的微生物群落组成的用途，其中，

例如，所述环境样品来源于自然环境，例如土壤环境，海洋环境和河流环境；或者，所述环境样品来源于体内环境，例如口腔环境和肠道环境。