一种基于宏基因组测序数据组装病原微生物基因组的方法

摘要

本发明涉及一种基于宏基因组测序数据组装病原微生物基因组的方法。该方法将原始数据过滤后与宿主数据库进行比对，以达到在数据层面进行去除宿主序列的目的；再使用soapdenovo对去宿主的reads进行组装，得到无参组装的contig序列，将病源数据库作为参考基因组，然后统计去宿主的reads比对情况中每个位点的测序深度，得到有参组装的contig，将无参组装的contig序列和有参组装的contig序列进行整合，得到合并后的contig序列，并进行病原的判别。本发明通过去宿主后，采取无参组装和有参组装相结合的方法进行病原微生物的组装，得到的宏基因组没有宿主的污染，并且在点突变以外加入了结构变异的信息，准确度更高。

权利要求

1.一种基于宏基因组测序数据组装病原微生物基因组的方法，其特征在于，具体步骤如下：

（1）使用fastp软件对二代下机数据进行质控，得到质控后高质量的二代数据，具体方法为：根据实验信息与数据的实际情况，使用fastp软件对原始下机数据进行过滤，除去接头序列，低质量序列以及长度不符合要求的序列，得到符合质控要求的下机数据；

（2）使用bwa软件将质控后高质量的二代数据比对到宿主数据库，得到比对情况，具体方法为：质控后的下机数据中同时存在宿主的DNA序列以及病原微生物的DNA序列，使用bwa软件，以宿主物种的基因组作为参考基因组，将质控后的下机数据与参考基因组进行比对，获取下机数据中所有序列的比对情况，即all_bam文件；

（3）使用samtools软件，根据步骤（2）中比对情况提取未比对上的数据，得到去宿主的reads，具体方法为：根据步骤（2）中得到的比对情况，使用samtools软件对all_bam文件中所有序列的比对情况进行过滤，对比对到参考基因组也就是宿主基因组的序列进行过滤，提取未比对上的所有序列，得到最终的去宿主的病原微生物的序列集合，即Clean reads；

（4）使用soapdenovo软件对去宿主的reads进行组装，得到无参组装的contig序列，具体方法为：使用soapdenovo软件对Clean reads进行组装，得到不依赖病原微生物参考基因组的contig序列；

（5）使用bwa软件将去宿主的reads比对到病原数据库，得到去宿主的reads比对情况，具体方法为：使用bwa软件，以病原微生物的基因组作为参考基因组，将去宿主的reads比对到参考基因组上，获得去宿主的reads在病原微生物参考基因组上的比对情况，即pathogen_bam文件；

（6）使用samtools软件统计去宿主的reads比对情况中每个位点的测序深度，得到有参组装的contig序列，具体方法为：使用samtools软件统计步骤（5）中pathogen_bam文件每个位点的测序深度，将每个位置上占比最高的碱基确定为该位置的碱基序列，并依次进行延伸；对于比对中出现的indel，根据位点两端序列以及位置信息，保留reads中的indel，得到有参组装的contig序列；

（7）将无参组装的contig序列和有参组装的contig序列进行整合，得到合并后的contig序列，具体方法为：使用自主编写脚本，以有参组装的contig序列为骨架，根据无参组装contig序列中的大片段易位或者倒位信息对骨架序列进行调整或重排，从而将无参组装的contig序列和有参组装的contig序列进行合并，得到合并后的contig序列；

（8）使用去宿主的reads，以合并后的contig序列为参考进行比对，并重复步骤（6），继续进行有参组装，得到最终版的contig序列，具体方法为，使用bwa软件，以步骤（7）中得到的合并的contig序列作为参考基因组，将Clean reads比对到参考基因组上，得到最终的比对情况，即final_bam文件，并参照步骤（6）中的方法，对得到的final_bam文件进行二次有参组装，得到最终版的contig序列；

（9）使用blast软件，将最终的contig序列与已有的病原微生物数据库进行比对，根据比对结果中序列之间的相似性与一致性参数，对contig序列进行病原属种和亚型鉴定。

说明书

一种基于宏基因组测序数据组装病原微生物基因组的方法

技术领域

[0001] 本发明涉及宏基因组测序与组装领域，具体涉及一种基于宏基因组测序数据组装病原微生物基因组的方法。

背景技术

[0002] 微生物群体几乎存在于这个世界每一个生态群落之中，从个体体表到肠道，从高原空气到深海海底淤泥，从冰川冻湖到火山岩浆都无处不在，并扮演着不可或缺的角色。对微生物的研究从列文虎克发明显微镜开始的数百年中，主要基于纯培养的研究方式。在数以万亿计的微生物种类中，仅0.1％～1％的物种可培养，极大地限制了对微生物多样性资源的研究和开发。

[0003] Metagenomics(翻译成元基因组学，或者翻译成宏基因组学)，是一种直接对微生物群体中包含的全部基因组信息进行研究的手段，即绕过对微生物个体进行分离培养，应用基因组学技术对自然环境中的微生物群落进行研究的学科。它规避了对样品中的微生物进行分离培养，提供了一种对不可分离培养的微生物进行研究的途径，更真实的反应样本中微生物组成、互作情况，同时在分子水平对其代谢通路、基因功能进行研究。

[0004] 近年来，随着测序技术和信息技术的快速发展，利用新一代测序技术(Next Generation Sequencing)，即二代测序技术研究Metagenomics，能快速准确的得到大量生物数据和丰富的微生物研究信息，从而成为研究微生物多样性和群落特征的重要手段。如致力于研究微生物与人类疾病健康关系的人体微生物组计划(HMP,Human Microbiome Project)，研究全球微生物组成和分布的全球微生物组计划(EMP,Earth Microbiome Project)都主要利用高通量测序技术进行研究。

[0005] 宏基因组的二代测序数据常用的组装软件为soapdenovo，将病原微生物组装到contig水平。样本中由于可能存在大量宿主细胞或宿主的游离核酸，病原微生物核酸在提取出的核酸样本中的丰度通常极低。因此宏基因组测序数据组装之后的序列会含有宿主的contig序列，导致组装准确度降低，影响宏基因组的后续分析。

发明内容

[0006] 本发明所要解决的技术问题就是提供一种基于宏基因组测序，通过对测序数据质控、与宿主数据库比对、去除宿主序列、进行有参和无参组装，进而得到较准确的宏基因组组装结果的方法,即提供一种基于宏基因组测序数据组装病原微生物基因组的方法。采用的技术方案为：

[0007] 一种基于宏基因组测序数据组装病原微生物基因组的方法，具体步骤如下：

[0008] (1)使用fastp对二代下机数据进行质控，得到质控后高质量的二代数据；

[0009] (2)使用bwa将质控后高质量的二代数据比对到宿主数据库，得到比对情况；

[0010] (3)使用samtools软件，根据步骤(2)中比对情况提取未比对上的数据，得到去宿主的reads，即Clean reads；

[0011] (4)使用soapdenovo对去宿主的reads进行组装，得到无参组装的contig序列；

[0012] (5)使用bwa将去宿主的reads比对到病原数据库，得到去宿主的reads比对情况；

[0013] (6)使用samtools软件统计去宿主的reads比对情况中每个位点的测序深度，得到有参组装的contig序列；

[0014] (7)将无参组装的contig序列和有参组装的contig序列进行整合，得到合并后的contig序列；

[0015] (8)使用去宿主的Clean reads，以合并后的contig序列为参考进行比对，并重复步骤6，继续进行有参组装，得到最终版的contig序列；

[0016] (9)使用blast软件，将最终的contig序列与已有的病原微生物数据库进行比对，根据比对结果中序列之间的相似性与一致性等参数，对contig序列进行病原属种和亚型鉴定。

[0017] 优选的，所述步骤(1)的具体方法为，根据实验信息与数据的实际情况，使用fastp软件对原始下机数据进行过滤，除去接头序列，低质量序列以及长度不符合要求的序列，得到符合质控要求的下机数据。

[0018] 优选的，所述步骤(2)的具体方法为，质控后的下机数据中同时存在宿主的DNA序列以及病原微生物的DNA序列，使用bwa软件，以宿主物种的基因组作为参考基因组，将质控后的下机数据与参考基因组进行比对，获取下机数据中所有序列的比对情况，即all_bam文件。

[0019] 优选的，所述步骤(3)的具体方法为，根据步骤(2)中得到的比对情况，使用samtools软件对all_bam文件中所有序列的比对情况进行过滤，对比对到参考基因组也就是宿主基因组的序列进行过滤，提取未比对上的所有序列，得到最终的去宿主的病原微生物的序列集合，即Clean reads。

[0020] 优选的，所述步骤(4)的具体方法为，使用soapdenovo软件对去宿主的reads进行组装，得到不依赖病原微生物参考基因组的contig序列。

[0021] 优选的，所述步骤(5)的具体方法为，使用bwa软件，以病原微生物的基因组作为参考基因组，将去宿主的reads比对到参考基因组上，获得去宿主的reads在病原微生物参考基因组上的比对情况，即pathogen_bam文件。

[0022] 优选的，所述步骤(6)的具体方法为，使用samtools软件统计步骤(5)中pathogen_bam文件每个位点的测序深度，将每个位置上占比最高的碱基确定为该位置的碱基序列，并依次进行延伸；对于比对中出现的indel，根据位点两端序列以及位置信息，保留reads中的indel，得到有参组装的contig序列。

[0023] 优选的，所述步骤(7)的具体方法为，使用自主编写脚本，以有参组装的contig序列为骨架，根据无参组装contig序列中的大片段易位或者倒位信息对骨架序列进行调整或重排，从而将无参组装的contig序列和有参组装的contig序列进行合并，得到合并后的contig序列。

[0024] 优选的，所述步骤(8)的具体方法为，使用bwa软件，以步骤(7)中得到的合并的contig序列作为参考基因组，将Clean reads比对到参考基因组上，得到最终的比对情况，即final_bam文件，并参照步骤(6)中的方法，对得到的final_bam文件进行二次有参组装，得到最终版的contig序列。

[0025] 与现有技术相比，本发明的有益之处在于：

[0026] 本发明将原始数据过滤后与宿主数据库进行比对，以达到在数据层面进行去除宿主序列的目的；

[0027] 本发明使用去宿主的序列分别进行无参组装和有参组装，并将无参和有参的结果进行合并，得到整合的基因组，然后将该整合基因组作为参考基因组，然后进行第二轮有参组装，通过去宿主后，采取无参组装和有参组装相结合的方法进行病原微生物的组装，组装的宏基因组没有宿主的污染，并且在点突变以外加入了结构变异的信息，准确度更高。

附图说明

[0028] 图1为本发明一种基于宏基因组测序数据组装病原微生物基因组的方法流程图。

具体实施方式

[0029] 下面结合附图、实施例、表格对本发明进行详细说明，附图仅用于示例性说明；应当理解，下面所提到的案例仅仅用来解释本发明，是为了便于描述本发明和简化描述，因此，不能理解为对本发明的限制。

[0030] 如图1所示，一种基于宏基因组测序数据组装病原微生物基因组的方法，具体步骤如下：

[0031] 步骤一.根据实验信息与数据的实际情况，使用fastp软件对含有甲型流感病毒的鸡胚尿囊液样本二代下机数据(Raw data)进行过滤，除去接头序列，低质量序列以及长度不符合要求的序列，得到符合质控要求的下机数据，即质控后高质量的二代数据。

[0032] 步骤二.质控后的下机数据中同时存在宿主的DNA序列以及病原微生物的DNA序列，使用bwa软件将质控后高质量的二代数据比对到宿主(鸡)基因组数据库，以宿主(鸡)的基因组作为参考基因组，将质控后的下机数据与参考基因组进行比对，获取下机数据中所有序列比对情况的all_bam文件。

[0033] 步骤三.根据步骤二中得到的比对情况，使用samtools软件对all_bam文件中所有序列的比对情况进行过滤，对比对到参考基因组(鸡基因组)的序列进行过滤，提取未比对上的所有序列，得到最终的去宿主的病原微生物的序列集合，即Clean reads。

[0034] 步骤四.使用soapdenovo软件对Clean reads进行组装，得到不依赖病原微生物参考基因组的重叠群(contig)序列。

[0035] 步骤五.使用bwa软件，以病原微生物的基因组作为参考基因组，将去宿主的Clean reads比对到参考基因组上，获得Clean reads在病原微生物参考基因组上的比对情况，即pathogen_bam文件。

[0036] 步骤六.使用samtools软件统计步骤五中pathogen_bam文件每个位点的测序深度，将每个位置上占比最高的碱基确定为该位置的碱基序列，并依次进行延伸；对于比对中出现的indel，根据位点两端序列以及位置信息，保留reads中的indel，得到有参组装的contig序列(即图1中的微生物基因组V1)。

[0037] 步骤七.使用自主编写脚本，以有参组装的contig序列为骨架，根据无参组装contig序列中的大片段易位或者倒位信息对骨架序列进行调整或重排，从而将无参组装的contig序列和有参组装的contig序列进行合并，得到合并后的contig序列(即图1中的微生物基因组V2)。

[0038] 步骤八.使用bwa软件，以步骤七中得到的合并的contig序列作为参考基因组，将Clean reads比对到参考基因组上，得到最终的比对情况(final_bam)，并参照步骤六中的方法，对得到的final_bam文件进行二次有参组装，得到最终版的contig序列(即图1中的微生物基因组V3)。

[0039] 步骤九.使用blast软件，将最终的contig序列与已有的病原微生物数据库进行比对，根据比对结果中序列之间的相似性与一致性等参数，对contig序列进行病原属种和亚型鉴定。

[0040] 为了设置生物学重复，避免实验中可能出现的错误并且减小误差，随机选取宿主鸡的三个样本进行甲型流感病毒数据统计，统计结果如表1所示。

[0041] 表1甲型流感病毒数据统计

[0042]

[0043]

[0044] 其中，HQ：高质量reads数量；

[0045] Host：宿主(鸡)的reads数量；

[0046] IAV：流感病毒的reads数量；

[0047] Others：除去宿主和流感病毒后剩余的其他reads数量；

[0048] Host％：宿主reads数占高质量reads数的百分比；

[0049] IAV％：流感病毒reads数占高质量reads数的百分比；

[0050] Others：其他reads数占高质量reads数的百分比。

[0051] 在对宿主鸡进行数据统计，将最终的contig序列与已有的病原微生物数据库进行比对，每一个样本分别选择8个目标片段，对这8个基因片段再进行分别的统计，统计结果如表2所示，表2是甲型流感病毒病毒基因组组装结果。

[0052] 表2甲型流感病毒病毒基因组组装结果

[0053]

[0054]

[0055] 其中，Gene：流感病毒基因节段；

[0056] Accession：参考序列GenBank号；

[0057] Length：组装序列长度；

[0058] Coverage：组装序列长度占参考序列长度的百分比；

[0059] Depth：组装序列中每个碱基的平均测序深度；

[0060] ≥10X：组装序列中测序深度≥10X的碱基百分比；

[0061] 根据表2的统计结果表明，采用本发明的方法，三个随机样本中，每个样本中的所有8个目标片段均被完整(100％)组装出来，并且目标片段的所有位点的覆盖度均在10X以上，准确度很高。

[0062] 本发明的组装方法适用于各类病原微生物，如细菌、真菌、病毒等的基因组组装。

[0063] 本发明采用的实施例只是较佳的实施例，用到的技术软件或物品如果没有明确说明，皆采用现有技术中常规使用标准的软件或物品，其他未述及的方法采用现有技术。

[0064] 本领域的普通技术人员将会意识到，这里所示的实施例是为了帮助读者理解本发明的原理，应理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其他各方面变换仍然在本发明的保护范围内。