一种利用细菌OTUs丰度预测河流水质理化指标的方法

一种利用细菌otus丰度预测河流水质理化指标的方法
技术领域
1.本发明涉及污水处理领域，具体涉及一种利用细菌otus丰度预测河流水质理化指标的方法。


背景技术：

2.河流水质评价是水环境管理的基础性工作，对污水处理工作的开展具有重要的指导意义。河流水质评价一般基于对多种河流水质指标进行检测和分析，主要分为物理水质指标、化学水质指标和生物水质指标，其中物理和化学指标(简称理化指标)是反映河流水质和需要常规监测的重要指标，例如传统的河流水质评价主要以碳、氮、磷、叶绿素等相关的理化指标为基础，通过对其进行定量分析以对水质物质导致的污染状态进行判断，并且通常需要多个理化指标的检测数据才能对污染状况进行综合评价。理化水质指标的检测通常费时费力，而且不同指标需要分别采用特定的检测方法来获得，针对不同的河流水样更需要一一进行检测，当河流水样的样本量较大时，水质评价工作面临挑战。
3.微生物能对外界营养物质及环境因素作出迅速而灵敏的反应，因此，长期以来，研究学者们一直致力于探索利用微生物指标评价水质状况的可行性，将水体污染所造成的微生物群落结构组成及功能的变化进行量化，并试图建立其与外界环境变化间的相关关系，从而利用微生物有效地评价和监测水污染状况。例如，有研究发现药物和表面活性剂对医院废水与城市污水的抗性组和微生物群具有显著影响；还有研究利用细菌16s rrna测序数据，对26个地球化学特征进行有效预测。然而，目前尚无采用微生物指标对河流水质理化指标进行准确预测的方法。


技术实现要素：

4.因此，本发明旨在提供一种利用细菌otus丰度预测河流水质理化指标的方法。
5.为实现上述目的，本发明提供如下技术方案：
6.本发明提供一种利用细菌otus丰度预测河流水质理化指标的方法，包括以下步骤：
7.(1)取若干河流水样进行16s rrna测序，确定细菌otus相对丰度，并对所述河流水样的水质理化指标进行检测；
8.(2)以细菌otus相对丰度作为输入变量，以水质理化指标作为输出变量，采用随机森林算法建立初始回归模型，并对输入变量进行重要度排序，从中选出重要细菌变量；
9.(3)以所述重要细菌变量作为输入变量，以水质理化指标作为输出变量，采用随机森林算法建立预测模型；
10.(4)利用所述预测模型对待测河流水样的水质理化指标进行预测。
11.进一步地，所述河流水质理化指标包括以下指标中的至少一种：总氮(tn)、硝态氮(no
3-‑
n)、氨氮(nh
4+-n)、总磷(tp)、可溶解性磷(srp)、叶绿素a(chl_a)、溶解性有机碳(doc)、溶解氧(do)、氧化还原电位(orp)、电导率(cond)、水温、tn/tp、no
3-‑
n/tp、nh
4+-n/
tn、no
3-‑
n/tn、no
3-‑
n/nh
4+-n、doc/tn、卡马西平(cbz)、红霉素(ery)、9-芴酮(9-fl)。
12.进一步地，步骤(3)中，分别以重要度排序前i名的变量作为重要细菌变量来建立预测模型，i＝10～50，比较不同预测模型之间的均方误差mse和/或拟合系数r2，选择mse值最小或者r2值最大的预测模型作为最终的预测模型，用来对待测河流水样的水质理化指标进行预测。
13.进一步地，步骤(2)中，采用％incmse评价输入变量的重要度；通过r软件中的“importance”和“varimpplot”获得输入变量的％incmse排序。
14.进一步地，通过r软件中的“randomforest”包建立所述初始回归模型和预测模型。
15.进一步地，步骤(4)中，利用所述预测模型对待测河流水样的水质理化指标进行预测的方法包括：
16.取所述待测河流水样进行16s rrna测序，确定细菌otus相对丰度，选取其中与所述预测模型中的重要细菌变量相应的变量作为输入变量，代入所述预测模型，得到的输出变量值即为所述待测河流水样的水质理化指标预测结果。
17.进一步地，所述的利用细菌otus丰度预测河流水质理化指标的方法，还包括：将所述河流水样的数据集随机划分为训练集和测试集，利用所述训练集建立所述初始回归模型和预测模型，利用所述训练集建立所述预测模型，优选的，以所述数据集的70％作为训练集，以所述数据集的30％作为测试集。
18.进一步地，以拟合系数r2衡量所述预测模型的预测能力：
19.当r2≤0.3时，预测值与观测值拟合效果差，所述预测模型的预测能力差；
20.当0.3＜r2≤0.4时，预测值与观测值拟合效果弱，所述预测模型的预测能力弱；
21.当0.4＜r2≤0.6时，预测值与观测值拟合效果中等，所述预测模型的预测能力中等；
22.当0.6＜r2≤1.0时，预测值与观测值拟合效果强，所述预测模型的预测能力强。
23.进一步地，步骤(1)中，所述确定细菌otus相对丰度的方法包括：
24.对所述河流水样进行dna提取，扩增16s rrna基因v4区；
25.对pcr产物进行纯化，并测量纯化pcr产物的dna浓度；
26.对纯化pcr产物进行基因测序，原始序列进行除杂，得到的高质量双端序列基于重复区域被连接成tags，过滤嵌合子，并将序列聚类成otus，计算细菌otus相对丰度。
27.进一步地，采用dna分离试剂盒进行dna提取。
28.进一步地，用条形码引物515f和806r扩增16s rrna基因v4区。
29.进一步地，使用dna凝胶提取试剂盒纯化pcr产物。
30.进一步地，使用荧光计测量纯化pcr产物的dna浓度。
31.进一步地，使用illumina hiseq 2500测序平台对纯化pcr产物进行基因测序。
32.进一步地，原始序列使用flash和trimmomatic软件进行去杂，去除引物错配数在2以上、拼接序列的重叠区错配率在0.2以上的序列。
33.进一步地，用uparse软件过滤嵌合子，并对序列进行聚类，根据97％的序列相似度将所有序列进行同源比对并聚类成otus。
34.本发明技术方案，具有如下优点：
35.本发明提供的利用细菌操作分类单元(operational taxonomic units，otus)丰
度预测河流水质理化指标的方法，通过预先取河流水样进行机器学习来建立预测模型，从而实现对待测河流水样理化指标的预测。具体是以细菌otus相对丰度作为输入变量，以水质理化指标作为输出变量，采用随机森林(random forest，rf)算法建立初始回归模型，并对输入变量进行重要度排序，从中选出重要细菌变量后再次采用随机森林算法建立预测模型。本发明提供的方法通过输入变量重要度排序优化输入集，使预测模型的预测结果更加精准；仅需要待测河流水样的16s rrna的测序数据即可对河流水质理化指标进行精准预测，尤其适用于大批量水样的预测，能够节约大量的时间、人力和检测试剂，具有重要的应用价值和现实意义；同时通过解释变量的重要度排序还可以识别出对其预测的水质理化指标敏感的细菌otus。
附图说明
36.为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
37.图1是本发明实施例中采用全部及重要度排序前10、20、30、40、50的细菌otus相对丰度数据对tp预测的模型训练效果(基于训练集)，其中，(a)为tp预测的均方误差mse，(b)为tp预测的决定系数r2；
38.图2是本发明实施例中采用全部及重要度排序前10、20、30、40、50的细菌otus相对丰度数据对tn/tp预测的模型训练效果(基于训练集)，其中，(a)为tn/tp预测的均方误差mse，(b)为tn/tp预测的决定系数r2；
39.图3是本发明实施例中基于全部及重要度排序前10、20、30、40、50的细菌otus相对丰度数据建立的tp预测模型的泛化性能(基于测试集)，其中，(a)为tp预测的均方误差mse，(b)为tp预测的决定系数r2；
40.图4是本发明实施例中基于全部及重要度排序前10、20、30、40、50的细菌otus相对丰度数据建立的tp预测模型的泛化性能(基于测试集)，其中，(a)为tn/tp预测的均方误差mse，(b)为tn/tp预测的决定系数r2；
41.图5是本发明实施例中基于重要度排序前20的细菌otus的最优预测效果，其中，(a)为tp的最优预测效果，(b)为tn/tp的最优预测效果；
42.图6是本发明实施例中tn和no
3-‑
n的最优预测效果，其中，(a)为基于重要度排序前20的细菌otus的tn最优预测效果，(b)为基于重要度排序前10的细菌otus的no
3-‑
n最优预测效果；
43.图7是本发明实施例中nh
4+-n和srp的最优预测效果，其中，(a)为基于重要度排序前10的细菌otus的nh
4+-n最优预测效果，(b)为基于重要度排序前20的细菌otus的srp最优预测效果；
44.图8是本发明实施例中doc和chl_a的最优预测效果，其中，(a)为基于重要度排序前30的细菌otus的doc最优预测效果，(b)为基于重要度排序前20的细菌otus的chl_a最优预测效果；
45.图9是本发明实施例中do和orp的最优预测效果，其中，(a)为基于重要度排序前20
的细菌otus的do最优预测效果，(b)为基于重要度排序前20的细菌otus的orp最优预测效果；
46.图10是本发明实施例中tds和cond的最优预测效果，其中，(a)为基于重要度排序前10的细菌otus的tds最优预测效果，(b)为基于重要度排序前20的细菌otus的cond最优预测效果；
47.图11是本发明实施例中ph和temp(水温)的最优预测效果，其中，(a)为基于重要度排序前30的细菌otus的ph最优预测效果，(b)为基于重要度排序前20的细菌otus的temp最优预测效果；
48.图12是本发明实施例中no
3-‑
n/tp和no
3-‑
n/tn的最优预测效果，其中，(a)为基于重要度排序前20的细菌otus的no
3-‑
n/tp最优预测效果，(b)为基于重要度排序前20的细菌otus的no
3-‑
n/tn最优预测效果；
49.图13是本发明实施例中nh
4+-n/tn和no
3-‑
n/nh
4+-n的最优预测效果，其中，(a)为基于重要度排序前30的细菌otus的nh
4+-n/tn最优预测效果，(b)为基于重要度排序前20的细菌otus的no
3-‑
n/nh
4+-n最优预测效果；
50.图14是本发明实施例中ery和cbz的最优预测效果，其中，(a)为基于重要度排序前20的细菌otus的ery最优预测效果，(b)为基于重要度排序前20的细菌otus的cbz最优预测效果。
51.图15是本发明实施例中doc/tn和9-fl的最优预测效果，其中，(a)为基于重要度排序前20的细菌otus的doc/tn最优预测效果，(b)为基于重要度排序前20的细菌otus的9-fl最优预测效果。
具体实施方式
52.提供下述实施例是为了更好地进一步理解本发明，并不局限于所述最佳实施方式，不对本发明的内容和保护范围构成限制，任何人在本发明的启示下或是将本发明与其他现有技术的特征进行组合而得出的任何与本发明相同或相近似的产品，均落在本发明的保护范围之内。
53.实施例中未注明具体实验步骤或条件者，按照本领域内的文献所描述的常规实验步骤的操作或条件即可进行。所用原料或仪器，均为可以通过市购获得的常规产品，包括但不限于本技术实施例中采用的原料或仪器。
54.实施例1
55.本实施例提供一种利用细菌otus丰度预测河流水质理化指标tp和tn/tp的方法，详述步骤如下：
56.(1)水样采集
57.选取中国海河水系中的潮白河(北纬39
°
至40
°
，东经116
°
至117
°
)流域为研究对象。于2016年12月(冬季)，2017年3月(春季)，2017年6月(夏季)，2017年9月(秋季)，在潮白河流域进行采样。一共34个采样点，采样点均位于潮白河的主河道及重要支流汇入口，分别为13个山区、7个城区及14个农业区样点。总共采集到126个水样。
58.(2)水样理化指标检测
59.对采集水样的理化指标的分析方法，主要参考中华人民共和国国家标准及《水和废水监测分析方法(第三版)》，总氮(tn)和总磷(tp)的检测方法、参考标准及主要仪器设备
如表1所示。
60.表1水质检测方法、参考标准及主要仪器设备
[0061][0062][0063]
注：各指标检测结果取三次检测结果的平均值。
[0064]
(3)水样微生物测序
[0065]
使用dna分离试剂盒(mobio laboratories inc.，carlsbad，ca，usa)对过滤后的水样进行dna提取。用条形码引物515f(5
′‑
gtgccagccggtaa-3
′
)和806r(5
′‑
ggactachvggwtctaat-3
′
)扩增16s rrna基因v4区。对于pcr，按照制造商的说明使用高保真pcr主混合物(美国马萨诸塞州伊普斯维奇新英格兰生物实验室)。循环条件为95℃持续3分钟，然后是30个循环，分别为95℃持续45秒、56℃持续45秒和72℃持续45秒，然后是72℃持续10分钟的延伸步骤。然后根据制造商的方案，使用axyprep dna凝胶提取试剂盒(美国axygen)纯化pcr产物。使用tbs-380荧光计(turner biosystems，ca，usa)测量纯化pcr产物的dna浓度。测序是在北京基因组研究所(bgi，中国深圳)的illumina hiseq2500测序平台上进行的。将重复的dna小份按相同浓度混合到一个dna池中进行测序。
[0066]
原始序列使用flash和trimmomatic软件对测序数据进行去杂，去除引物错配数在2以上、拼接序列的重叠区错配率在0.2以上的序列。剩下的高质量双端序列基于重复区域被连接成tags，用uparse软件(uparse v7，http://drive5.com/uparse/)过滤嵌合子，并对序列进行聚类，根据97％的序列相似度将所有序列进行同源比对并聚类成otus，计算细菌otus相对丰度。总共确定1045个细菌otus相对丰度。
[0067]
(4)采用随机森林算法建立初始回归模型
[0068]
以1045个细菌otus相对丰度作为输入变量，分别以tp和tn/tp作为输出变量，采用随机森林算法建立初始回归模型。输出变量进行对数转换，输出变量值为0时取检出限的1/3进行替代。在每个初始回归模型建立之前，通过10次交叉验证确定关键参数mtry的最佳值，在计算结果达到稳定时确定参数ntree。对于每个初始回归模型，70％的数据集(86个样本)被随机分成训练集，剩下的30％(37个样本)被用作测试集。考虑到数据集的不同划分会影响预测结果的稳定性，设置了10个随机种子，形成10次不同的数据集划分进行10次计算，计算结果以箱线图的形式进行统计分析，以此避免有限数据量和数据集不平衡造成的预测结果不稳定的影响。此外，同时采用基于训练集和基于测试集预测的均方误差(mean-square error,mse)、拟合系数r2两个性能参数对模型进行评估。由于不同指标数值范围不同，mse作为绝对数值的指标，仅适用于同一指标间预测结果的比较。另外，通常mse越小则r2越大，故将r2作为主要评估参数。当r2≤0.3时，预测值与观测值拟合效果差，所述预测模型的预测能力差；当0.3＜r2≤0.4时，预测值与观测值拟合效果弱；当0.4＜r2≤0.6时，预测值与观测值拟合效果中等；当0.6＜r2≤1.0时，预测值与观测值拟合效果强。在随机森林算法通过r软件(v.3.5.2)中的“randomforest”包实现，采用r软件(v.3.5.2)
createdatapartition函数或sample函数进行数据集随机划分。
[0069]
(5)重要变量选择及模型优化
[0070]
由于细菌otus数据集特征很多、样本较少，在建立基于所有细菌otus的模型时，对输入变量进行重要度排序，以选出重要输入变量进行进一步的模型优化。采用平均精度下降(increase in mse，％incmse)评价输入变量重要性。％incmse值越大，变量就越重要。通过r软件(v.3.5.2)中的“importance”和“varimpplot”获得变量％incmse排序。分别选取重要度排序前10、20、30、40、50名的变量作为输入变量，即重要细菌变量，再次按照(4)中的方法进行建模，得到5个优化模型。
[0071]
(6)基于训练集的模型训练效果
[0072]
基于训练集预测tp和tn/tp的均方误差mse和r2统计结果见表2，统计箱线图见图1和图2。
[0073]
表2基于训练集的模型训练效果
[0074][0075]
如表2所示，5种优化模型均能显著降低mse、提高r2，特别是基于重要度排序前20名变量的方法表现最好：tp预测的mse范围在0.83至1.00之间，中位数为0.90，r2范围在0.62至0.73之间，中位数为0.67；tn/tp预测的mse范围在0.67至0.78之间，中位数为0.73，r2范围在0.70至0.76之间，中位数为0.73，说明预测均可以达到“强”拟合效果，且预测结果
稳健。
[0076]
预测tp和tn/tp的前20名变量的具体注释信息分别见表3和表4。
[0077]
表3tp预测模型中重要度排序前20名变量的注释信息
[0078]
[0079]
表4tn/tp预测模型中重要度排序前20名变量的注释信息
[0080][0081]
(7)检验预测模型的预测能力
[0082]
将前述基于训练集得到的初始回归模型和5个优化模型分别用于测试集，以验证模型对于新数据集的泛化能力，采用均方误差mse和决定系数r2进行评估。对测试集进行多
次预测的r2统计结果如表5所示，统计箱线图见图3和图4。
[0083]
表5基于测试集的模型预测效果
[0084][0085]
由表5可以看出，仍然是基于重要度排序前20名的变量所建立的模型预测性能最好，其中tp最好的预测结果如图5(a)所示，基于测试集预测的均方误差mse为0.57，决定系数r2为0.75；tn/tp最好的预测结果如图5(b)所示，基于测试集预测的均方误差mse为0.63，决定系数r2为0.79。
[0086]
(8)预测待测河流水样的tp和tn/tp值
[0087]
分别选择重要度排序前20名变量作为输入变量所建立的tp预测模型或tn/tp预测模型(基于训练集)来预测待测河流水样的tp和tn/tp值：
[0088]
取待测河流水样进行16s rrna测序，确定细菌otus相对丰度；以待测河流水样的重要细菌变量作为输入变量，代入tp预测模型或tn/tp预测模型，得到的输出变量值即为待测河流水样的tp值或tn/tp值。
[0089]
实施例2
[0090]
本实施例提供一种利用细菌otus丰度预测河流水质理化指标tn、no
3-‑
n、nh
4+-n、srp、doc、chl_a、do、orp、tds、cond、ph、temp、no
3-‑
n/tp、no
3-‑
n/tn、nh
4+-n/tn、no
3-‑
n/nh
4+-n、doc/tn、cbz、ery、9-fl的方法，水样同实施例1，各水质理化指标的检测方法、参考标准及
主要仪器设备如表6所示。
[0091]
表6水质检测方法、参考标准及主要仪器设备
[0092][0093][0094]
注：各指标检测结果取三次检测结果的平均值。
[0095]
本实施例的建模方法参照实施例1，得到各水质理化指标的最优预测效果如图6～15所示，由图可以看出，tn、no
3-‑
n、nh
4+-n、doc、chl_a、orp、tds、cond、ph、temp、no
3-‑
n/tp、no
3-‑
n/tn、nh
4+-n/tn、doc/tn、cbz、ery、9-fl对测试集预测的决定系数r2均大于0.6，为强拟合效果；srp、do、no
3-‑
n/nh
4+-n对测试集预测的决定系数r2接近强拟合效果。由此表明本实施例基于重要度排序前10～30名的变量所建立的预测模型可以准确预测以上水质指标。
[0096]
显然，上述实施例仅仅是为清楚地说明所作的举例，而并非对实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或
变动。这里无需也无法对所有的实施方式予以穷举。而由此所引伸出的显而易见的变化或变动仍处于本发明创造的保护范围之中。