辽河口沉积物反硝化过程研究
——以反硝化功能基因丰度及nirK 型细菌群落结构分析为例

明红霞，陈泉睿,3，史银银,，苏洁，于颖，樊景凤*

( 1. 国家海洋环境监测中心 国家生态环境保护近岸海域生态环境重点实验室，辽宁 大连 116023；2. 大连海事大学 环境科学与工程学院，辽宁 大连 116026；3. 厦门大学 海洋与地球学院，福建 厦门 361005)

1 引言

氮元素作为生物体所必需的营养元素，在生态系统中发挥着不可替代的作用，是限制海洋生态系统生产力的关键因子之一[1]。氮的转化途径复杂多变，反硝化过程是全球近岸海洋沉积物氮循环脱氮的关键反应过程，异养细菌利用硝酸盐作为呼吸末端的电子受体，将海洋环境中的硝酸盐和亚硝酸盐还原为N2O 和N2两种气态产物，主导着海洋生态系统的脱氮过程，进而去除水体中的氮元素，控制富营养化发生速率[2]。主导反硝化过程的主要功能基因有：调控硝酸盐向亚硝酸盐还原的narG 基因、将亚硝酸盐还原为NO 的nirK 基因、控制NO 还原为N2O 的norB基因、使得N2O 生成N2的nosZ 基因。

Helen 等[3]从反硝化过程的角度认为，参与反硝化的还原剂至少应当能够完成NO2−到气态NO 的转化，并从中获取能量，该过程离不开亚硝酸还原酶的作用。亚硝酸还原酶将亚硝酸盐还原为绝大多数生物都不可利用的气态NO。常见的亚硝酸盐还原酶有两种，分别是由nirS 功能基因编码的细胞色素cd1 亚硝酸还原酶和由nirK 基因编码的Cu 依赖性亚硝酸还原酶，它们功能相同但结构不同[4]。通常认为nirS 型反硝化菌具有更高的基因丰度，在环境中占据了主导地位，而nirK 型反硝化菌则存在于更为广泛的微生物种群中，涵盖了更为多样的分类单元，nirK功能基因所主导的反硝化过程在环境中的地位也一直被低估[4–5]。研究表明，nirK 型反硝化菌在对不同环境梯度的响应中表现出了更为强大的生境选择[6–7]。对nirK 型反硝化细菌的研究，可以帮助我们进一步认识和了解nirK 型反硝化细菌在环境中的独特地位，进而对近岸海洋氮循环过程做出更加正确、合理的推断。

河口是连接河流与海洋的过渡地带，是水圈、岩石圈和大气圈的重要交汇地区，也是微生物地球化学循环最为活跃的地区[8]。河口沉积物是海洋微生物反硝化作用的重要发生地和富集地，也是氮素及营养盐重要的源和汇[8]。辽河口位于辽宁省盘锦市河流入海的过渡地带，属暖温带季风气候区，降水多集中于7−9 月[9]。辽河口河流含沙量较高，水文及理化条件复杂，环境因子多变。辽河口区域为当地经济生产活动提供了丰富的海洋和生物资源，而大量陆源营养物质随河流搬运作用进入辽东湾海洋生态系统，引发了辽河口水质下降、水体富营养化[10]等，这些污染对河口氮循环微生物群落造成了不同程度的干扰。

本研究以荧光定量PCR 的方法测定了辽河口沉积物反硝化过程关键功能基因的丰度，采用高通量测序技术对nirK 型反硝化细菌群落结构和多样性进行解析，并进一步分析了细菌群落与环境因子之间的相关性，可以帮助我们进一步探究反硝化细菌在河口生态环境中的分布及其生态策略，以期为河口或近岸海域水资源保护和生态环境治理提供科学的理论指导。

2 方法

2.1 样品采集

实验样品采集于2016 年7 月，在辽河口区域按照河口向海延伸方向设置12 个采样站位，按照盐度梯度划分为3 组（盐度0.8～7.0 为近岸淡水组，图1 中C 组；盐度13.1～20.7 为混合盐度组，图1 中B组；盐度25.5～31.5 为高盐度组，图1 中A 组），采样站位布设如图1 所示。使用有机玻璃采水器收集各站位的底层水，保存于无菌玻璃瓶中，在4℃条件下保存并运回实验室，用以测量其他理化参数。在各个采样站位随机抓取3 斗海洋沉积物，将3 斗沉积物样品混匀均质化形成复合样品，用无菌密封袋保存，在−20℃冷藏运输回实验室，分装后部分样品于−80℃长期保存。

图1 辽河口12 个采样点站位布设及分组情况Fig. 1 Distribution and grouping of 12 stations in Liaohe Estuary

2.2 环境参数测定

将YSI 水质测量仪（SG23-FK-ISM，上海）直接放进水层底部，实时测量水温、盐度、pH 和电导率等数据。实验室利用碘量法测定底层水溶氧值，萘乙二胺分光光度法测定 NO−2含量，锌镉还原法测定 NO−3含量，次溴酸盐氧化法测定 NH+4含量，磷钼蓝分光光度法测定 PO34−含量。使用激光粒度分析仪（LS13320SW，美国）测定沉积物粒度，重铬酸钾氧化还原容量法测定沉积物中总有机碳含量，分光光度法测定总磷含量，凯式滴定法测定总氮含量。本研究的环境参数与文献[11]一致（详细参数见表1，表2）。

2.3 DNA 提取及qPCR 技术

在各采样站点取0.35 g 海洋沉积物，使用DNeasy PowerSoil Kit 试剂盒（Qiagen，德国）提取沉积物中的DNA，各站位取3 组平行样品，操作过程严格按照试剂盒提取步骤进行。使用微量紫外分光光度计（Scandrop 200，德国）测定DNA 的纯度和浓度，使用Goldview 核酸染料（Solarbio，北京）对提取到的DNA染色，以1.5%的TAE-琼脂糖凝胶电泳检测DNA 完整性。分别构建各功能基因的质粒，采用荧光定量PCR 仪（ABI-7 500 型号，APPLIED BIOSYSTEM）依次测定组成反硝化过程的4 个主要功能基因：硝酸盐还原酶基因narG、亚硝酸盐还原酶基因nirK、一氧化氮还原酶基因norB 和一氧化二氮还原酶基因nosZ 的丰度，功能基因的引物见表3。

qPCR 扩 增体系为：10 μL 的TB Green Premix Ex Taq II（宝生物，大连），上下游引物（10 mmol/L）各0.4 μL，ROXII（宝生物，大连）0.4 μL，DNA 模板2 μL 和6.8 μL的ddH2O（宝生物，大连），设置阴性对照和3 组平行。梯度PCR 确定narG、nirK、norB、nosZ 4 个反硝化功能基因引物的最佳退火温度分别为51℃、57℃、55℃和59℃。选择插入所需目的基因DNA 片段的质粒建立标准曲线，通过观察熔解曲线检验扩增产物的特异性。qPCR 技术测定目的基因的扩增效率应在80%～105%之间，标准曲线的相关系数大于0.99。

2.4 高通量测序技术

将nirK 型反硝化功能基因作为反硝化细菌的代表功能基因，对其进行高通量测序。使用与qPCR 相同的特异性引物对nirK 基因进行PCR 扩增。扩增产物使用AxyPrepDNA 凝胶回收试剂盒（宝生物，大连）回收目的条带。委托上海微基生物科技有限公司使用2×300 bp 末端碱基配对，通过Illumina MiSeq 平台对nirK 基因文库进行高通量测序。测序完成后，根据barcode 标签筛选有效序列并区分样品来源，对有效序列进行定量质量控制和序列拼接，优化基因序列，去除低质量序列和小于50 bp 的短序列。

2.5 数据统计分析

使用软件 USEARCH 将有效序列确定为分类操作单元（Operational Taxonomic Units，OTUs），序列相似性在97%以上的OTUs 进行聚类，分析细菌群落结构的α 多样性。Chao1 指数、ACE 指数、Shannon 指数、Simpson 指数、Coverage 指数通过Mothur v.1.32.1对OTU 进行计算，分析环境样品中OTU 数目组成的相似性和重叠情况。使用R 语言建立韦恩图、花瓣图以及RDA 冗余分析图。采用SPSS 23.0 对环境参数与基因丰度进行Pearson 相关分析，显著性水平设定为p<0.05。通过R 语言及语言包完成分析多样性指数和环境因子之间的关系，LEfSe 分析采用线型判别分析。通过LEfSe 软件进行组间群落差异分析。将Illumina 测序原始序列提交至GenBank。

表1 辽河口底层水中环境因子理化参数[11]Table 1 Physical and chemical parameters of environmental factors in bottom water of the Liaohe Estuary[11]

表2 辽河口表层沉积物中环境因子理化参数[11]Table 2 Physical and chemical parameters of environmental factors in surface sediments of the Liaohe Estuary[11]

表3 反硝化功能基因的引物Table 3 Primers for denitrification genes

2.6 基因序列登录号

本研究中高通量测序得到的nirK 型反硝化功能基因的原始序列已经登陆至NCBI 数据库，基因序列登 录 号 为SRR8907784、SRR8907785、SRR8907786、SRR8907787、SRR8907789、SRR8907790、SRR8907791、SRR8907792、SRR8907793、SRR8907794、SRR8907795和SRR8907796。

3 结果

3.1 反硝化功能基因在辽河口表层沉积物中的丰度分布

采用qPCR 技术的方法测定narG、nirK、norB、nosZ 4 个功能基因在辽河口沉积物中的丰度，结果如图2所示。narG 型反硝化功能基因的丰度分布在1.63×108～1.13×109copies/g，平均丰度为5.92×108copies/g，在B1 站位反硝化功能基因丰度最高，在C4 站位基因丰度最低。nirK 型反硝化功能基因的平均丰度为2.33×107copies/g，在B1 站位丰度最高，为5.58×107copies/g，在C4 站位丰度最低，为2.12×106copies/g。norB 型反硝化功能基因的平均丰度为9.13×106copies/g，C3 站位基因丰度最高，为2.72×107copies/g；A2 站位基因丰度最低，为2.60×106copies/g。nosZ 型反硝化功能基因平均丰度为1.35×108copies/g，在B1 站位丰度最高，为1.03×108copies/g；C4 站位丰度最低，为6.48×107copies/g。

从分组数据可以看出，在辽河口12 个站位中narG基因丰度最高，norB 基因丰度最低。narG、nirK、nosZ 型反硝化功能基因的相关性较好，其丰度的最高值都出现在淡咸水混合组，最低值出现在低盐度组。norB 型反硝化功能基因丰度表现出随盐度升高而丰度降低的趋势，nosZ 型反硝化功能基因的丰度与narG、nirK 相比分布更为均匀，在按盐度梯度划分的分组中差异最小。

对辽河口12 个站位的环境因子与功能基因丰度进行皮尔逊相关性分析，图3 中红色表示两个变量呈正相关关系，绿色为负相关关系。在功能基因丰度上，narG 与nirK、narG 与nosZ 以及nirK 与nosZ 之间均具有显著相关性（p＜0.05），其中narG 与nirK 型基因丰度与沉积物粒径也表现出了显著相关性，与其他环境因子无显著相关性，而调控着NO 向N2O 转化的norB 基因丰度与盐度、pH、 NO−2含量、 PO34−含量具有显著相关性，与N H+4含量、N O−3含量具有极显著相关性。

3.2 nirK 型反硝化细菌多样性

对辽河口12 个采样点站位表层沉积物进行nirK型反硝化功能基因高通量测序，测出有效序列数、数据覆盖度，对所测得的数据进行OTUs 分类，并计算其丰度和多样性指数，结果如表4 所示。

图2 辽河口narG、nirK、norB、nosZ 功能基因丰度Fig. 2 Functional gene abundance of narG, nirK, norB, and nosZ in the Liaohe Estuary

图3 辽河口环境因子与功能基因丰度的皮尔逊相关性Fig. 3 Pearson correlation between environmental factors and the abundance of functional genes in the Liaohe Estuary

样品中的Coverage 指数最小值为99.89%，表明高通量测序的结果基本覆盖了所有的样本信息，可以客观地反映样品中nirK 型反硝化细菌群落的真实信息。Chao1 指数和ACE 指数可以反映出基因丰富度。由表4 可以看出，在低盐度C 组中菌群丰富度较高，高盐度A 组中菌群丰富度较低。Shannon 指数和Simpson指数则可以估算样品中细菌群落的多样性，二者皆表明低盐度组具有更为丰富的物种多样性，高盐度组微生物物种群落结构较为单一。本次高通量测序在12 个站位中共得到417 272 个有效序列，将这些序列进行聚类分析，一致性在97%以上的聚为1 个OTU。对辽河口nirK 型反硝化功能基因的测序共得到619 个OTUs 序列，从韦恩图（图4a）中可以比较直观的看出各采样点中OTUs 数目组成的相似性和重叠情况。高盐度A 组、混合区B 组、低盐度C 组分别含有219 个、328 个、480 个OTUs，3 组样品共有的OTUs 是126 个。低盐度组中OTUs 数最多，其次是淡咸水混合组，高盐度组的OTUs 最少，原因可能是盐度的胁迫筛选了特定的菌属。花瓣图（图4b）可以表示出12 个采样站位中独有和共有的OTUs 数目，其中最大OTUs 值为305 个，出现在C2 采样点，最小OTUs 为96 个，出现在A2 站位；12 个采样点共有的核心OTUs 为30 个。

表4 辽河口nirK 型反硝化细菌丰度及多样性Table 4 Abundance and diversity of nirK-type denitrifying bacteria in the Liaohe Estuary

3.3 nirK 型反硝化细菌群落结构组成及差异分析

根据高通量测序结果对OTU 进行注释，发现在辽河口nirK 型反硝化细菌共有7 个门、9 个纲、23 个目、36 个科、63 个属、121 个种。其中变形菌门（Proteobacteria）占绝对优势地位，在各站位菌门中最低比例为99.37%，平均值为99.83%。在属水平上，优势菌主要为德沃斯氏菌（Devosia）、海洋嗜酸杆菌（Phaeobacter）、产碱杆菌（Alcaligenes）、假单胞菌（Pseudomonas），约占总菌属的68.97%，详见图5。

图4 nirK 型功能基因OTUs 分布的韦恩图（a）和花瓣图（b）Fig. 4 Venn diagram (a) and petal map (b) of the nirK functional gene OTUs distribution

从图5 中可以看出，德沃斯氏菌和产碱杆菌在高盐度组向低盐度组的过渡中丰度逐渐下降，海洋嗜酸杆菌在混合盐度区拥有最大丰度，而假单胞菌则在近岸淡水组呈现最大丰度。

按盐度梯度的划分，对nirK 型反硝化细菌群落进行线性判别分析，构建了细菌的进化分支树图和线性判别分析（LDA）值分布柱状图（图6）。进化分支树图由内向外的每层圆弧分别代表了门、纲、目、科、属、种，共6 个水平层次上的细菌群落结构分布。发育树各树枝节点的颜色代表细菌的分组情况，节点的直径大小与该物种的相对丰度正相关。LDA 值分布柱状图显示了大于预设值2.0 的具有显著差异性的物种。将图6 各分组中的生物标志物对应到属水平，可以得出产碱杆菌、德沃斯氏菌、假单胞菌、根瘤菌（Rhizobium）、无色杆菌（Achromobacter）为A 组的生物标志物，在高盐度中发挥着重要的作用。而低盐度C 组中起着重要作用的菌属为假单胞菌、Defluviimonas、螯 台 球 菌 属（Chelatococcus）、 中 慢 生 根 瘤 菌（Mesorhizobium）、胺杆菌属（Aminobacter）、副球菌属（Paracoccus）、褐杆菌属（Phaeobacter）等。我们在淡咸水混合的B 组没有检出具有统计学差异的生物标志物。”

选择冗余分析（Redundancy Analysis，RDA）图（图7）分析不同类别OTUs 的组成来反映样本间的差异，可以得到辽河口nirK 型反硝化细菌群落组成与张慧珍等[11]测得的环境因子之间的相关性，其中RDA图的总贡献度为80.93%，在很大程度上可以解释环境因子变量与采样站位细菌群落结构之间的关系。盐度、pH、溶解氧和亚硝酸盐等环境因子对细菌微生物群落有更明显的影响作用。德沃斯氏菌、产碱杆菌、Breoghania的群落丰度与盐度、pH 正相关，与NH+4、 PO34−、 NO−2含量负相关。海洋嗜酸杆菌与 NO−3含量正相关，与DO 含量呈现出明显的负相关，假单胞菌与盐度具有显著的负相关关系。其中，N O−2含量、pH、NH+4对nirK 型细菌群落结构的影响具有统计学意义（p<0.05，蒙特卡洛检验1 000 倍）。

图5 nirK 型反硝化细菌属水平的群落结构柱状图Fig. 5 Histogram of community structure at the genus level of nirK-type denitrifying bacterias

图6 基于LEfSe 分析的辽河口沉积物中nirK 型反硝化细菌的丰度柱状图(a)和进化分支图(b)Fig. 6 Abundance histogram (a) and evolutionary branch map (b) of nir K-type denitrifying bacteria in the sediments of the Liaohe Estuary based on LEfSe analysis

图7 nirK 型反硝化细菌群落与环境因子变量的冗余分析图Fig. 7 Redundancy analysis of nirK-type denitrifying bacteria community and environmental factor variables

4 讨论

4.1 辽河口反硝化基因丰度分析

功能基因的丰度可以反映辽河口潜在的反硝化能力。narG 基因丰度最高，表明硝酸盐还原作用在辽河口反硝化过程中的重要地位，较高的基因丰度也可能是由部分细菌基因组内存在多个基因拷贝数所引起的；norB 基因随河口盐度升高而丰度降低，呈现出明显的规律变化，我们推测在辽河口区域NO 的还原过程主要发生在淡水输入比较突出的低盐度区域。narG、nirK 和nosZ 3 个功能基因的丰度均在盐水−淡水的过渡混合区沉积物中最高，其中nosZ 在各盐度组中的丰度变化最小，其原因可能为混合区复杂的环境条件为其提供了更加多样化的生存条件。对反硝化功能基因的丰度进行相关性分析，结果表明narG、nirK 和nosZ 功能基因之间表现出两两显著相关，进一步说明了辽河口 NO−3还原和N2生成的潜力较高，与环境因子的相关性较弱，受河口环境因子变化的影响较小。皮尔逊相关性分析中也指出narG 和nirK功能基因的丰度与沉积物粒径呈显著相关性（图3）。沉积物粒径与原位水文条件密切相关，可以通过河流径流、淡咸水混合强度，对其他物理化学特性的动态控制，进而影响反硝化细菌群落的丰度和多样性[15]。Jiang 等[16]在中国长江口的研究结果也表明沉积物类型是预测nirK 型反硝化细菌丰度和群落组成的关键因子。虽然反硝化是由多个不同的功能基因所调控，但是反硝化作用的发生可能是由环境中某一类优势菌群所主导，从高通量测序的注释结果可以看出在辽河口反硝化菌群中99.37%以上的菌群均隶属于变形菌门，这也可能是反硝化功能基因丰度之间具有较高相关性的原因。

研究中，个别站点的局部差异影响功能基因丰度的分布，但该差异并没有造成大范围区域性的影响。这提醒我们在今后的研究中设置更多的采样站位，以消除局部空间误差对实验带来的影响，从而更好地解释河口地区复杂生境下功能微生物的分布和群落特征。

4.2 辽河口nirK 型反硝化细菌群落生态功能分析

对nirK 型反硝化功能基因的高通量测序发现，变形菌门（Proteobacteria）在辽河口亚硝酸盐还原菌中的丰度和多样性最高，在辽河口沉积物中变形菌门为起到重要作用的优势菌群，这与Yang 等[17]在辽河口的研究一致。此外，本研究中所使用的nirK 引物主要扩增出α-变形菌纲的基因片段，但是 Helen 等[3]也指出，nirK 型反硝化细菌还存在于更广泛的细菌分类中，包括以放线菌门、厚壁菌门等为主的Clade II 类群。辽河口环境中也可能存在着更高的反硝化菌丰度和更丰富的群落多样性，其中德沃斯氏菌、海洋嗜酸杆菌、产碱杆菌和假单胞菌占据较大的丰度。

有研究表明，在突尼斯北部海岸[18]、智利中部金特罗湾（Quintero Bay）[19]被石油污染的海洋沉积物中同样分离出了产碱杆菌，且该菌种可利用石油、原油以及其他碳氢化合物作为唯一碳源，产生生物表面活性剂，同时具有重金属抗性基因，可以适应有毒化合物、高盐、氧化胁迫等环境条件，其生物降解能力可用于原位生物修复等[18]。辽东湾北部区域蕴藏着丰富的石油和天然气资源，陆源营养物质的输入对沉积物中的菌属进行了特定的筛选，这些特定菌属的大量存在可能对海洋环境中石油和天然气等资源的分布具有潜在的指向作用。在太平洋东海岸[20]、韩国所测的东海的淡咸水交界处[21]以及波兰地区的东南部[22]也都检测并分离出了德沃斯氏菌，且证实了该菌属与水生植物的固氮和结瘤共生有着密切的联系[23]。德沃斯氏菌在辽河口沉积物中除了参与氮循环的反硝化过程，也极有可能参与了微生物固氮的过程，并在两种过程中发挥重要作用，对该菌属下的细菌还有待更深入的研究。海洋嗜酸杆菌在南海[24]、丹麦港口[25]等海域中也被检测出，它可以利用环境中的铵盐，增加生物体内的氮积累，并分泌某种物质从而保证自己免受其他微生物的捕食，其生命活动与海洋及河口中还原性氮的分布密切相关[26]。

4.3 环境因子对反硝化细菌的影响

环境因子与反硝化细菌群落的相关性表明，细菌群落对环境因子的响应不是一成不变的，而是多种复杂因素共同作用的结果。从张慧珍等[11]的研究中可以看出，在B1 采样点沉积物主要为黏土和粉砂，其粒径较小，而C4 采样点沉积物中含有84.88%的砂粒，有文献指出，按照沉积物粒径划分的黏土、粉砂、砂在沉积物中所占的比例可以影响基因的拷贝数，黏土和粉砂中的细菌丰度更高[27–28]，较高的含砂比是造成采样站位功能基因丰度较低的主要原因。

nirK 型反硝化功能基因的丰度指数和多样性指数表明，nirK 型功能基因在淡水沉积物中具有更高的OTU 数目和更为丰富的物种多样性。有研究指出，反硝化细菌的多样性及其分布受盐度的影响，在盐度发生变化时，反硝化细菌表现出随盐度波动的抗性或阻力[29]。我们的研究结果也进一步证实了盐度在nirK 型反硝化细菌群落组成和多样性中的重要性。Liu 等[30]在对长江口沉积物的研究中指出水的理化性质和养分组成是预测河口沉积物中反硝化作用和N2O 生成速率的主要影响要素。同时，Zhou 等[31]的研究提出温度是影响细菌反硝化群落结构变化的关键环境因子之一。nirK 型反硝化菌群主导着反硝化过程中 NO−2向NO 的转化，其中 NO−2是反硝化过程中的电子受体，并且 NO−3与 NO−2浓度的比值影响反硝化菌群的多样性[32]。以上观点也佐证了盐度、DO 以及NO−2含量与nirK 型反硝化微生物群落结构具有较大影响的观点。

本研究探究了在辽河口沉积物反硝化作用中主要的功能基因和反硝化优势菌，同时分析了环境因子对于功能基因丰度和细菌群落结构的影响，旨在为河口环境富营养化中氮元素的循环提供理论依据。

5 结论

(1) 功能基因的丰度在不同反硝化过程中存在差异。narG 基因占据最高基因丰度，而norB 基因随河口向海延伸而丰度下降，且与多数环境因子具有显著相关性。功能基因narG、nirK 和nosZ 之间具有显著相关性。

(2) 辽河口的盐度梯度对nirK 型反硝化细菌群落的影响显示出：高盐度组中主要作用菌属为产碱杆菌属、德沃斯氏菌属、假单胞菌属、根瘤菌属等；近岸低盐度组中主要作用菌属为假单胞菌属、Defluviimonas、螯台球菌属、中慢生根瘤菌属、胺杆菌属等；在淡水−咸水的混合交界组中没有筛选出有差别作用的重要菌属。

(3) 辽河口沉积物中反硝化细菌的丰度主要受沉积物粒径（泥砂比）的影响。nirK 型反硝化细菌多样性随河口向海延伸逐渐下降，沉积物的盐度、pH、DO 以及 NO−2含量是影响nirK 型反硝化微生物群落结构的主要因素。