北运河浮游细菌集合群落空间变化的环境解释

靳 燕,邱 莹,董 志,邸琰茗,赵栋梁,郭逍宇*（.首都师范大学资源环境与旅游学院,北京 00048；.北京大学生命科学学院,北京 0087；.北京市北运河管理处,北京 000）

浮游细菌是水生生态系统的基本组成部分,具有极高的遗传多样性,并在生物地球化学循环中发挥重要作用[1].作为浮游细菌群落结构多样性构建的两个基本竞争过程,环境选择和扩散过程对浮游细菌在空间变化的贡献率随研究尺度的变化,结论有所差异[2-6].受水流持续扩散和污染物源汇动力学同质化效应影响,扩散作用在较小地理尺度下的自然流动河流中占主导.在高度人工化的水生态系统中,再生水通过补给口与上游来水汇流,形成小尺度地理空间上明显的环境梯度,物种分类作用增强[7].再生水污水处理厂受纳水体微生物群落和多样性的空间变异已成为热点研究主题[7-8].

环境选择作用导致浮游细菌群落中复杂共存的物种间存在互惠、共生、竞争等关系[9-10],随机扩散运动进一步影响物种对特有生境的占有频率[11].二者综合作用,群落中具有较高竞争能力、生态幅度较宽、适应多样生境的物种成为共有核心微生物类群,而具有较强抗逆能力的物种则适应特定的生境类型,成为群落中的特有微生物类群[12-14].集合群落的提出为更好的理解区域尺度范围内的微生物群落组成及其结构动态变化提供了强有力的方法.Jack & Master生态策略假说认为:共有的核心微生物类群在胁迫环境和有利环境下都具有较强的表型可塑性,最终在区域和本底尺度下均表现为较高的丰度和频度.Pandit等[12]发现岩石池中环核心微生物类群变异主要由扩散作用主导,而特有微生物类群变异则主要受环境选择的作用.Sriswasdi等[15]研究显示共有核心微生物类群物种形成率是特有微生物类群的19倍,在维持分类学多样性和系统发育多样性方面发挥着关键作用.Graham等[16]提出扩散和选择作用在多样性构建过程中的相对贡献率对群落生物地化循环功能具有较大影响,高度的扩散速率增加与本地群落中不具有高度适应性的类群比例,从而导致本地群落生物地球化学功能多样性降低.Chen等[17]研究表明在金沙江梯级水库空间和环境因素共同支配着核心和特有分类群的组装和生物地理学,环境因素对核心丰富群落解释度更高.Cao等[18]发现环境因素和地理距离在控制两个连通湖泊区的特有及核心亚群落中起着关键作用.但在环境主导的浮游微生物群落中,具有不同扩散能力和不同生物相互作用模式的共有核心微生物类群是否与特有微生物类群具有相同的环境选择模式未受到关注.基于集合从群落分析再生水污水处理厂受纳水体不同局域群落的空间变化机制尚未见报道.基于此,本文拟探讨高度人工化的再生水补给河道核心浮游细菌类群和特有浮游细菌类群空间变异的环境解释,为揭示高度人工化的再生水补给河道浮游细菌在污染胁迫下群落变化的机制提供科学依据.

1 材料与方法

1.1 研究区域概况

北运河源于北京市昌平区及海淀区一带,后向东南流入通州区,流经河北,在天津市汇入海河,全长 260km.流域总面积为 6166km2,其中,山区面积952km2,平原面积 5214km2.研究区多年平均降雨量 643mm,降雨主要集中在 6～9月,约占全年的84%[19].作为京杭大运河的重要组成部分,北运河主要流经城市化水平高且人口、产业密集的北京、天津,是海河流域骨干行洪排涝河道之一[19].北运河流域土地利用类型和污染源差异导致河流地理空间尺度上存在潜在环境梯度,其次,上、中游密布71座再生水处理厂以排污口形式补给河道,形成小尺度环境梯度.其中,上游流经北京城区,93%的水源来源于城市人口消耗的经处理废水,4%是未经处理的废水;位于京津之间的中游地区主要流经农田,农业化肥和农药污染是中游地区主要的污染物;下游地区进入天津地区,接受来自天津经处理和未处理的废水[20].由于地处平原地区,水流流速缓慢致使水流中裹挟的浊物逐步沉淀,最终导致河道淤积严重[21].

1.2 采样点的布设和样品采集

因非汛期北运河河道断流频繁,水体流通性较差,且汛期湿地植物水质净化特征和非点源污染特征均较为典型.因而,于2019年 9月 3～4日在北运河京津冀段结合周边土地利用类型和均匀布点的基本原则设置样点15个,如图 1所示.每个采样点采集平行的 2L混合均匀的水样 2个,冰袋低温保存并快运回实验室.同时,通过全球定位系统（GPS,Garmin Legend,Garmin USA）定位每个样点的坐标,并采用多参数水质探测仪（AP-7000）测定 pH值、溶解性总固体（TDS）、温度（T）、氧化还原电位（ORP）、电导率（EC）、浊度（Turbidity）.带回室内的水样,定量1L水样充分摇匀（0.25L）通过0.22mm混合纤维素膜过滤,滤膜密封于-80℃保存用于微生物测定.

图1 北运河样点分布Fig.1 Sampling locations in the North Canal River

1.3 水质理化分析

带回室内的水质 24h内完成水质理化分析.水中高锰酸盐指数（CODMn）采用滴定法测定;氨氮（NH4+-N）采用纳氏试剂分光光度法测定;硝酸盐氮（NO3--N）采用紫外分光光度法测定;亚硝酸盐氮（NO2--N）采用N-（1-萘基）-乙二胺分光光度法测定;总氮（TN）采用过硫酸钾消解紫外分光光度法测定;总磷（TP）采用钼酸盐分光光度法测定;总有机碳（TOC）和总碳（TC）采用 TOC 测定仪测定,阴阳离子（Na+、NH4+、K+、Mg2+、Ca2+、Cl-、NO3-、SO42-）采用离子色谱仪测定.

1.4 微生物高通量测序

将已采集好的滤膜送至诺禾致源科技有限公司进行高通量测序分析.扩增引物采用（CCTAYGGGRBGCASCAG,GGACTACNNGGGTATCTAAT）,扩增16S rRNA基因V3-V4区.

1.5 数据预处理与分析

将测序得到的分类操作单元（Operational Taxonomic Unit,OTU）序列按 97%的序列相似度对有效序列进行比对,对 OTU序列进行聚类,同时去除仅出现一次的 OTU以消除假阳性影响,样点均出现的 OTU定义为核心 OTUs（538OTUs）,各样点出现频率＜40%的 OTU定义为特有 OTUs（439OTUs）.利用 Primer 5进行聚类分析和非度量多维尺度分析（NMDS）.基于 Galaxy微生物数据分析平台计算α多样性指数（Shannon指数、Chao指数、Simpson指数、Pielou指数）,利用Origin 2017绘制箱图.利用 Canoco 5.0进行典范对应分析（CCA）.利用 SPSS 24进行非参数检验（Kruskal-Wallis test）和主成分分析（PCA）.

2 结果与讨论

2.1 京津冀北运河水环境理化指标分析

2.1.1 水环境理化指标的统计特征分析 京津冀北运河段各采样点水质参数的描述统计见表1.将水质理化指标与对应的地表水环境质量标准对比后可知,北运河京津冀段 TN均值（3.82mg/L）高于地表水Ⅴ类标准 1.96倍,是北运河应重点加以控制的主要因子;其次,水质 pH值均值为 8.97,接近于地表水Ⅴ类标准 6～9的上限值,最高值高达 9.48,表明北运河京津冀段水质呈现碱性,需重点加以监测控制.结合各阳离子均值看,Na+＞Mg2+＞Ca2+＞K+,而阴离子均值则呈现 Cl-＞SO42-＞NO3-,典型再生水消毒过程导致北运河河道水质以 Na+和 Cl-污染为主,及其诱发的二次有机污染是再生水河道需要重点关注的另一关键因子.结合变异系数可知,NO3--N、SO42-和 Cl-变异系数分别高达 104.74%、111.52%和 110.36%,表明这三个指标具有极强的空间变异性,NH4+-N、NO2--N 和 K+变异系数均大于 75%,属于强变异,表明其离散程度较高,空间分布不均匀.

表1 水质指标的统计特征及水环境质量标准Table 1 Summary statistics of measured parameters and the national standards for surface water quality

2.1.2 水环境理化指标的主成分分析 水质理化指标的主成分分析共提取了 5个旋转因子（VF）,累计解释水质变异 90.23%（表 2）.VF1（方差贡献率为37.22%）与 TN、NO3--N、K+、Ca2+、Na+、Mg2+正相关,表征了再生水补给河道中具有典型再生水水质特征的无机氮和钙、镁盐离子污染.VF2（方差贡献率为24.84%）与T、pH值、TDS、浊度、高锰酸盐指数、TP、TOC、TC、Cl－、SO42-正相关,主要表征再生水补给口物理扰动对水文水质的间接影响.再生水的汇入改变河道水温和酸碱度,并通过再生水补给口水流扰动导致表层底泥悬浮,间接改变水体中浊度和溶解氧,进而改变水质高锰酸盐指数含量.VF3（方差贡献率为 11.76%）和 VF4（方差贡献率为 9.14%）分别与 NH4+-N、NO2--N正相关.再生水中NH4+-N和NO2--N不稳定,短暂存在后转化为较为稳定的NO3--N.VF5（方差贡献率为7.27%）与ORP正相关,再生水汇入引起温度、pH值以及水力学变化,引起表层底泥扰动改变高锰酸盐指数,进而改变ORP.

表2 水质指标主成分分析结果Table 2 Principal component analysis of water quality index

2.2 核心OTUs和特有OTUs局域群落空间变化模式

基于高通量测序的微生物数据经过预处理之后,采用 Primer 进行聚类分析,聚类结果如图 2所示.从图可以看出,不同类型OTUs类型的局域群落结构在空间上均发生较明显的变化,但其空间分布格局模式存在较大差异.以相似性为63%为标准,核心OTUs和特有OTUs局域群落均被划分为两大类.核心OTUs局域群落类型的第一类包含:样点1～11,第二类包含:样点 12～15,核心 OTUs的局域群落聚类结果与地理空间分布格局具有高度一致性.特有OTUs局域群落类型的第一类包含:样点1,3,5,7,8,9,11,13,14,15;第二类包含:样点 2,4,6,10,12;特有OTUs的局域群落聚类结果完全有别于地理空间距离模式.其次,以相似性为25%为标准对不同类型的 OTUs局域群落类型分类结构可以看出,特有OTUs局域群落具有更多的群落类型,即核心OTUs局域群落结构的相似性明显高于特有 OTUs局域群落结构.

图2 基于高通量测序的核心OTUs和特有OTUs局域群落结构聚类分析Fig.2 Dendrogram of hierarchical cluster analysis of bacterioplankton communities based on the high-throughput sequencing

由NMDS结果（图3）可知,核心OTUs局域群落随地理距离的变化,群落类型间存在明显差异（Stress=0.11,Global R=0.887,P=0.001）.特有 OTUs局域群落分类模式与地理空间分布模式无关,基于聚类的不同类群间差异显著（Stress=0.21,Global R=0.473,P=0.003）.其次,基于聚类的核心 OTUs排序图中表现为类内样点间分布集中,类间样点间具有清晰的分布界限和范围,特有 OTUs排序图中类内样点间具有更广泛的分布格局和范围,类间分布界线较弱的特点,从而进一步验证了核心 OTUs局域群落结构的相似性明显高于特有OTUs局域群落结构这一结论.

图3 基于NMDS的核心OTUs和特有OTUs局域群落排序Fig.3 NMDS sequence diagram of the core OTUs and unique OTUs local community

图4 基于核心OTUs和特有OTUs局域群落多样性指数Fig.4 Boxplots of the core OTUs and unique OTUs local community alpha diversity indexes

本研究中高度人工化的北运河浮游集合群落中,核心OTUs和特有OTUs两个局域群落物种个体迁移能力存在显著差异是导致二者地理分布模式存在显著差异的主要原因.一方面核心OTUs局域群落中 OTUs片段的丰度显著高于特有 OTUs局域群落中的OTUs片段丰度,较高的片段丰度具有更高的个体扩散概率[11],从而使得核心OTUs局域群落具有显著的地理空间分布模式.另一方面,核心 OTUs局域群落中较高扩散概率促进其具有更高的生境频率[12,14,22],从而在污染物源汇动力学同质化效应影响下进一步正反馈调节核心 OTUs局域群落较特有OTUs局域群落具有更高的生境频率[7].反之,特有OTUs局域群落分类模式与地理空间分布模式无关,并表现为排污口尺度的变异规律,究其原因一方面与特有OTUs片段丰度较低,个体数不易扩散有关,另一方面,具有较低扩散能力的OTUs片段对源汇动力学同质化效应敏感性较低,北运河河道人工闸坝的存在进一步限制了浮游细菌的扩散能力,尤其是限制特有OTUs局域群落的扩散.

2.3 核心OTUs和特有OTUs局域群落空间变化的环境解释

北运河15个采样点两大类在CCA轴上得到了很好的分化（图5）.基于核心OTUs聚类的第1类和第2类样点沿x轴从左到右的变化趋势.这种变化趋势与北运河水系地理空间距离的变化趋势一致,同时能够很好地表征再生水补给对河道水的空间差异性的影响.与群落分布相关程度最高的环境因子是 TN、NO3--N、K+、Ca2+、Na+、Mg2+.已有大量文献证明再生水中富含 Na+、K+、Ca2+、Mg2+、Cl-、TDS等盐分离子和 TN、NO3--N等多形态无机氮[8,23-27].北运河北京段（上游和部分中游样点）密布再生水补给口,在上、中游集中补给后随径流径向渐变汇入下游海河（部分中游和下游样点）的过程中水中浮游微生物对再生水的净化主要体现在无机氮和盐分离子的降解转化.与群落分布具有较高相关程度的环境因子还包含了高锰酸盐指数和浊度,表征再生水补给对河道水力学特征的空间差异性的影响.城市河道水流缓慢导致沉积物大量淤积,而污染物质的沉降、吸附使得沉积物成为内源污染源[28].再生水汇入影响河流水动力条件,引起流量、流速变化[8].上中游频繁的再生水补给过程扰动河道底部含有大量污染物底泥,使表层污泥中的污染物向上覆水体释放[29],引起河道水体浊度增加.同时,富含丰富有机成分的底泥颗粒使得上中游高锰酸盐指数显著高于下游.作为再生水水质特性中的另外一类关键环境因子,氨氮与亚硝酸盐氮受氧化还原环境的限制,在水中存在形态不稳定[30-32],因而与群落分布相关程度较低,但与 x具有密切的关系.基于特有OTUs聚类的第1类和第2类样点沿y轴主导的对角线方向上从左下到右上的变化趋势.与 y轴相关性最高的环境因子为温度、氧化还原电位.其次,pH值与对角线的变化趋势较为一致.表明这三个环境因子与特有 OTUs群落分类具有密切的关系.取样点2,4,6,10,12靠近再生水补给口是导致第二类群落样点温度显著高于第一类群的直接原因,温度对水体的氧化还原体系有较大影响[33-34],影响水中不稳定变价元素进而进一步反馈促进氧化还原体系空间变异.再生水富含的无机氮磷和盐分离子与 y轴负半轴具有较小的夹角,表征该类环境因仅影响第二类群群落内部空间变异.

图5 再生水补水河道采样点核心OTUs和特有OTUs局域群落与环境因子的排序Fig.5 Sequence diagram of the sample point of the core OTUs and unique OTUs local community and environmental factors on Reclaimed water in river

进一步结合 Kruskal-Wallis检验分析核心OTUs和特有OTUs局域群落内类间显著差异性环境因子.核心OTUs局域群落筛选出TDS、Turbidity、高锰酸盐指数、硝酸盐氮、TN、TOC、TC、Na+、K+、Mg2+、Ca2+等 11个环境因子,这些因子在核心OTUs第1类和第2类样点间具有显著差异（P＜0.05）,并呈下降趋势（表3）;特有OTUs局域群落筛选出T.类间显著差异性环境因子筛选结果与 CCA排序结果一致,且CCA排序较Kruskal-Wallis检验对类间差异性更为敏感.

表3 Kruskal-Wallis检验具有显著差异的水质指标的统计特征Table 3 Kruskal-Wallis test of environmental variables

2.4 核心OTUs和特有OTUs局域群落环境选择模式差异机制

与核心OTUs群落显著的地理分布模式具有较高相关程度的环境因子分别是TN、NO3--N、K+、Ca2+、Na+、Mg2+,这些不同形态的氮素和盐分离子不仅随水流方向呈现显著下降趋势,同时与地理扩散距离表现出较强的协同变化规律.以 TN、NO3--N、K+、Ca2+、Na+、Mg2+主导的环境选择和地理扩散共同作用形成核心OTUs局域群落典型地理空间分布模式.同时表征核心 OTUs局域群落中物种具为宽广的生态幅度,是能适应北运河上、中、下游不同污染特征、不同水力学特征、不同生境特征条件下的具有较高生境频度的广布种.特有OTUs局域群落不具有显著的地理空间分布模式,表明其空间分布格局主要受环境因素主导.主导特有OTUs群落空间变异的主要环境因子为 T,其次为 ORP和pH值.再生水处理厂较高水温和 pH值出水对水体浮游细菌类群形成温度选择压力,进而对浮游细菌群落产生强烈温度和酸碱度选择作用.其次,再生水补给口较强的水力学特征影响河道底泥颗粒物的二次沉降,颗粒物及其富含的有机物质导致出水口氧化还原点位有别于其他样点,进而对浮游细菌群落产生二次选择作用.结合水环境因子主成分源解析表明 TN、NO3--N、K+、Ca2+、Na+、Mg2+与 T和pH值不具有相似的空间变化规律,表征这两类环境因子不同源,从而间接证明了核心 OTUs局域群落和特有OTUs局域群落空间变化具有不同的环境驱动模式.

3 结论

3.1 基于水环境理化指标的主成分分析的污染源源解析结果表明再生水补给带来的TN、NO3--N等无机氮和 K+、Ca2+、Na+、Mg2+盐离子污染是主要污染源,同时再生水汇入引起温度、pH值以及水力学变化,引起表层底泥扰动改变高锰酸盐指数,进而改变 ORP,是再生水补给口物理扰动对水文水质产生间接影响形成的有机污染污染源.

3.2 核心OTUs局域群落和特有OTUs局域群落具有不同的空间变化规律.较高的片段丰度、丰富度和污染物源汇动力学同质化效应共同作用促进核心OTUs具有较强扩散能力,并与地理空间分布格局高度一致.较低的片段丰度限制特有 OTUs扩散能力,从而导致特有OTUs局域群落空间分布模式表现为与距离排污口远近有关的排污口尺度变异.

3.3 核心OTUs局域群落和特有OTUs局域群落空间变化具有不同的环境驱动模式.核心 OTUs局域群落空间分布相关程度最高的环境因子是 TN、NO3--N、K+、Ca2+、Na+、Mg2+,随水流方向显著下降并与地理扩散距离协同变化.特有 OTUs局域群落相关程度最高的环境因子为温度、氧化还原电位和pH值.空间分布格局主要受再生水补给口水力学特征影响所致.