太湖上游丘陵区典型水塘反硝化潜力变化特征及其影响因素*

2022-11-09王慧亮张汪寿李恒鹏夏天宇陈永娟

湖泊科学 2022年6期

王慧亮，曹恒,，张汪寿，李恒鹏，何鹏，夏天宇,，陈永娟,3**

(1：郑州大学水利科学与工程学院，郑州 450001) (2：中国科学院南京地理与湖泊研究所，中国科学院流域地理学重点实验室，南京 210008) (3：滁州学院土木与建筑工程学院，滁州 239099)

随着人口的增长和粮食需求的增加，在追求高产量的同时，更多的农业化肥被投入使用. 由于其利用率较低，全球超过一半的施用化肥通过水塘等汇入河流、湖泊和海洋，并造成严重的生态环境问题[1]. 丘陵区的农业发展模式主要以“坡地-沟塘”为主[2]，近年来，丘陵区的土地开发呈现快速增长的趋势[3]，其中坡地种植以茶树为主，施氮量可达959 kg/(hm2·a)，约是水田的2倍[4]，造成了大量的氮素流失，严重威胁水质安全，并导致水体富营养化、生物多样性降低等一系列环境问题[5]. 水塘作为氮汇的重要场所，同时也多分布在流域上游，因此研究水塘反硝化作用对于流域源头的水质保护具有重要意义.

丘陵区水塘是一种常见的连接农业区与河流的重要水文通道[6]，其能够截留、过滤和净化大量的污染物，是氮磷等污染物发生迁移转化的重要场所，不仅保障了农业用水，也促进了养分的消纳和循环利用[7]. 水塘中氮循环过程主要包括：氨化作用、矿化作用、植物吸收、硝化作用和反硝化作用等，其中氮素去除的机制主要有吸附沉降、植物吸收、反硝化作用，前两者都难以从根本上将氮素去除，而后者使活性氮最终以惰性氮(N2)形式回归大气，将氮素永久去除[8]. 反硝化是氮循环的重要过程，也是水塘氮素去除的主要机制[9]，对于减少氮素输出、改善丘陵区水质具有重要作用.

反硝化作用受诸多因素的影响，如温度、pH、硝态氮浓度、有机碳和溶解氧(DO)等[1]. 由于水塘水体受污染程度和污染类型不同，导致反硝化作用的主要影响因素也存在差异. 反硝化作用主要发生在淹水环境中(如稻田、水塘等)，是氮素有效去除的最佳途径[10]. 然而，反硝化具有强烈的时空异质性，精确、简便测定反硝化速率一直是一个世界性难题[11]. 这是由于反硝化主要产物N2在大气中的背景值高达79%，在如此高的背景环境中一般难以直接测定[12]. 当前对原位状态下水塘反硝化过程的研究报道较少. 已经建立的反硝化测定方法多以间接测定为主，如乙炔抑制法、15N同位素示踪法、化学计量法等[9]，但是这些方法普遍存在人为扰动大、操作繁琐、误差大等不足，无法精确测定水土界面反硝化速率. 近年来，基于N2∶Ar的膜进样质谱法(MIMS)被广泛用于直接测定淹水环境反硝化速率[13]，其测定速度快(<3 min)、精度高(<0.03%)、进样少(<7 mL)、操作简便等，被广泛用于不同水体环境反硝化作用研究[9]. 目前，MIMS不仅已大量应用于河道[14]、河口以及稻田土壤反硝化速率直接测定[14-15]，还可以与反硝化培养装置串联，模拟室内环境条件控制下反硝化速率的变化[16-17].

天目湖位于太湖流域上游，由沙河和大溪两大水库组成，是溧阳市重要的饮用水源地. 天目湖流域内分布着1000多个大小不一的水塘，每1 km2约有7个[18]. 这些水塘主要由降雨冲刷和坡面产流汇集而成，一般多具有小流域结构. 由于其数量众多，这些水塘构成了丘陵坡地污染拦截和消纳的第一道防线，其反硝化脱氮潜力大小将直接影响下游天目湖湖体输入的氮量及水生态健康. 本文根据天目湖水塘集水区土地利用及水塘类型的特点，选择天目湖流域内4类(茶园塘、村塘、养殖塘、林塘)共14个典型水塘，采用膜进样质谱法(MIMS)测定水体中溶存的氮气(N2)浓度，并结合水气交换通量模型估算水体反硝化潜力，分析水塘反硝化潜力变化特征及其主要影响因素，揭示其时空变化规律，以期为太湖上游丘陵区氮污染源头控制和水质改善提供科学支撑.

1 材料与方法

1.1 研究区域概况

天目湖流域位于江苏溧阳境内，为太湖上游低山丘陵区，流域总面积246 km2，由沙河和大溪两大水库流域组成，面积分别为154和92 km2. 根据2017年天目湖流域土地利用遥感数据统计分析，天目湖流域土地类型主要为：林地、耕地、茶园、水域等，其中林地是该流域最主要的土地类型，面积约110 km2，占流域总面积44.6%；耕地面积仅次于林地，约为47 km2，占比18.9%；茶园面积29 km2，位居第三，占比11.9%. 天目湖流域属亚热带季风气候，四季分明，夏季和秋季温度高，降水量大. 流域内农业种植包括茶树、水稻、油菜、小麦、稻-麦和稻-油，氮肥施用量分别约为959、445、315、247、692和760 kg/(hm2·a)[4]，农业施肥也主要集中在夏季和秋季，受气候条件和施肥的影响，反硝化作用在夏季和秋季较为剧烈.

1.2 样品采集与沉积物培养实验

本文选择天目湖流域内14个水塘，根据其集水区土地利用类型及用途不同划分为4种类型水塘，其中茶园塘3个、村塘4个、林塘4个、养殖塘3个(图1). 在反硝化作用剧烈的夏季和秋季，即2020年6-11月期间，每月月末采集1次各水塘表层水样，测定各水塘水体中溶存的N2浓度. 同时，于2020年10月采集1次各水塘的柱状沉积物，用于室内反硝化的测定. 其中，由于林塘沉积物厚度较薄，无法采集完整柱状沉积物，因此林塘的反硝化速率未进行测定.

图1 天目湖流域采样点位分布及不同类型水塘Fig.1 Distribution of sampling points and diagram of different types of ponds in Lake Tianmu Basin

采集样品用于N2测定时，应全程动作缓慢，避免气泡产生，并减少与空气接触. 首先，将采样器缓慢没入水塘水体表面约10 cm后，缓慢提起采样器. 打开采样器底部止水阀，将采样器中的水样通过硅胶管缓慢引流至带橡胶垫片的螺口瓶(7 mL)的底部，当水样装满整个样品瓶后应继续引流，让其溢出样品瓶体积3倍以上. 随后，在样品瓶液面下注射0.1 mL饱和ZnCl2溶液，以抑制微生物的活动. 最后，迅速拧紧瓶盖，并将样品瓶倒置观察其中是否有气泡产生，如果有气泡则需要重新采集，每次样品均取3个重复. 同时采集水塘表层水样约550 mL，用于测定水体中各项水质指标.

在样品采集时的采样现场，利用便携式多参数水质检测仪(YSI 6000，USA)对采样点位的各项环境要素进行现场测定，包括水温、pH、溶解氧等指标.

此外，还在2020年10月采集水塘无扰动沉积物柱状(10 cm)样品，在各水塘不同区域采集3个平行沉积物柱状样品，用于室内反硝化培养实验，采样结束后并采集各水塘原位水样作为反硝化培养的上覆水. 室内培养时，首先将上覆水倒入培养装置中，并将柱状沉积物垂直放入培养装置中，并控制上覆水水位高出培养柱约3～5 cm，调节温度与采样时的水温一致，同时利用通气泵保持上覆水通气6 h以上，使得沉积物达到稳定条件，并使得柱样淹水环境溶解氧浓度与原位水体基本一致. 培养开始前，首先关闭通气阀，将柱盖在水下拧紧，并连接柱盖上的进出水管；然后打开进出水管的阀门，当进水管的原位水样可在重力作用下进入培养柱，同时出水管有水样流出时，说明整个培养系统是密闭的. 同时，要观察进出水管和密闭柱样是否有气泡残留，若有气泡则应排空气泡. 培养开始后，开启培养装置顶部的可调速电机，转速设置为75转/min，以保证沉积物培养柱中水体的充分混匀，采样第一个样品作为0时刻样品，同时依次在2、4、6、8、10、12、22、24 h时刻分别取样，每个时空均采集3个5 mL的平行样品，以减少误差，采集的样品均于48h内测定完毕. 具体的室内培养实验方法详见李晓波等[9].

1.3 理化指标的室内测定

对于N2浓度，利用膜进样质谱仪(MIMS)测定水样(水塘表层水和室内沉积物培养上覆水)中N2/Ar的比值来计算其氮气浓度.

1.4 反硝化速率计算

水塘反硝化作用的估算包括2种方法，一是对水体的直接测定，采用MIMS测定水体N2浓度，计算N2的过饱和浓度，并结合水气交换通量模型计算水体反硝化速率；二是沉积物室内培养测定的反硝化速率，根据不同培养时间反硝化产物N2浓度和时间变化进行线性拟合，计算反硝化速率. 第一种方法操作简便，对原位条件无扰动，能够快速获取野外原位条件下反硝化速率. 但同时存在较大的不确定性，这主要是由于当前有多个备选的水气交换通量估算模型，当这些模型应用至丘陵区沟塘反硝化速率测算中可能存在误差不一的问题；第二种方法是目前最为常用的方法，但一般在实验室控制条件下进行，操作繁琐且不可避免地改变了原位条件(比如光照、风速扰动等)，对估算精度会产生一定的影响. 尽管如此，众多学者认为该方法估算的反硝化速率相对可靠，与真实的反硝化速率较为接近[8,13-14,17]. 本研究根据沉积物测算的反硝化速率来约束模型估算的结果，并优选适合于丘陵沟塘的水气交换通量模型，来进一步提升反硝化测算的可靠性和精度.

上述2种方法N2浓度的测定均采用MIMS法. MIMS法是利用气体过滤膜系统直接将水样中溶存的气体吸入质谱仪真空分析室，经过离子化，由质谱仪测出气体组成和含量，对N2、Ar的测定精度分别小于0.5%、0.05%，但对于N2/Ar的测定精度小于0.03%，水环境中Ar溶解度一般只受温度和盐度的影响，其浓度非常稳定.

通过水气交换通量模型来计算水体反硝化速率时，首先需计算出过饱和浓度(Δ[N2]，μmol/L)，具体可通过下述公式得出：

Δ[N2]=[N2]/[Ar]×[Ar]*-[N2]*

(1)

式中，[N2]/[Ar]为经过质谱校准过的水样中N2和Ar的浓度比值，[N2]*、[Ar]*为特定温盐条件下N2、Ar的理论平衡浓度(μmol/L).

依据水-气界面的分子扩散模型和亨利定律，估算不同单位面积水塘单位时间反硝化N2的增量(F，mmol/(m2·d))，计算公式如下：

F=k·Δ[N2]

(2)

式中，k为气体扩散系数，由现场测定的温度、风速等参数结合气体Sc数(schmidt number)计算. 考虑到风速变化对估算k值的影响，本研究采用Cole和Caraco提出的CC98模型[19]，该模型开发时的风力条件与本研究的水塘最为相似，估算k值的计算公式如下：

CC98:k=[2.07+(0.215U101.7)](Sc/600)-2/3

(3)

式中，U10为水面上方10 m处的风速(m/s)，由流域内安装的气象站获得. 而Sc数为水的动力黏度与待测气体分子扩散速率之比，对于特定气体，Sc数与水温、盐度等物理参数有关，Wanninkhof[20]提出N2气体Sc数与水温的关系式如下：

Sc=1970.7-131.45t+4.1390t2-0.052106t3

(4)

式中，t为水温(℃).

除了CC98模型中的k值计算方法外，还有许多其他的方法来计算k值. 本文列举了另外2种计算k值的经验模型，将不同模型计算的结果进行比较，以探究不同k值计算方法的选择所带来的不确定性. 模型方程如下：

1)LM86模型: Liss and Merlivat (1986)[21]：

k=0.17U10(Sc/600)-2/3(0k=(2.85U10-9.65)(Sc/600)-1/2(3.6k=(5.9U10-49.3)(Sc/600)-1/2(U10>13)

(5)

2)CW03模型: Crusius and Wanninkhof (2003)[22]：

(6)

在室内沉积物培养实验中，根据反硝化产生的N2浓度变化与培养时间的关系来估算反硝化速率. 具体计算公式如下：

Y=mX+b

(7)

式中，Y是不同培养时间对应的反硝化产物N2浓度(μmol/L)；斜率m是沉积物上覆水中N2浓度变化速率(μmol/(L·h))；X是培养时间(h)；b是初始上覆水中N2浓度(μmol/L)；结合上覆水体积V(L)和沉积物柱样横截面积S(m2)，即得反硝化速率(F，mmol/(m2·d))，计算公式如下：

F=0.024mV/S

(8)

式中，0.024为单位换算系数.

2 结果与分析

2.1 上覆水与沉积物理化特性

表1 夏季和秋季各个水塘水体水质参数(平均值±标准偏差)*

表2 各水塘沉积物理化参数(平均值±标准偏差)*

2.2 水塘水体N2过饱和浓度变化

用膜进样质谱法测定水塘水体反硝化作用，结果表明，不同类型水塘反硝化作用存在一定的时空差异. 图2a为不同类型水塘夏季和秋季N2过饱和浓度(ΔN2)的变化. 14个水塘上覆水溶解N2过饱和浓度介于1.36～28.35 μmol/L之间，平均值为(8.23±6.04) μmol/L. 夏季4种类型水塘N2过饱和浓度平均为(8.81±4.08) μmol/L，其中茶园塘N2过饱和浓度在所有类型水塘中最高，达到(13.91±4.79) μmol/L；其次是村塘，为(8.28±1.36) μmol/L；林塘和养殖塘较为接近，分别为(6.54±1.34)、(6.51±1.79) μmol/L. 夏季水塘上覆水溶解N2过饱和浓度大小依次为：茶园塘>村塘>林塘>养殖塘.

秋季4种类型水塘N2过饱和浓度均值为(7.64±7.46) μmol/L，其中茶园塘ΔN2最高，为(19.39±5.48)μmol/L；养殖塘次之，为(4.61±1.57) μmol/L；村塘为(4.57±1.29) μmol/L；林塘ΔN2最低，为(2.00±0.42)μmol/L，秋季水塘上覆水ΔN2大小依次为：茶园塘>养殖塘>村塘>林塘. 总体上，水塘水体夏季ΔN2高于秋季，村塘、林塘和养殖塘均呈现相同的规律，而茶园塘夏季ΔN2低于秋季，可能与秋季茶园高强度施肥有关.

2.3 水塘水体反硝化潜力

采用CC98模型中的水气交换通量模型来估算不同类型水塘水体反硝化潜力，由图2b可知水塘反硝化潜力在时间和空间上均存在差异. 14个水塘的反硝化速率均值为(4.75±3.27) mmol/(m2·d). 其中茶园塘反硝化速率最高，为(9.36±1.44) mmol/(m2·d)，村塘和养殖塘分别为(5.47±1.52)和(3.09±0.52)mmol/(m2·d)，林塘最低为(1.09±0.48)mmol/(m2·d)，茶园塘反硝化速率大约是村塘的1.7倍、林塘的8.6倍、养殖塘的3倍，其反硝化速率均值大小顺序为：茶园塘>村塘>养殖塘>林塘. 反硝化速率与氮气过饱和浓度的季节变化规律一致.

室内沉积物培养实验结果表明，不同类型水塘沉积物的反硝化速率表现出一定的差异性(图2b). 9个水塘的沉积物反硝化速率在2.21～10.34 mmol/(m2·d)之间，平均值为(5.73±2.65)mmol/(m2·d). 其中村塘和养殖塘反硝化速率分别为5.24和2.75 mmol/(m2·d)，茶园塘的反硝化脱氮能力最强，反硝化速率为9.19 mmol/(m2·d). 总体上，室内培养试验测得的反硝化速率与CC98模型估算结果较为接近.

图2 夏季和秋季不同类型水塘N2过饱和浓度(a)和水体反硝化速率(b)(b图中“-”表示水塘沉积物反硝化速率)Fig.2 The excess concentrations of N2(a) and denitrification rate in water(b)in different types of ponds in summer and autumn(“-” in the figure b represents the measured denitrification rates from pond sediments)

2.4 水塘反硝化作用与环境因子的关系

图3 水塘水体N2过饱和浓度与水体(a)和沉积物(b)理化指标的相关性分析Fig.3 Correlation analysis of N2 supersaturated concentration in pond water with physicochemical indexes of water body (a) and sediment (b)

沉积物反硝化速率与硝态氮和总氮浓度呈现显著的正相关性(P<0.01，图3b)，而与容重存在显著的负相关性(P<0.05)，pH、含水率和氨氮与沉积物反硝化速率的相关性不显著. 水体和沉积物反硝化作用均与硝态氮浓度呈显著正相关关系.

3 讨论

3.1 水塘沉积物反硝化估算及其不确定性

通过膜进样质谱法测定上覆水N2过饱和浓度，结合水气交换通量模型估算上覆水反硝化速率在0.53～12.65 mmol/(m2·d)之间，平均值为(4.75±3.27)mmol/(m2·d). 与其他研究相比，本研究报道的水塘反硝化速率范围明显高于其他不同类型水体的范围. 比如，张波等对太湖流域典型河流沉积物反硝化速率研究发现，利用15N同位素配对法测得25℃时反硝化速率为0.79 mmol/(m2·d)，10℃时为0.26 mmol/(m2·d)，明显低于水塘沉积物反硝化速率[23]. Bonnett等利用乙炔抑制法测定英国西南部Tamar河流沉积物反硝化速率为0.64 mmol/(m2·d)，低于本研究水塘反硝化速率[24]. Lazar等利用15N同位素配对技术测得美国东北部农村海狸池塘沉积物反硝化速率为1.73～4.21 mmol/(m2·d)，也总体低于本研究的结果[25]. 然而，一些河流的反硝化速率要高于本研究的结果. 比如，陈能汪等利用膜进样质谱仪对九龙江河口区溶解N2浓度测定，并采用水气交换模型估算了其反硝化速率，发现其均值为21.7 mmol/(m2·d)，明显高于本研究水塘反硝化速率[15]. 龙虹竹等利用改进的乙炔抑制法对川中丘陵区农业源头沟渠沉积物反硝化速率测定，其反硝化速率均值为8.96 mmol/(m2·d)，高于本研究水塘沉积物反硝化速率[26]. 李佳霖等同样利用15N同位素配对法测得弥河沉积物反硝化速率为16.34 mmol/(m2·d)，明显高于天目湖流域水塘沉积物反硝化速率[27]. 通过这些对比发现，本研究报道的水塘反硝化速率波动范围在全球不同类型的水体中均有报道. 由于丘陵区水塘与坡地相连并直接受纳坡地产生的污染物，本研究的结果证实了其相对较强的反硝化速率在削减氮负荷及减少下游污染方面发挥了重要的作用.

本研究通过对水样中氮气的直接测定，并通过水气交换模型测算了反硝化速率. 与其他传统的测定方法相比，膜进样质谱法直接测定水体反硝化作用具有精度高、操作简便、结果可靠等特点[28-30]，并且适用于无沉积物或少沉积物的水体. 本研究对不同类型水塘水体反硝化速率直接测定，同时采集各水塘沉积物柱样进行室内培养测定反硝化速率，两种测定方法的结果具有可比性(图4)，并且前者实验周期更短.

然而，本研究的估算存在一定的不确定性. 首先，反硝化估算受采样频率的影响. 例如每月采集一次水塘上覆水存在很大的不确定性，因为N2浓度时刻在变化，低频的观测难以准确获得N2的动态变化. 此外，本研究仅估算了N2水气通量，未考虑不完全反硝化产生的N2O排放. 尽管对于大多数水体而言，N2O水气通量一般只占总通量的0.1%～1.0%，比如九龙江河口区反硝化作用产生的N2O水气交换通量只占0.25%[15]. 但在某些情况下，最高可达到6%[31]，这说明我们在一定程度上低估了水塘反硝化潜力，需要对N2和N2O水气通量同时评估. 此外，反硝化估算是否准确与水气交换模型的选择有很大的关系，这是因为当前有很多不同的估算模型(见本文2.4节部分). 其中，Cole和Caraco[19]提出的水气交换模型被应用最广泛、精度最高，且适用于丘陵区水文气象条件. 该模型采用风速和温度作为模型的输入，使其更接近原位环境条件，能更精准地刻画水气交换过程. 本研究还进一步采用沉积物柱状培养测定反硝化速率来验证水气交换通量的结果，研究进一步证实CC98模型估算的反硝化速率与传统室内反硝化模拟试验测得的速率非常接近. 然而，对于其他2个常用的水气交换模型，估算的结果远小于室内沉积物培养测得的反硝化速率. 比如，当采用LM86和CW03模型时，估算的反硝化速率仅为实验室测得的速率的1/10～1/3，这说明这些模型不适用于丘陵区水塘反硝化潜力的估算. 然而，水气交换过程不仅受到大气扰动，也受到水文和微气象等因素的影响，在未来的工作中，应原位测定k值，以期更准确地估算反硝化速率.

图4 不同k值计算模型估算的水塘反硝化速率差异(LM86、CW03和CC98指采用相应的公式计算得出的反硝化速率，具体公式见本文1.4节部分；图中“-”代表水塘沉积物反硝化速率)Fig.4 The difference of denitrification rates of ponds estimated by different k calculation models

3.2 影响水塘反硝化作用的主要因素

溶解氧水平是影响反硝化脱氮的重要因素. 这是因为反硝化一般为厌氧过程，高浓度溶解氧的存在会显著抑制反硝化功能微生物的活性，进而抑制反硝化脱氮过程[37]. 本研究发现，水塘反硝化速率与溶解氧浓度之间呈显著负相关关系(r=-0.47，P<0.05). 在不同O2浓度水平下检测土壤反硝化速率，其最高值出现在无氧条件下[38]；丘陵区水塘水体溶解氧浓度夏季低于秋季，其反硝化速率夏季大于秋季；在不同类型水塘上，茶园塘和村塘水体溶解氧浓度低于林塘和养殖塘，低溶解氧浓度有利于反硝化细菌的生长，从而促进反硝化作用[39]；在沉积物中较高的含水率会阻止氧气的传输，从而形成厌氧条件，茶园塘沉积物含水率较高，其反硝化速率也较高，二者有一定的相关关系，但相关性不显著.