技术领域
本发明涉及生物技术领域，特别涉及一种水体微生物定性与定量的检测方法。
背景技术
水体微生物是水体质量的指标，也是人畜病原菌的重要来源，因此，水体微生物精 确地定性与定量检测与监控是十分必要的。
现有水体微生物定性与定量检测技术包括形态学计数、芯片检测、16SrRNA测序、 宏基因组测序和实时定量PCR(PolymeraseChainReaction，聚合酶链式反应)。形态学计 数检测需要对微生物进行预培养，耗时长，不可培养微生物不可检测，一次仅能够检测一种 微生物，通量低，在计数时抽样量有限，且结果粗糙，无法对种以下的分类单元进行区分。芯 片检测所需的待测样品的DNA量大，需要对微生物进行预培养及富集处理，检测结果不准 确，且无法做定量检测。16SrRNA测序无法对种以下的分类单元进行区分。宏基因组测序深 度有限，对于低含量的微生物的定量检测准确度很差。实时定量PCR一次只能检测一种微生 物，通量低。另外，已有方法共有缺陷是，无法计算微生物定性与定量的可靠性，使得结论实 用性差。
发明内容
为了解决现有技术中微生物定性与定量检测不准确的问题，本发明实施例提供了 一种水体微生物定性与定量的检测方法。所述技术方案如下：
本发明实施例提供了一种水体微生物定性与定量的检测方法，所述方法包括：
确定待测样品中的目标微生物类群、目标微生物和非目标生物、以及不存在于所 述待测样品中的参考微生物，所述待测样品为水体；
根据所述目标微生物类群、所述目标微生物、所述参考微生物和所述非目标生物 的参考基因组序列，获得所述目标微生物类群的特征区域、所述目标微生物的特征区域和 所述参考微生物的特征区域；
制备扩增所述目标微生物类群的特征区域的第一多重扩增引物、扩增所述目标微 生物的特征区域的第二多重扩增引物和扩增所述参考微生物的特征区域的第三多重扩增 引物，将所述第一多重扩增引物、所述第二多重扩增引物和所述第三多重扩增引物混合得 到混合多重扩增引物；
向所述待测样品中加入所述参考微生物，获得混合样品；
提取所述混合样品的核酸；
利用所述混合多重扩增引物和所述混合样品的核酸进行扩增反应，获得扩增产 物；
利用所述扩增产物进行高通量测序，获得高通量测序片段；
根据所述高通量测序片段，对所述目标微生物类群和所述目标微生物进行定性和 定量分析。
具体地，所述目标微生物类群的数目≥1个，且每个所述目标微生物类群包括≥0 种所述目标微生物；
所述目标微生物为细菌、病毒、真菌、放线菌、立克次体、支原体、衣原体、螺旋体、 藻类和原生动物中的至少一种；
所述参考微生物为细菌、病毒、真菌、放线菌、立克次体、支原体、衣原体、螺旋体、 藻类和原生动物中的至少一种。
具体地，所述确定待测样品中的非目标生物的方法包括：将所述非目标生物确定 为除所述目标微生物类群之外的所有生物，若能获得所述目标微生物类群的特征区域，则 所述非目标生物指除所述目标微生物类群之外的所有生物；若不能获得所述目标微生物类 群的特征区域，则所述非目标生物指所述混合样品中，除所述目标微生物类群之外的其它 生物。
具体地，所述目标微生物类群的特征区域为所述目标微生物类群内的微生物的参 考基因组上的核酸序列；所述目标微生物类群的特征区域的两侧的序列在所述参考基因组 中为单一序列；所述目标微生物类群的特征区域的两侧的序列在所述目标微生物类群内不 同微生物间保守；所述目标微生物类群的特征区域的区分度≥3；
所述目标微生物的特征区域与所述目标微生物类群的特征区域同源；所述目标微 生物的特征区域的m2值≥2，其中，m2值为所述目标微生物的特征区域与所述目标微生物类 群内除所述目标微生物外的其它所述微生物间的差异碱基数的最小值；
所述参考微生物的特征区域为所述参考微生物的参考基因组上的核酸序列；所述 参考微生物的特征区域的两侧的序列在所述参考微生物的参考基因组中为单一序列；所述 参考微生物的特征区域的两侧的序列在除所述参考微生物外的其它生物中不具有同源性。
进一步地，所述区分度是指由同一所述混合多重扩增引物扩增的任一所述目标微 生物类群的特征区域与任一非特征区域间的差异碱基数的最小值，其中，所述非特征区域 是所述混合多重扩增引物以所述混合样品的核酸为模板的扩增产物，且所述非特征区域不 为所述目标微生物类群的特征区域，若无所述非特征区域，则所述区分度＝3×L1/4，其中， L1为所述目标微生物类群的特征区域的核酸序列长度。
具体地，在提取所述混合样品的核酸时，若所述待测样品中核酸的含量过低，则在 提取所述混合样品的核酸的过程中，加入所述混合多重扩增引物不能扩增的外源核酸。
具体地，所述目标微生物类群和所述目标微生物的定性分析方法如下：
将所述高通量测序片段与每种所述目标微生物类群的特征区域进行比对，当差异 碱基数≤n1时，则比对成功，相应的所述高通量测序片段为所述目标微生物类群的特征区 域，其中，n1为所述目标微生物类群的特征测序片段的最大容错碱基数；若比对成功的所述 目标微生物类群的特征区域≥1种时，则判断所述高通量测序片段为所述目标微生物类群 的特征测序片段；
将所述目标微生物的特征区域与每种同源的所述目标微生物类群的特征区域进 行比对，在所述目标微生物的特征区域中提取差异碱基组成所述目标微生物的标准基因 型；在所述目标微生物类群的特征测序片段上，提取所述目标微生物的标准基因型所对应 的碱基，组成所述目标微生物的测试基因型；若所述目标微生物的测试基因型与所述目标 微生物的标准基因型的差异碱基数≤n2，其中，n2为所述目标微生物的特征测序片段的最 大容错碱基数，则所述目标微生物的测试基因型所在的所述高通量测序片段为所述目标微 生物的特征测序片段；
将所述参考微生物作为仅包含一个所述目标微生物的所述目标微生物类群，计算 获得的所述目标微生物的特征测序片段，即为所述参考微生物的特征测序片段；
若所述目标微生物类群的特征测序片段存在的概率P5≥α5，则判断所述待测样品 中存在所述目标微生物类群，其中，α5为概率保障；若所述目标微生物类群的特征测序片段 存在的概率P5<α5，则判断所述待测样品中不存在所述目标微生物类群；
若所述目标微生物的特征测序片段存在的概率P6≥α6，则判断所述待测样品中存 在所述目标微生物，其中，α6为概率保障；若所述目标微生物的特征测序片段存在的概率P6 <α6，则判断所述待测样品中不存在所述目标微生物；
n1使得P1≤α1且P3≤α3，其中，P1为一条不是所述目标微生物类群的特征测序片 段的所述高通量测序片段被误判为所述目标微生物类群的特征测序片段而产生的假阳性 的概率；P3为一条所述目标微生物类群的特征测序片段被误判为不是所述目标微生物类群 的特征测序片段而产生的假阴性的概率；α1和α3为判断阈值；
n2使得P2≤α2且P4≤α4，其中，P2为一条不是所述目标微生物的特征测序片段的 所述高通量测序片段被误判为所述目标微生物的特征测序片段而产生的假阳性的概率；P4 为一条所述目标微生物的特征测序片段被误判为不是所述目标微生物的特征测序片段而 产生的假阴性的概率；α2和α4为判断阈值；
P5＝1-BINOM.DIST(S1,S1,P1,FALSE)，P6＝1-BINOM.DIST(S3,S3,P2,FALSE)，S1 为所有的所述目标微生物类群的特征区域的所述目标微生物类群的特征测序片段的数量 的中位数；S3为所有的所述目标微生物的特征区域的所述目标微生物的特征测序片段的数 量的中位数，FALSE为参数值；BINOM.DIST函数返回一元二项式分布的概率。
进一步地，所述目标微生物类群和所述目标微生物的定量分析方法如下：
所述目标微生物类群的量M1＝Mr×S1/S2，所述目标微生物类群的量的置信区间 为[M11，M12]，其中，Mr为加入所述待测样品中的所述参考微生物的量；S2为所有的所述参 考微生物的特征区域的所述参考微生物的特征测序片段的数量的中位数；M11和M12分别为 M1值的置信区间的下限与上限；
所述目标微生物的量M2＝M1×S3/S1，所述目标微生物的量的置信区间为[M21， M22]，M21和M22分别为M2值的置信区间的下限与上限；
M11＝M1×(1-S4/S1)，M12＝M1×(1+S5/S1)，M21＝M2×(1-S6/S3)，M22＝M2×(1+ S7/S3)；其中，S4为假阳性的所述目标微生物类群的特征测序片段的数量且S4＝CRITBINOM (nS,P1,α9)，其中，nS为计算S1的所述目标微生物类群的特征区域的所述多重扩增引物所 扩增的所述非特征区域的所述高通量测序片段的数量；S5为假阴性的所述目标微生物类群 的特征测序片段的数量且S5＝CRITBINOM(S1,P3,α9)，其中，α9为概率保障；S6为假阳性的 所述目标微生物的特征测序片段的数量且S6＝CRITBINOM(S1,P2,α10)，S7为假阴性的所述 目标微生物的特征测序片段的数量且S7＝CRITBINOM(S3,P4,α10)，其中，α10为概率保障； CRITBINOM函数返回使累积二项式分布大于等于临界值的最小值。
进一步地，P1＝BINOM.DIST(n1,m1,1-E,TRUE)，P2＝BINOM.DIST(n2,m2,1-E, TRUE)，P3＝1-BINOM.DIST(n1,L1,E,TRUE)，P4＝1-BINOM.DIST(n2,L2,E,TRUE)，其中，m1为 所述区分度；所述m2为所述目标微生物的特征区域与所述目标微生物类群的其它所述微生 物间差异碱基的最小值；L1为所述目标微生物类群的特征区域的长度；L2为所述目标微生 物的标准基因型的长度；E为碱基错误率。
本发明实施例提供的技术方案带来的有益效果是：本发明提供的方法不需要对微 生物进行预培养与增殖，耗时短，可同时检测多种微生物，通量高，计数时抽样量大，检测结 果精细，能够对分类单元进行区分，无需大量的DNA并避免了富集培养，检测结构无噪音且 准确，对于低含量微生物的定量准确度高，且对于微生物定性和定量的检测结果准确、分辨 率高、灵敏度高、有概率保障，检测过程简单、快速且流程规范。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面将对本发明实施方式作进一 步地详细描述。本发明中未标注说明的试剂均为常用市售试剂，在大多数生物技术公司均 可购买到且效果几乎无差别。
实施例、三角湖水体微生物的鉴定
水体可以包括江、河、湖、海、冰川、积雪、水库、池塘等中的水，也可以包括地下水 和大气中的水汽，本实施例中的待测样品为湖北武汉三角湖中的水体，检测水体中的微生 物是为了通过微生物的量监控水体受人畜粪便等污染的情况。
步骤一、确定待测样品中的目标微生物类群、目标微生物和非目标生物、以及不存 在于待测样品中的参考微生物，具体方法如下：
目标微生物类群的数目≥1个，且每个目标微生物类群包括≥0种目标微生物；目 标微生物可以为细菌、病毒、真菌、放线菌、立克次体、支原体、衣原体、螺旋体、藻类和原生 动物中的至少一种。本实施例的目的是鉴定待测样品中的肠道沙门氏菌，该菌的主要来源 为人畜粪便。三角湖紧邻江汉大学和多个生活小区，受粪便污染的可能性较大。肠道沙门氏 菌的拉丁学名为Salmonellaenterica，在NCBI(Nationalcenterforbiotechnology information，国家生物技术信息中心)上，已知参考基因组的肠道沙门氏菌生理小种共33 个(截止时间2015年6月2日)，具体见http://www.ncbi.nlm.nih.gov/genome/ genomegroups/152，这些生理小种构成本实施例的目标微生物类群。在这些生理小种中， Salmonellaentericasubsp.houtenaestr.ATCCBAA-1581致病性较强，作为本实施例的 目标微生物。
参考微生物可以为细菌、病毒、真菌、放线菌、立克次体、支原体、衣原体、螺旋体、 藻类和原生动物中的至少一种。参考微生物不存在于待测样品中。参考微生物的作用是为 待测样品中的目标微生物类群和目标微生物的定量提供一个参照。由于根癌农杆菌存在于 植物根中，所以不存在于待测样品中，因此，在本实施例中选取根癌农杆菌作为参考微生 物，其拉丁学名为AgrobacteriumtumefaciensK84。
具体地，确定待测样品中的非目标生物的方法包括：将非目标生物确定为除目标 微生物类群之外的所有生物，若能获得目标微生物类群的特征区域，则非目标生物指除目 标微生物类群之外的所有生物，其中，所有生物是指具有参考基因组的生物，是非目标生物 的最严格的标准。在本实施例中，将非目标生物确定为目标微生物类群之外的所有己知生 物时，可以找到目标微生物类群的特征区域(特征区域的获取过程见后，结果见表1)，因此， 本实施例中的非目标生物为除目标微生物类群之外的所有生物的集合。
将非目标生物确定为除目标微生物类群之外的所有生物，若不能获得目标微生物 类群的特征区域，则非目标生物指混合样品中，除目标微生物类群之外的其它生物，以缩小 非目标生物的范围，增加找到目标微生物类群的特征区域的可能性。在混合样品中，除目标 微生物类群之外的其它生物可根据经验确定，例如，本实施例中，混合样品包括水体和参考 微生物，则混合样品中不可能存在陆生植物成份以及专性寄生于陆生植物的微生物，因此， 若将本实施例中非目标生物确定为目标微生物类群之外的所有己知生物时，无法获得目标 微生物的特征区域，则非目标微生物可以确定为除目标微生物、陆生植物、专性寄生于陆生 植物的微生物之外的生物的集合。
步骤二、根据目标微生物类群的参考基因组序列、目标微生物的参考基因组序列、 参考微生物的参考基因组序列和非目标生物的参考基因组序列，获得目标微生物类群的特 征区域、目标微生物的特征区域和参考微生物的特征区域，具体方法如下：
目标微生物类群的特征区域为目标微生物类群内的微生物的参考基因组上的核 酸序列；目标微生物类群的特征区域的两侧的序列在参考基因组中为单一序列；目标微生 物类群的特征区域的两侧的序列在目标微生物类群内不同微生物间保守；目标微生物类群 的特征区域的区分度≥3。非特征区域不为目标微生物类群的特征区域，非特征区域是指混 合多重扩增引物以混合样品的核酸为模板的扩增产物；区分度是指由同一混合多重扩增引 物扩增的任一目标微生物类群的特征区域与任一非特征区域间的差异碱基数的最小值，若 无非特征区域，则区分度＝3×L1/4，其中，L1为目标微生物类群的特征区域的核酸序列长 度。
具体地，目标微生物类群的特征区域用于代表目标微生物类群，目标微生物类群 的特征区域存在，则代表目标微生物类群存在，目标微生物类群的特征区域的测序片段的 数量，代表目标微生物类群的数量。理想的目标微生物类群的特征区域的多重引物只扩增 目标微生物类群的特征区域，不能扩增非目标生物。这就要求目标微生物类群的特征区域 的两侧序列，即引物设计区域在非目标生物中不同源，那么，非目标生物不能被扩增，也不 能产生非特征区域。此时，特征区域与非特征区域间只能随机产生相同碱基，碱基共4种，相 同与不同的概率分别为1/4和3/4，因此，区分度为3×L1/4。目标微生物类群的特征区域的 区分度≥3是为了保证目标微生物类群的特征测序片段判断的假阳性率与假阴性率都较 低，其原理见表2。目标微生物类群的特征区域的两侧的序列在目标微生物类群内不同微生 物间保守，就可以用相同的引物扩增目标微生物类群内的不同微生物，以排除扩增效率对 目标微生物类群不同微生物间的相对定量的影响。
目标微生物的特征区域与目标微生物类群的特征区域同源；目标微生物的特征区 域的m2值≥2，其中，m2值为目标微生物的特征区域与目标微生物类群内除目标微生物外的 其它微生物间的差异碱基数的最小值。本实施例中的其它微生物是指除目标微生物外目标 微生物类群中的其它生理小种，m2值为目标微生物的特征区域分别与目标微生物类群中的 其它生理小种同源区域比较，所获得的差异碱基数目中的最小值。目标微生物定性与定量 分析时，重点是与目标微生物类群内其它微生物进行区分。目标微生物与目标微生物类群 亲缘关系往往较近，序列间相似性高，因此，难以区分。在目标微生物定性与定量分析时，只 关注了扩增子中与目标微生物类群内其它微生物间有差异的标准基因型，减少了误差的来 源，从而可以更好地将目标微生物从目标微生物类群内区分出来。当m2≥2时，判断测序片 段为目标微生物的特征测序片段的假阳性率与假阴性率都较低，因此，可以将目标微生物 从目标微生物类群中区分出来，其原理见表2。
参考微生物的特征区域为参考微生物的参考基因组上的核酸序列；参考微生物的 特征区域的两侧的序列在参考微生物的参考基因组中为单一序列；参考微生物的特征区域 的两侧的序列在除参考微生物外的其它生物中不具有同源性。
本实施例中，区分度是目标微生物类群的特征区域的唯一选择标准，根据检测的 目的的不同，也可以将具有特定基因序列的微生物作为目标微生物类群，并且将特定基因 序列作为目标微生物类群的特征区域。例如，可以将具有特定致病基因的微生物作为目标 微生物类群，并将该致病基因作为目标微生物的特征区域，以便根据致病基因的类型，指导 用药。同样，耐药性基因作为特定基因序列也可以指导用药。
步骤三、制备扩增目标微生物类群的特征区域的第一多重扩增引物、扩增目标微 生物的特征区域的第二多重扩增引物和扩增参考微生物的特征区域的第三多重扩增引物， 将第一多重扩增引物、第二多重扩增引物和第三多重扩增引物混合得到混合多重扩增引 物。
结合步骤二和步骤三具体方法如下：
在ftp://ftp.ncbi.nlm.nih.gov/genomes/中下载目标微生物类群中的不同生理 小种的基因组序列，采用软件Megablast(版本2.2.26)将它们的基因组与query序列(参考 序列)进行比对分析，本实施例中，query序列为NCBI上接收号为CM001471的基因组序列。 Megablast软件比对的各参数设置为：参数-e设置为1e-5；参数-p设置为0；参数–v设置为 5000；参数-m设置为1。比对完成后，获得在目标微生物类群的所有微生物间的同源序列，从 中挑选仅在query序列中出现1次的同源序列。以110bp为窗口大小，以10bp为步长，在挑选 的同源序列内作窗口平移。对于每一次平移获得的窗口，比较获得至少在目标微生物类群 内两种微生物间存在差异的碱基，截取该窗口中从第一个差异碱基起至最后一个差异碱基 止的区域作为特征区域，并统计该特征区域内的差异碱基的数量。向特征区域的两侧各延 伸长度为160bp-特征区域长度的区域作为引物搜索区，在引物搜索区内，搜索存在长度大 于20bp且在目标微生物类群内所有微生物间没有任何碱基差异的区域，作为特征区域的引 物设计区，放弃缺乏引物设计区的特征区域。
登录多重扩增引物在线设计网页https://ampliseq.com，在“Applicationtype” 选项选择“DNAHotspotdesigns(single-pool)”。若在本实施例中选择multi-pool，则多重 PCR将分多管进行，成本会有所增加。而选择single-pool的引物只需要一次多重PCR即可， 节省成本，但缺点是某些特征区域的引物设计可能失败，但由于基因组上的特征区域的数 目较多，少数特征区域引物设计失败不影响结果，所以，本实施例选择single-pool。将以上 获得的所有目标微生物类群的特征区域及其对应的引物设计区用100个碱基N(N代表A、T、C 和G四种碱基中的任意一种)连接起来，生成为一个引物设计的参考基因组。在“Selectthe genomeyouwishtouse”选项中选择“Custom”后，上传生成的引物设计的参考基因组。 DNAType选项选择“StandardDNA”，在AddHotspot选项中，填写特征区域在生成的引物设 计的参考基因组中的起始和终止位置。最后点击“Submittargets”按钮提交并获得目标微 生物类群的特征区域的多重扩增引物序列。
利用设计的多重扩增引物对目标微生物类群利用BLASTN(BasicLocal AlignmentSearchTool，基本局部比对搜索工具)(version2.2.26)做比对分析，从中挑 选正反向引物中至少有一个具有特异性的引物。将挑选出来的引物再与非目标生物的基因 组做BLASTN比对分析，检查它们是否可以扩增非目标生物的基因组。本实施例中，非目标生 物为除目标微生物类群外的所有生物，非目标生物的基因组为NCBI的NT/NR库。判断引物可 以扩增的标准为：扩增区长度不超过200bp，引物匹配长度大于15bp且引物3’端的5个碱基 以内没有碱基缺失或错配。若引物不能扩增任何非目标生物，此时，引物所对应的目标微生 物的特征区域的区分度m1＝3×L1/4，若引物可扩增部分非目标生物，则将将该引物扩增的 任一非目标生物的扩增产物与任一目标微生物类群的特征区域进行比对，获得所有比对 中，差异碱基数目的最小值为区分度m1，保留m1≥3的目标微生物类群的特征区域，进一步 去掉含有简单重复序列或在基因组上为多拷贝的特征区域。从保留的目标微生物类群的特 征区域中，进一步优选目标微生物类群的特征区域并选择目标微生物的特征区域。
进一步地，目标微生物类群的特征区域的优选方法如下：将特征区域与非目标生 物的参考基因组做BLASTN比对，去掉与非目标生物存在95％以上同源性的特征区域，将剩 余的特征区域在目标微生物与目标微生物类群内其它微生物间利用软件muscle(版本： V3.6)按其默认参数进行比对，获得差异碱基数的最小值，即m2值。保留m2≥2的目标微生物 类群的特征区域，从保留的特征区域中，任意挑选区分度m1与m2均较大的2个及2个以上的 特征区域同时作为目标微生物类群的特征区域和目标微生物的特征区域，其对应的多重扩 增引物同时作为第一多重扩增引物和第二多重扩增引物。
按与寻找目标微生物类群的特征区域相类似的方法，获取参考微生物特征区域及 其对应的第三多重扩增引物，下面重点描述不同之处，相同之处不再重复描述。同样采用软 件Megablast(版本2.2.26)对参考微生物基因组与query序列(参考序列)进行比对分析， query序列为AgrobacteriumtumefaciensK84的基因组序列。比对完成后，获得参考微生 物基因组中，仅在query序列中出现1次的单一序列。将单一序列与NCBI的NT/NR库比对，放 弃在非目标生物中存在同源序列的单一序列，从单一序列中随机挑选不重叠的110bp长度 的作为特征区域，其两侧的序列作为引物设计区域。于多重扩增引物在线设计网页 https://ampliseq.com上设计特征区域的多重扩增引物，进一步筛选成功设计了多重扩增 引物的特征区域，具体方法如下：去掉含有简单重复序列或在基因组上为多拷贝的特征区 域，将剩下的特征区域与非目标生物的参考基因组做BLASTN比对，去掉与非目标生物存在 95％以上同源性的特征区域。从保留下来的特征区域中，随机挑选2个及2个以上特征区域 作为参考微生物类群的特征区域，其对应的多重扩增引物作为第三多重扩增引物。
由生工生物工程(上海)股份有限公司逐一合成以上获得的第一多重扩增引物、第 二多重扩增引物和第三多重扩增引物中的每一重扩增引物、以及由每个多重扩增引物对应 的扩增的模板序列，模板序列是指每个多重扩增引物在填入AddHotspot选项的扩增区域。 按照美国赛默飞世尔公司的StepOne实时定量PCR仪的操作手册(PartNumber 4376784Rev.E)检测每个多重扩增引物的扩增效率，仅保留扩增效率在95％～105％的多重 扩增引物，以减少扩增效率的差异对微生物定性与定量的影响。由于扩增效率影响较少，因 此，目标微生物类群与目标微生物的特征区域也可以不同，以方便更容易分别找到各自的 特征区域。将以上获得的第一多重扩增引物、第二多重扩增引物和第三多重扩增引物保留 下来的多重扩增引物按多重扩增引物在线设计网页https://ampliseq.com上的合并程序 进行合并，获得混合多重扩增引物，混合多重扩增引物由美国赛默飞世尔公司合成后，以液 体形式提供。本实施例最终获得的特征区域相关信息见表1。表1中的起始位置与终止位置 是指特征区域在query序列上的参考基因组上的起始和终止位置。
表1本实施例提供的引物相关信息
步骤四、向待测样品中加入参考微生物，获得混合样品，具体方法如下：
参考微生物不存在于待测样品中，所以，可以将参考微生物作为内部参照，并与待 测样品中的微生物进行平行操作，对待测样品中的目标微生物类群与目标微生物进行定 量。参考微生物的加入量控制为大约可以提取10ng的混合样品的核酸(DNA)，以正常构建高 通量测序文库，同时，参考微生物的加入量又不至于使得参考微生物所占的比例过大，占用 过多的高通量测序数据量。本实施例混合样品的获取方法如下：从三角湖中取水样10公斤， 选4000rpm离心20分钟后，再用12000rpm离心10分钟，收集菌体深沉。将浓度为2OD(OD为菌 液最大吸光度值)的参考微生物的菌液0.2mL置于1.5mL的离心管中真空冷冻离心干燥后， 加入菌体深沉中，混匀，即得到待测样品与参考微生物的混合样品。通过血球板计数，计算 获得加入混合样品的参考微生物的量见表2。
步骤五、提取混合样品的核酸，具体方法如下：
在提取所述混合样品的核酸时，若待测样品中核酸的含量过低(低于1ug)，将影响 混合样品的核酸的提取，则可以在混合样品的核酸的提取过程中，加入混合多重扩增引物 不能扩增的外源核酸。所加入的外源核酸不存在与自然界中，因而不干扰微生物检测。外部 RNA对照协会设计了并验证了一套核酸序列，它们在自然界中不存在，可以作为本发明实施 例中的外源核酸，其序列可参考https://tools.lifetechnologies.com/content/sfs/ manuals/cms_095047.txt。外源核酸的加入量为1ug左右，该加入量可以保证混合样品的核 酸能够正常提取。在本实施例中，待测样品为水体，其核酸含量较低，因此，向混合样品中加 入外源核酸，即1ug外部RNA对照协会设计的ERCC-00014基因。利用细菌基因组DNA提取试剂 盒(货号：DP302，生产公司：天根生化科技(北京)有限公司)按其操作手册提供的方法提取 获得的混合样品的核酸。
步骤六、利用混合多重扩增引物和混合样品的核酸进行扩增反应，获得扩增产物， 具体方法如下：
利用文库构建试剂盒2.0(由美国LifeTechnology公司生产，货号为4475345)多重PCR扩增混合样品的核酸后，利用扩增产物构建高通量测序文库。该试剂盒包括以下试剂：5×IonAmpliSeqTMHiFiMix、FuPa试剂、转换试剂、测序接头溶液和DNA连接酶。文库构建的方法按该试剂盒的操作手册《IonAmpliSeqTMLibraryPreparation》(出版号：MAN0006735，版本：A.0)进行。多重PCR的扩增体系如下：5×IonAmpliSeqTMHiFiMix4μl、合成的混合多重扩增引物4μl、提取的混合样品的核酸10ng和无酶水11μl。多重PCR的扩增程序如下：99℃，2分钟；(99℃，15秒；60℃，4分钟)×25个循环；10℃保温。利用FuPa试剂消化掉多重PCR扩增产物中多余的引物后，再进行磷酸化，具体方法为：向多重PCR的扩增产物中加入2μLFuPa试剂，混匀后，在PCR仪上按如下程序反应：50℃，10分钟；55℃，10分钟；60℃，10分钟；10℃保存，得到混合物a，混合物a为含有经过磷酸化的扩增产物溶液。将磷酸化的扩增产物连接上测序接头，具体方法为：向混合物a中加入转换试剂4μL、测序接头溶液2μL和DNA连接酶2μL，混匀后，在PCR仪上按如下程序反应：22℃，30分钟；72℃，10分钟；10℃保存，得到混合液b。利用标准的乙醇沉淀方法纯化混合液b后溶解于10μL无酶水中。利用美国Invitrigen公司生产的dsDNAHSAssayKit(货号为Q32852)并按照其说明书进行测定，获得混合液b的质量浓度后，将纯化后混合液b稀释至15ng/ml，得到浓度约100pM的高通量测序文库。
步骤七、利用扩增产物进行高通量测序，获得高通量测序片段，具体方法如下：
利用获得的高通量测序文库和试剂盒IonPITemplateOT2200Kitv2(美国 invirtrigen公司生产，货号为4485146)进行测序前的ePCR(EmulsionPCR，乳化聚合酶链 反应)扩增，操作方法按该试剂盒的操作手册进行。利用ePCR产物和试剂盒IonPI Sequencing200Kitv2(美国invirtrigen公司生产，货号为4485149)在Proton二代高通 量测序仪上进行高通量测序，操作方法按该试剂盒的操作手册进行。在本实施例中，高通量 测序量设置为1M测序片段(1M＝100万)。
根据测序片段的引物，将高通量测序片段比对到对应的目标微生物类群的特征区 域、目标微生物的特征区域和参考微生物的特征区域。去掉比对不成功和特征区域不完整 的测序片段，比对不成功的测序片段多为非特异扩增产物，特征区域不完整的测序片段是 指没能将表1中的特征区域的起始位置到终止位置的序列检测完整。
步骤八、根据高通量测序片段，对目标微生物类群和目标微生物进行定性和定量 分析，具体方法如下：
本发明提供的微生物定性定量分析的基本原理是：特征区域代表了目标微生物类 群和目标微生物，若存在特征区域的测序片段，表明目标微生物类群或目标微生物存在，而 特征区域的测序片段的数量也代表了目标微生物类群和目标微生物的数量。与其它微生物 定性与定量检测不同，本发明实施例计算了微生物定性与定量的可靠性，同时，增强了结论 的实用性。本发明实施例需要先理清参数间复杂的关系，才能实现任意微生物的定性、定量 检测，并获得可靠的结论，本发明的具体参数及其推算原理见表2。表2中单元格、符号与公 式的定义与Excel2010相同，其中，单元格“基本参数”为A1，其它单元格参照A1按Excel 2010的规则进行定义。
定性分析方法如下：将高通量测序片段与每种目标微生物类群的特征区域进行比 对，当差异碱基数≤n1时，则比对成功，相应的高通量测序片段为目标微生物类群的特征区 域，其中，n1为目标微生物类群的特征测序片段的最大容错碱基数；若比对成功的目标微生 物类群的特征区域≥1种时，则判断高通量测序片段为目标微生物类群的特征测序片段。
将目标微生物的特征区域与每种同源的目标微生物类群的特征区域进行比对，在 目标微生物的特征区域中提取差异碱基组成目标微生物的标准基因型，这里的差异碱基是 指目标微生物的特征区域与任何一个目标微生物类群内的微生物比较，存在差异的碱基的 总和。在目标微生物类群的特征测序片段上，提取目标微生物的标准基因型所对应的碱基， 组成目标微生物的测试基因型；若目标微生物的测试基因型与目标微生物的标准基因型的 差异碱基数≤n2，其中，n2为目标微生物的特征测序片段的最大容错碱基数，则目标微生物 的测试基因型所在的高通量测序片段为目标微生物的特征测序片段。特别地，当目标微生 物类群中仅包含了一个目标微生物时，此时的标准基因型和测试基因型的碱基数均为0个， 因此，它们之间的差异碱基数也为0个，则不论n2为多大，均将目标微生物的测试基因型所 在的高通量测序片段判定为目标微生物的特征测序片段。按以上方法，分别获得了目标微 生物类群和目标微生物的特征区域的特征片段数，其结果列于表1。在本实施例中，n1和n2 的值见表2，其推算过程见后。
n1使得P1≤α1且P3≤α3，其中，P1为一条不是目标微生物类群的特征测序片段的 高通量测序片段被误判为目标微生物类群的特征测序片段而产生的假阳性的概率；P3为一 条目标微生物类群的特征测序片段被误判为不是目标微生物类群的特征测序片段而产生 的假阴性的概率；α1和α3为判断阈值。
n2使得P2≤α2且P4≤α4，其中，P2为一条不是目标微生物的特征测序片段被误判 为目标微生物的特征测序片段而产生的假阳性的概率；P4为一条目标微生物的特征测序片 段被误判为不是目标微生物的特征测序片段而产生的假阴性的概率；α2和α4为判断阈值； 本发明实施例中的各种阈值的大小由现实需要确定，例如，某些病菌危害性极大，漏检(假 阴性)将引起严重的后果，那么，就要控制假阴性，α2和α4值要低。若无特殊要求，则采用较 低假阳性与假阴性为原则，本实施例子属于后者，α1和α3取值为0.01％，即大约1万条特征 序列出现1条假阳性或假阴性，其准确性是很高的，之所以可以控制如此高的准确性，是因 为特征序列中的m1值较大，很容易与其它非目标生物区分开，从而将假阳性率与假阴性率 都控制在一个很低的水平。α2和α4的取值为0.5％，即大约1千条特征序列出现5条假阳性或 假阴性，可见其准确性很高。P1＝BINOM.DIST(n1,m1,1-E,TRUE)，P2＝BINOM.DIST(n2,m2, 1-E,TRUE)，P3＝1-BINOM.DIST(n1,L1,E,TRUE)，P4＝1-BINOM.DIST(n2,L2,E,TRUE)，其中， m1为区分度，具体指用于计算S1的目标微生物类群的特征区域对应的区分度，本实施例中， m1的值见表1和表2；m2为目标微生物的特征区域与目标微生物类群的其它微生物间差异碱 基的最小值，具体指用于计算S3的目标微生物对应的特征区域的m2的值，本实施例中，m2的 值见表1和表2；L1为目标微生物类群的特征区域的长度，本实施例子中，L1的值见表2；L2为 目标微生物的标准基因型的长度，本实施例子中，L2的值见表2；E为碱基错误率，其由测序 错误率E1和自然突变率E2组成，本实施例中，PROTON高通量测序仪的测序错误率E1≤1％， 根据我们的调查，微生物小种(如P1-P6白叶枯小种)的参考基因组之间的变异率一般小于 0.5％，而自然突变率是低于小种间的变异率的，因此，自然突变率E2≤0.5％，为了本发明 的方法适应性更广，取E2≤1％，则本实施例中，E≤2％，为了使得本实施例中微生物的定性 与定量的结论正确率的概率更可靠，取E值的最大值2％进行计算。将以上参数值代入P1和 P3的公式中后，将n1的值从0开始逐渐增加，计算得P1和P3的值，当n1＝13时，计算得P1≤α1 且P3≤α3，因此，本实施例中，n1＝13(见表2)，n1＝13对应的P1和P3的值为本实施例中P1和 P3的值。按类似的方法，将以上参数值代入P2和P4的公式中后，将n2的值从0开始逐渐增加， 计算得P2和P4的值，当n2＝2时，P2≤α2，P4≤α4，因此，本实施例中，n2＝2(见表2)，n2＝2对 应的P2和P4的值为本实施例中P2和P4的值。
将参考微生物作为仅包含一个目标微生物的目标微生物类群，计算获得的目标微 生物的特征测序片段，即为参考微生物的特征测序片段。参考微生物的特征区域的特征片 段数见表1和表2。
若目标微生物类群的特征测序片段存在的概率P5≥α5，则判断待测样品中存在目 标微生物类群；若目标微生物类群的特征测序片段存在的概率P5<α5，则判断待测样品中不 存在目标微生物类群，其中，α5为概率保障，本实施例中，α5取值为99.99％。P5＝1- BINOM.DIST(S1,S1,P1,FALSE)，S1为所有的目标微生物类群的特征区域的目标微生物类群 的特征测序片段的数量的中位数，在本实施例中，目标微生物类群的第2个特征测序片段的 数量为所有目标微生物类群的特征测序片段的数量的中位数，所以本实施例S1的值见表1 和表2，将本实施例中S1的值和P1的值代入P5的计算公式计算获得P5≥α5，因此，判断本实 施例中，待测样品中存在目标微生物类群，FALSE为参数值，BINOM.DIST函数返回一元二项 式分布的概率。
若目标微生物的特征测序片段存在的概率P6≥α6，则判断待测样品中存在目标微 生物；若目标微生物的特征测序片段存在的概率P6<α6，则判断待测样品中不存在目标微生 物；α6为概率保障。本实施例中，α6取值为99.99％。P6＝1-BINOM.DIST(S3,S3,P2,FALSE)， BINOM.DIST函数返回一元二项式分布的概率，S3为所有的目标微生物的特征区域的目标微 生物的特征测序片段的数量的中位数，在本实施例中，目标微生物的第2个特征测序片段的 数量为所有目标微生物的特征测序片段的数量的中位数，其对应的S3的值见表1和表2，将 本实施例中S3的值和P2的值代入P6的计算公式计算获得P6≥α6，因此，判断本实施例中，待 测样品中存在目标微生物。
此外，α5和α6均是人们根据实际需要定的，α5和α6均可以相同也可以不同，其区别 取决于实际需要，当要严格控制某种微生物时，α5和α6的取值均较大，反之，α5和α6的取值 均较小。此外，本发明实施例中所有的a值的取值均遵循该原理。
定量分析方法如下：目标微生物类群的量M1＝Mr×S1/S2，其中，Mr为加入待测样 品中的参考微生物的量，本实施例中，Mr的值见表2；S2为所有的参考微生物的特征区域的 参考微生物的特征测序片段的数量的中位数，在本实施例中，参考微生物的第2个特征测序 片段的数量为所有参考微生物的特征测序片段的数量的中位数，其对应的S2的值见表1和 表2；将以上参数和通过定性分析获得的S1的值代入M1的计算公式中，计算获得M1值，即待 测样品中，目标微生物类群中的微生物的量为M1＝126522个。
目标微生物类群的量的置信区间为[M11，M12]，M11和M12分别为M1值的置信区间 的下限与上限。M11＝M1×(1-S4/S1)，M12＝M1×(1+S5/S1)，其中，S4为假阳性的目标微生 物类群的特征测序片段的数量且S4＝CRITBINOM(nS,P1,α9)，S5为假阴性的目标微生物类 群的特征测序片段的数量且S5＝CRITBINOM(S1,P3,α9)，其中，α9为概率保障，本实施例中， α9取值为99.50％，CRITBINOM函数返回使累积二项式分布大于等于临界值的最小值，nS为 计算S1的目标微生物类群的特征区域的多重扩增引物所扩增的非特征区域的高通量测序 片段的数量，即是指多重引物所扩增的除目标微生物的特征测序片段之外的其它高通量测 序片段。在本实施例中，nS为目标微生物类群中第2个特征区域的多重扩增引物扩增所产生 的非特征区域的高通量测序片段的数量，本实施例中，nS的值见表2。将nS的值和P1的值代 入S4的公式计算获得S4的值，将本实施例S1的值和P3的值代入S5的公式计算获得S5的值。 获得M11和M12公式中所有参数的值后，计算获得本实施例中M11和M12的值，进而获得M1的 置信区间，即目标微生物类群的量的置信区间为[126522，126547]。
目标微生物的量M2＝M1×S3/S1，将M1、S3和S1的值代入上述公式，获得目标微生 物的量M2＝16192。
目标微生物的量的置信区间为[M21，M22]，M21和M22分别为M2值的置信区间的下 限与上限；M21＝M2×(1-S6/S3)，M22＝M2×(1+S7/S3)；其中，S6为假阳性的目标微生物的 特征测序片段的数量且S6＝CRITBINOM(S1,P2,α10)，S7为假阴性的目标微生物的特征测序 片段的数量且S7＝CRITBINOM(S3,P4,α10)，其中，α10为概率保障；CRITBINOM函数返回使累 积二项式分布大于等于临界值的最小值。本实施例中，α10取值为99.50％，将本实施例S1和 S3的值，以及P2和P4的值代入S6和S7的计算公式，计算得到S6和S7的值。进一步将S6、S7、M1 和S3的值代入M21和M22的计算公式，计算得到M21和M22的值，进而得到目标微生物的量的 置信区间为[15809，16213]。
表2为本实施例微生物定性与定量分析参数及其推算原理
本发明改变了已有方法中一次只能检测少数几种微生物、只能将微生物区分到 种、定量不准、检测结果无概率保障、需要预培养、检测周期长、某些微生物不可培养因而不 可检测、微生物可培养性不同而导致的定量失真、定量粗糙等诸多问题，为微生物检测提供 了一种全面、准确、快速、精细的定性与定量检测新方法。