我国工业化和城市化的快速发展, 使得河流等水体普遍受到了不同程度的污染, 水体受到污染后, 河道水流缓慢导致水中污染物部分沉积或通过吸附作用在底泥中富集, 长期累积下使得底泥成为内源污染源, 而底泥中微生物可以通过同化或异化作用对污染物进行降解^[1]。另外, 由于微生物对外界水环境变化较为敏感, 河流水环境质量的变化可使底泥微生物在丰度和群落结构上发生一定的改变, 在微生物作为环境现状指示生物的同时, 河流底泥微生物可以反过来影响外界水环境, 从而使得河流水体生态修复成为可能, 因此研究其群落结构特征具有重要意义。目前高通量测序技术已广泛应用于研究环境微生物的群落结构、多样性及进化关系^[2], 国内外学者也越来越关注河流微生物群落结构与水质的响应关系^[3-5]。

妫水河为延庆的母亲河, 是官厅水库三大入库河系之一。针对妫水河冬季水流动性差、水体富营养化较为严重的问题, 本文采用MiSeq高通量测序方法研究了妫水河不同断面底泥微生物群落结构分布特征, 分析底泥微生物群落与环境因子之间的关系, 研究结果对妫水河水生态环境保护和水质管理提供数据支撑及理论指导意义。

1 材料与方法 1.1 研究区域与样品采集

妫水河流域分布在北京市延庆城区内, 属于永定河水系, 是北京市重要的供水水源河道和水源保护区, 地理坐标为115°49′21″E—116°16′16.75″E, 40°24′52.18″ N—40°33′25.91″N, 发源于延庆城区东北13公里处, 横贯延庆盆地, 在下屯乡大路村北入官厅水库后入永定河, 沿途有支流古城河、三里河和蔡家河汇入, 河长18.5 km, 流域面积1064.66 km², 占延庆区总面积的52%。妫水河流域属大陆季风气候区, 冬季寒冷干燥, 在农场橡胶坝至南关桥段间水体流动性差, 河岸及浅滩区水生植物稀少且分布不均, 水体富营养化较为严重, 为官厅水库水环境改善带来了沉重的负担, 也极大的影响了延庆城区周边的生态环境^[6]。

根据《水质采样方案设计技术规定》(HJ 495-2009)的断面设置原则在妫水河干流设置9个断面(D1—D9), 古城河(DG)、三里河(DS)、蔡家河(DC)三个支流分别各设置1个断面, 共12个, 如图 1。于2017年11月采用多点混合采样法对每个断面进行采样, 每个断面设置距离为10 m的3个采样点, 每个采样点设置3—5个重复。使用有机玻璃采水器采集36个采样点上覆水样500 mL, 利用便携式柱状采泥器采集河流表层底泥(0—10 cm), 底泥样品在去除杂物后, 将各断面样品分别各自混合, 用无菌袋封装并4℃冷藏运送至实验室, 放入-80℃冰柜保存至DNA提取。

1.2 样品理化分析

采用便携式水质分析仪(Thermo Orion A329, USA)现场测定上覆水pH、水温T和溶解氧DO。其他理化指标(包括化学需氧量COD、总磷TP、总氮TN、氨氮NH₄⁺-N和硝态氮NO₃^--N)则在水样采样后24 h之内于实验室按国家标准方法^[7]检测完成。

底泥样品经冷冻干燥后研磨, 称取0.1 g样品于50 mL比色管中, 加入25 mL碱性过硫酸钾溶液, 在0.15—0.16 MPa压力下保持120—124℃的温度30 min, 自然冷却后过滤, 滤液定容到100 mL, 总氮采用紫外分光光度法、总磷采用钼酸铵分光光度法测定^[8]。氨氮和硝态氮按照《HJ 634-2012土壤氨氮、亚硝酸盐氮、硝酸盐氮的测定氯化钾溶液提取-分光光度法》国家标准方法测定, 并取1 g粉末原样用2 mL的1 mol/L盐酸超声波3 h除去无机碳后, 用德国Elemantar vario EL cube元素分析仪进行总TOC/TN的测定^[9]。

1.3 DNA的提取及高通量测定

根据E.Z.N.A.^® soil试剂盒(Omega Bio-tek, Norcross, GA, U.S.)说明书进行总DNA抽提, DNA浓度和纯度利用NanoDrop2000进行检测, 利用1%琼脂糖凝胶电泳检测DNA提取质量, PCR扩增及其高通量测序采用细菌引物为338F：5′-ACTCCTACGGGAGGCAGCA-3′和806R：5′-GGACTACHVGGGTWTCTAAT-3′。扩增程序为：95℃预变性3 min, 27个循环(95℃变性30 s, 55℃退火30 s, 72℃延伸30 s), 最后72℃延伸10 min。所有样本反应结束后, 每个样本3个重复, 将同一样本的PCR产物浓度进行等量混样后使用2%琼脂糖凝胶回收PCR产物, 利用AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, USA)进行纯化, Tris-HCl洗脱, 2%琼脂糖电泳检测。利用QuantiFluor^TM-ST (Promega, USA)进行检测定量。利用Illumina公司的Miseq PE300平台进行测序, 测序委托北京中科晶云科技有限公司进行, 对原始序列进行预处理, 根据97%的相似度对序列进行OTU聚类, 选取高质量的序列进行分析。

1.4 数据分析

高通量测序所得的序列通过QIIME进行处理, 主要包括OUT分类、多样性分析和PCoA聚类分析；采用CANOCO 4.5进行线性回归的冗余分析(RDA)分析, 获取底泥微生物群落对环境因子之间的相关性。

2 结果与分析 2.1 底泥及上覆水的理化性质分析

如表 1所示, pH整体呈弱碱性, 各样点底泥的pH值基本上低于上覆水。妫水河中下游COD、TN含量未满足“水十条”标准, D5和DS点NH₄⁺-N均为Ⅳ类水标准。其中, 支流三里河DS点污染严重的水体流入妫水河主干流是导致D8点COD、NH₄⁺-N、TN污染超标的重要原因之一, TN含量整体偏高可能是由于冬季水温低, 微生物的活性较低, 降低了生物作用对氮磷的去除效果^[10], 且秋冬季降雨补给减少, 污染稀释作用降低。由表 1、表 2对比可知上腹水TN含量主要超标点位由高到低依次为：DS＞D5＞D4＞D6＞D7, 底泥中TN含量变化顺序为：DS＞D5＞D4＞D6＞D7, 由相关性分析可得, 各点位之间上覆水TN含量与底泥TN含量呈极显著正相关(P=0.914)。通过C/N比示踪环境与气候演变的过程, 可以得出沉积物有机质的来源信息, 通常有机质中C/N比介于10—13, 表示来源为水生植物, 小于10为藻类^[11], 因此根据表 2中C/N比可得出妫水河底泥污染主要来源为水生植物和藻类。

2.2 优势微生物菌群分布

利用Miseq测序技术研究妫水河底泥微生物群落结构, 共检出70门228纲1168属, 其中丰度最大的10个菌门为变形菌门(Proteobacteria)、拟杆菌门(Bacteroidetes)、酸杆菌门(Acidobacteria)、放线菌门(Actinobacteria)、绿弯菌门(Chloroflexi)、厚壁菌门(Firmicutes)、蓝细菌门(Cyanobacteria)、绿菌门(Chlorobi)、疣微菌门(Verrucomicrobia)和硝化螺旋菌门(Nitrospirae), 如图 2所示。以上菌门之和在各个样品中均占84%以上, 但这些菌门的分布有所不同。变形菌门(Proteobacteria)在各个样品中均是第一优势门, 而在丰度上也存在差异, 其中D7样品最高, 占69.1%, 而D5样品相对较低, 占45.3%。

对12个底泥样品的变形菌门微生物的分布特征进行分析, 结果见表 3。可以看出, β-变形菌纲(β-proteobacteria)和δ-变形菌纲(δ-proteobacteria)是变形菌门丰度最大的两个菌群, β-变形菌纲约占细菌总量的11.01%—31.50%, δ-变形菌纲约占细菌总量的11.6%—20.24%, 除此之外, γ-变形菌纲(γ-proteobacteria)和α-变形菌纲(α-Proteobacteria)也占有较高的丰度比例, 而ε-变形菌纲(ε-proteobacteria)的丰度最低。DC样品中α-变形菌纲(α-proteobacteria)含量偏高, 达到12.8%。

从属水平上对微生物的占比情况进行统计, 其中相对丰度大于1%的包括20个属(平均丰度), 如图 3所示, 分别是Bacteroidales_norank、Sinobacteraceae_norank、脱氯单胞菌属(Dechloromonas)、Comamonadaceae_norank、Ellin6067_norank、地杆菌属(Geobacter)、Stramenopiles_norank、iii1-15_norank、SC-I-84_norank、MND1_norank、Syntrophobacteraceae_norank、Myxococcales_norank、硫杆菌属(Thiobacillus)、Dok59、Ignavibacteriaceae_norank、Betaproteobacteria_norank、Gemm-1_norank、Syntrophaceae_norank、Deltaproteobacteria_norank、Xanthomona-daceae_norank。其中, 脱氮单孢菌属、Comamonadaceae_norank、SC-I-84_norank、硫杆状菌属、Ellin6067_norank、Dok59、MND1_norank均隶属于β-变形菌纲。

不同断面样品的优势种属有所差别。其中DG、D5的优势菌属为Bacteroidales_norank, D3、D6的优势菌属集中在Sinobacteraceae_norank, D4、D7含有较高比例的脱氮单孢菌属微生物, 但上游点位D1、D2和下游点位DS、D8、DC、D9的优势菌属比较分散, 分别为BPC076_norank、Comamonadaceae_norank、Dok59、硫杆状菌属、iii1-15_norank、MND1_norank。

2.3 微生物群落多样性及差异性分析

如表 4所示, 本研究将相似度为97%的序列划分为1个OTU, 并进行多样性指数计算。由Ace丰富度指数可知, D7样品的丰富度较其他样品偏低, 低至2673, 而DG古城河样品的Ace指数最高, 达4869, D5样品其次, 为4826, 表明DG、D5样品中物种总数也相对最高, Chao丰富度指数所表现趋势和Ace指数相同。各样品Shannon指数均大于6.0, 其中DC蔡家河样品的Shannon指数值最高, 为8.72, 其多样性相对较高, D1、D2、D3、D6、D7点Shannon指数均小于8.0, 其中D7点最小, 为6.56。

为了进一步分析样本间在群落物种组成上的差异性, 基于weighted Unifrac^[12]采用PCoA分析衡量样本间群落物种组成的相似度, 结果表明12个样品主要分布在两个不同的位置。如图 4所示PC1是造成样品差异性最大的主坐标成分, 解释度为35.08%, 其次为PC2, 解释度为25.83%, 对物种分布的总解释量为60.91%。在PC1和PC2维度上, 将采样地分为两组, 妫水河中游点位D2、D3、D4、D6、D7、D8和DS微生物群落结构相似度较高为一组, 上游点位D1、DG及下游点位D9、DC相似度较高为一组, 进而将两个分组的Shannon指数进行方差分析, 结果也表明两组多样性存在显著性差异(P=0.037)。

2.4 微生物群落结构与环境因子关系的冗余分析

基于12个样品优势属微生物群落信息和底泥主要环境因子进行冗余分析(RDA), 考察微生物群落因子与环境因子变量之间的关系。图 5中RDA模型分析发现, 两主坐标分别解释了43.1%和21.3%的底泥微生物群落因子与环境因子关系, 共解释了微生物属水平64.4%的分布差异。NH₄⁺-N和T对妫水河底泥微生物群落结构影响最大(F=2.92, P=0.038；F=2.81, P=0.014)。值得注意的是, NH₄⁺-N与Bacteroidales-norank相关性最大, 脱氮单孢菌属和硫杆状菌属均与NH₄⁺-N、T呈正相关, 同时与DO呈负相关。

3 讨论 3.1 妫水河不同断面微生物群落结构特征

妫水河第一优势菌门为变形菌门(Proteobacteria), 在水环境中普遍存在。相关研究表明其在湿地、养殖废水和深海沉积物等环境中均以第一优势门存在^[13-15], 大多数在生物脱氮、生物除磷及诸多污染物降解过程中起重要作用的微生物均归属于变形菌门^[16]。在纲水平上, β-变形菌纲占比最多, 达到11.01%—31.50%, 其广泛存在河流和湖泊中(如松花江^[17]、海宁长山河^[18]、东江^[19]、Prealpine湖^[20]等), 并在氮磷等污染物的去除中起到了重要的作用^[21]。

妫水河各样品Shannon多样性指数均较高, 与微生物多样性指数在水体底泥表层通常最高有关^[22]。D7点Shannon多样性指数较其他点位低, 而D7点变形菌门丰度最高, 经实地调研这可能是由于D7点位河面较宽, 水流速极小(小于0.02 m³/s), 微生物迁移速度慢导致菌种单一稳定, 多样性小。妫水河微生物群落结构在12个样品间存在断面差异, 妫水河中游点位D2、D3、D4、D6、D7、D8和DS微生物群落结构相似度较高, 这与中游点位均分布在延庆城区内, 受人类活动影响大有关。下游D9入官厅水库点位和支流DC蔡家河点位差异性较小, 与两个点位地理位置相近、附近环境相似有关, 而在三条支流中古城河DG点和蔡家河DC点具有一定的相似性, 是由于两条支流附近均有村庄分布, 与城区有一定距离, 受农业面源污染影响大, 因此生态环境相似度高, 微生物群落结构差异性较小。

3.2 妫水河微生物群落生态功能及对环境因子的相关性

底泥微生物群落结构会随着环境因子的变化呈现出一定的空间分布特征, 而微生物群落结构的变化又会在一定程度上反映环境因子的变化情况^[3]。王娜等^[23]在研究太湖沉积物中微生物生物量与碳、氧、磷等主要元素关系时发现, 太湖沉积物TOC/TN的变化显著影响微生物群落结构。鲍林林等^[24]通过高通量测序研究得出氨氮是影响河流底泥氨氮化微生物群落结构特征的主要因子, 与本研究结论一致。刘吉文等^[25]在研究典型海域微生物群落结构中得出DO对细菌属水平上具有重要显著影响, 本研究中Bacteroidales-norank与DO呈正相关, 而Bacteroidales-norank多存在与生物粪便内, 常被作为粪便来源标识物^[26], 作为DG、D5点的优势菌属, 与DG点农村动物粪便、D5点城区生活污水排放相关。DC蔡家河点位优势菌属iii1-15_norank隶属于酸杆菌门(Acidobacteria), 与pH相关性最大, 并呈正相关, 这与酸杆菌是生长较缓慢的寡营养菌, 其丰度受pH影响较大相符, 且pH对不同亚群相关性影响不同^[27], 有研究表明氮肥会降低酸杆菌细菌的丰度^[28], 本研究中iii1-15_norank与pH正相关, 这可能与施肥导致pH下降有关^[29]。

不同环境因子对微生物群落丰度具有不同的影响。pH可以通过H⁺浓度改变土壤中营养元素的形态从而影响自养微生物类群^[30], 妫水河不同断面之间pH值差异性不大, 呈弱碱性。妫水河上游古城河支流DG点硝态氮含量较其他点位高, 但其微生物群落丰富度高, 可能与上游农业面源污染有关, 在长期施用氮肥累积下硝态氮含量升高, 从而使微生物的数量和生物量增多^[31]。经实地调研妫水河D8谷家营国控断面点位河面很宽, 水位较高, DO相对较小, D8点优势菌属为硫杆状菌属(Thiobacillus), 占细菌总数达到3.10%, 同为浑河底泥优势菌属, 在厌氧条件下以NO₃^-为电子受体, 在氧化硫化物的过程中将其转化为N₂, 可实现脱硫反硝化^[32]。并且本研究中硫杆状菌属(Thiobacillus)与温度相关性较大, 当水温升高时, 营养物质丰富的环境中微生物生长较旺盛, 而硫杆状菌属主要参与反硝化过程^[33], 因此反硝化速率会有所升高。Zhao等^[34]在太湖流域研究中得出沉积物反硝化速率随NH₄⁺含量的增加而升高, 而且当高含氮污水排入河流时, 氮素的去除效率也相应的提高, 本研究中反硝化菌^[35]脱氮单孢菌属(Dechloromonas)与氨氮含量呈正相关, 同时与DO呈负相关, 这与反硝化菌群在高DO条件下的反硝化能力会受到抑制有关^[36]。作为妫水河底泥中主要反硝化菌属, 硫杆状菌属和脱氮单孢菌属在不同断面中的平均占比和达到4.55%, 妫水河水体流速缓慢, 为底泥中反硝化细菌提供了较好的缺氧环境, 而在冬季温度较低条件下底泥微生物反硝化速率减弱。因此, 通过分析微生物群落及丰度变化与环境因子的相关性, 反映水环境因子的变化, 体现水生态环境状况的好坏, 并为水质改善、水生态保护和管理提供指导意见。但基于河流的季节性变化, 妫水河底泥微生物群落结构与环境因子相互作用的具体机制较为复杂, 还需今后进行更深入全面的探究。

4 结论

(1) 妫水河中下游水质偏差, 上覆水中营养物质对底泥污染存在一定的影响, 基本呈正相关性；支流三里河水质和底泥污染较严重, 对妫水河干流水质影响较大, 因此, 要进一步加大对支流水源的监督与控制。

(2) 妫水河底泥微生物共检出70门228纲1168属, 作为第一优势菌门, 变形菌门(Proteobacteria)在各样品中占比高达45.3%—69.1%, 其中β-变形菌纲(β-proteobacteria)和δ-变形菌纲(δ-proteobacteria)是变形菌门丰度最大的两个菌群, 占比分别为11.01%—31.50%和11.6%—20.24%；妫水河底泥微生物群落丰度总体较高, 多样性也相对较高, 但不同样品之间存在差异性, 其中世园段D7点Ace丰富度指数和Shannon多样性指数均较其他点位低。

(3) 妫水河底泥微生物群落受NH₄⁺-N和温度的影响较大, 妫水河优势反硝化菌属脱氮单孢菌属和硫杆状菌属的丰度与NH₄⁺-N、温度均呈正相关, 与DO呈负相关。在冬季温度较低条件下底泥反硝化速率减弱, 通过分析微生物群落及丰度变化与环境因子的相关性, 加强妫水河冬季对水质的管理与监测, 可为水环境治理和水质改善提供指导意见。

致谢: 感谢罗明科、李阳、王源在野外采样中给予的帮助。