ACTA Scientiarum Naturalium Universitatis Pekinensis

富营养湖泊稳态转换的恢复时间及影响因素模拟研究

吴思枫梁中耀刘永

2018-09-20 - 北京大学环境科学与工程学院, 水沙科学教育部重点实验室, 北京 100871; † 通信作者, E-mail: yongliu@pku.edu.cn

摘要为了揭示富营养湖泊恢复时间和负荷削减强度之间的非线性关系及其影响因素, 以典型湖泊为理论分析对象, 基于经典磷循环模型, 通过数值模拟, 计算得到不同负荷削减强度下的湖泊恢复时间, 即不同负荷下, 湖泊从初始富营养状态, P浓度逐渐降低回到清水稳态所需要的时间。通过计算不同参数取值下湖泊恢复时间的变化, 探究湖泊不同方面的特征对湖泊恢复时间的影响及可能的管理契机, 得到如下结果。1) 湖泊由浊水稳态向清水稳态的恢复时间和负荷削减强度之间存在非线性关系。若将外源负荷强度控制在略低于“浊水–清水”转换阈值, 则恢复时间在40年以上, 若加大削减强度, 则可缩短恢复时间; 但随着削减强度的持续增加, 其边际效应将逐渐减弱。2) 湖泊形态和状态会对恢复时间产生明显的影响。在同样的负荷削减强度下, 寒冷地区较深的湖泊恢复时间更短; 沉积物释放较快的湖泊恢复时间更长; 水力停留时间短的湖泊恢复时间也更短。因此, 从理论上讲, 通过生态修复工程措施来降低沉积物P的释放或改善水动力条件, 能够缩短富营养湖泊的恢复时间; 湖泊状态的改善还可降低“浊水–清水”稳态转换阈值, 进而降低治理难度。关键词富营养化; 稳态转换; 恢复时间; 磷循环模型中图分类号 X37

湖泊富营养化是当前全球环境领域面临的突出问题[1–2],也是典型的稳态转换现象[3–4]。湖泊存在两种稳态状态: 清水稳态和浊水稳态。在清水稳态下, 湖泊藻类生物量低, 沉水植被覆盖率、透明度和溶解氧含量高; 在浊水稳态下, 湖泊藻类生物量高, 水质差[3,5–6]。导致湖泊发生稳态转换的驱动因素较多, 其中最重要、最常见的驱动因素是过量营

[7–9]养物质的输入, 尤其是磷的输入。由于生态系统内部复杂的动力学过程, 湖泊在不同稳态下都存在负反馈(stabilizing feedback)机制, 使其能够承受一定的外部压力, 维持当前的稳态而不发生明显变化[3,10]。对于清水稳态的湖泊, 主要的负反馈机制包括沉积物的缓冲作用、沉水植被对沉积物释放的抑制、浮游动物对浮游植物的捕食压力、硅藻在低营养盐浓度下的竞争优势等[7,11–12]。当环境压力逐渐增加, 湖泊内部的负反馈作用逐渐减弱, 直到环境压力超过一定的阈值(图1中T1)时, 湖泊变为正反馈(destabilizing feedback)的系统。此时湖泊沉水植被减少, 沉积物释放迅速增加, 内源和外源负荷之和大于流出和沉降之和, 水体P浓度迅速增加。P浓度的增加进一步促进藻类生长, 水体浊度增加, 底层光照减少, 加剧沉水植被的丧失。在这一正反馈作用机制下, 湖泊迅速偏离原有状态, 进入浊水稳态[3,12–13]。

当湖泊进入浊水稳态后, 系统演变呈现另一套不同的负反馈机制。例如, 浊水稳态下湖泊溶解氧含量低, 沉积物释放量大; 水体浊度高, 水底光照条件差, 沉水植被消失, 加剧沉积物再悬浮和释放;营养盐浓度高, 蓝藻由于对光照竞争的优势而成为优势种, 进一步降低水体透明度[5,11–12]。由于这些负反馈机制的存在, 使得湖泊长期维持在浊水稳态[13]。由于湖泊在不同稳态下的负反馈机制不同,因此退化和恢复过程的变化轨迹并不重合(hysteresis),导致要将湖泊扭转回到清水稳态, 必须将外部压力减小到比发生稳态转换前的负荷(图1中T1)更低的状态(低于T2)[3,14], 且需要维持长时间的削减。发生稳态转换的本质是湖泊内部反馈机制的重建, 打破维持现存稳态的负反馈机制, 使得湖泊从浊水稳态回到清水稳态。因此, 稳态转换的前提是生态系统内部结构的变化, 而该变化通常需要很长的时间,一般在几年到几十年之间[3,15]。在恢复过程中, 在稳态转换之前, 可能有较长时间湖泊对外源负荷削减的响应并不明显, 使得水质改善效果缓慢。以往对湖泊稳态转换的研究大多关注稳态转换中的反馈机制和过程及其驱动因素的识别、阈值和预警[13,16–22],而对恢复时间的讨论相对较少。本文定义湖泊由浊水稳态向清水稳态转换所需要的时间为富营养湖泊生态系统的恢复时间。对于富营养湖泊的治理, 驱动因素、阈值和预警都有重要的管理意义, 然而, 这些都不能回答湖泊恢复所需要的时间。理论上, 将外源负荷控制在阈值(T2)以下, 系统均实现稳态转换, 并回到清水稳态。同一湖泊在不同负荷削减强度下的恢复时间不一样; 不同类型的湖泊需要的恢复时间也不相同。恢复时间也是需要考虑的一个重要问题: 负荷削减强度和维持的时间影响治理所需要的经济投入, 但对于一些严重污染的湖泊, 较短的恢复时间可能比较低的经济成本更重要。所以, 对恢复时间进行定量并探究其影响因素具有重要意义。本研究基于经典的湖泊磷循环模型, 探讨稳态转换下富营养湖泊的恢复时间及其影响因素。该模型被广泛应用于磷限制湖泊稳态转换阈值和预警的研究中; 也有一些复杂的水质–水生态模型(如Pclake), 将磷循环模型耦合入复杂模型中以提高准确度[15,23–25]。

1模型构建1.1模型方程及参数

本文选取的磷循环模型方程为

其中, P为水体中p的浓度(g/(m2·a)); lp为p的外源负荷(g/(m2·a)); s为p的损失速率(a1), 包括沉降和流出

速率, 与水体P浓度成正比; r为P的最大释放速率

释放的半饱和常数(g/m2), 即沉积物P释放速率达到1/2最大速率时对应的水体P浓度; q是无量纲常数,反映沉积物释放随水体P浓度变化的快慢。

模型中考虑P的外源负荷、损失(流出和沉降)及沉积物释放这些主要过程。模型假设: 对于分层湖泊而言, 沉积物P释放主要为溶解态P的扩散, 而底层缺氧状况决定溶解态P的浓度; 对于P限制湖泊, 水体P浓度决定初级生产力, 进而决定有机物沉降和底层缺氧状况, 且底层缺氧与水体P浓度之间的关系呈“S”型曲线(sigmoid function), 因此沉积

的适用范围为存在夏季热分层的P限制湖泊。本文的参数取值(表1)参考符合模型假设且被广泛研究的Lake Mendota, 该湖泊平均深度为12.5 m[15,26]。1.2 恢复时间及影响因素分析方法为了获取不同入湖负荷下湖泊 P浓度的变化过程, 需遵循如下步骤: 1)对不同负荷下的湖泊源汇关系进行分析, 其中, 源包括外源负荷和沉积物释

衡点附近 P 浓度的变化趋势, 判断平衡点是否为稳定平衡点, 若平衡点受到干扰偏离之后有回到平衡点的趋势, 则为稳定平衡点; 3) 通过求解不同外源负荷下式(1)的平衡点, 得到湖泊稳态转换曲线, 并定量地得到外源负荷阈值 T1和 T2以及阈值负荷下对应的 P 浓度P 和P 。T1 T2由于本研究关注的是富营养湖泊的恢复过程,因此起始P浓度应该设定在富营养状态, 即大于P。T2

本研究将起始P浓度设定为0.4 mg/l (劣V类), 为便于分析, 湖泊深度为12.5 m, 将体积浓度转换为以面积计算的P浓度(面积浓度)为5 g/m2。根据前面的界定, 恢复时间可细化为: 在一定的负荷强度削减下, P浓度从特定富营养稳态降低到P 所需的时

T1间。初始负荷为该负荷下稳态浓度等于初始P浓度

荷 lp , 如式(2)所示。本文假设起始P浓度为5 g/m2,

0因此起始外源负荷为2.265 g/(m2·a)。削减比例为负荷削减前后的差值与削减前负荷的比值。

通过对式(1)积分, 得到P浓度随时间的变化曲线, 以及不同负荷下的恢复时间。与湖泊类型和形态相关的参数lp, s, r和q会对恢复时间产生影响。本文通过改变参数取值, 得到不同参数取值下恢复时间与负荷削减之间关系的变化, 进而分析得到湖泊类型和形态对恢复时间的影响。本研究的模拟通过R语言[27]实现, 使用的工具包主要有desolve (获取微分方程的时间变化曲线)[28]、rootsolve (求解微分方程平衡点)[29]和phaser (相图)[30]。

2研究结果与讨论2.1稳态转换过程及阈值

由模拟结果可知, 当外源负荷较低时, 湖泊生态系统只存在1个稳定平衡点, 为清水稳态; 当外源负荷较高时, 湖泊也只存在1个稳定平衡点, 为浊水稳态; 当外源负荷适中时, 湖泊可能存在3个平衡点, 其中2个为稳定平衡点, 中间为不稳定平衡点(图2(a))。在湖泊只存在1个稳定平衡点的负荷条件下, 无论湖泊初始状态如何, 最终都逐渐趋于稳定平衡点(图2(b)和(d))。在多稳态情况下, 湖泊最终状态与起始浓度有关, 若起始P浓度低于不稳定平衡点P浓度, 则最终趋向于清水稳态; 若高于不稳定平衡点P浓度, 则最终趋向于浊水稳态。计算得到不同负荷输入条件下的平衡点P浓度(图3), 湖泊存在多稳态的外源负荷范围为0.825~ 1.35 g/(m2·a), 由浊水稳态变为清水稳态的阈值(T2)为0.825 g/(m2·a), 由清水稳态变为浊水稳态的阈值(T1)为1.35 g/(m2·a)。

2.2 不同流域负荷削减强度下湖泊的恢复时间

由于湖泊初始的P面积浓度为5 g/m2, 为达到清水稳态, 外源负荷必须低于T 2, 即0.825 g/(m2·a)。为了得到不同流域负荷削减强度下的恢复时间, 在lp不同取值下对式(1)积分, 得到P浓度随时间的变化曲线(图4)。可以看出, 若负荷削减至稍低于T2,湖泊的恢复时间将会很长, 此时加大削减力度对缩短恢复时间有明显效果。例如, 当外源负荷削减至0.82 g/(m2·a)时, 恢复时间为46年, 而将外源负荷降低到0.80 g/(m2·a)时, 恢复时间则为21年。继续降低外源负荷, 对缩短恢复时间的作用逐渐减弱。例如,将削减后的负荷从0.80 g/(m2·a)降低到0.10 g/(m2·a),恢复时间仅缩短约17年, 效果不如将负荷从0.82 g/(m2·a)降低到0.80 g/(m2·a)。综上所述, 可知恢复时间与流域负荷削减以及削减比例之间具有非线性关系(图5)。

2.3 湖泊类型和形态相关参数对恢复时间的影响

模型参数沉积释放半饱和常数m、指数q、最大沉积物释放速率r和p损失(沉降和流出)速率s与湖泊类型和形态相关。通过模拟不同参数取值下湖泊恢复时间的变化, 可以得到不同湖泊类型恢复时间的差异, 参数取值范围都在存在多稳态解的范围内。图6显示, 参数m, q, r和s均对恢复时间有明显的影响; m和s值越大, 在同等削减比例下恢复时间越短; q和r值越小, 在同等削减比例下恢复时间

越短。q值主要与湖泊深浅和温度有关, 一般较浅的、温暖的湖泊q值较大, 而较深的、寒冷的湖泊q值较小; m值主要与湖泊均温层深度、厚度和温度有关[26],而均温层的深度和厚度与湖泊深度直接相关[6]。一般而言, q值较大的湖泊m值相对较小。q值越小, m值越大, 恢复时间越短。这说明在同样的负荷削减条件下, 较深的、寒冷地区的湖泊相对容易治理。r值越小, 恢复时间越短。r值与沉积物表面性质和沉积物P浓度有关。对于分层湖泊, 沉积物P释放过程的决速步骤一般是间隙水中溶解态P的浓度梯度扩散[5–6,26], 所以沉积物中可溶性P浓度越低,间隙水中溶解态P浓度也越低, r值也越小, 湖泊恢复时间越短。s值越大, 恢复时间越短。损失速率s包括流出和沉降。一般而言, 流出速率越大, s值越大, 恢复时间越短; 湖泊受风力扰动越小, 净沉降速率越大, s值也越大,恢复时间越短[26,31]。

2.4 讨论

富营养湖泊恢复时间和负荷削减强度之间具有非线性关系。当水体P浓度在P 附近时, 变化速率T2明显变慢, 而在阈值P 附近时则明显变快, 且负荷T1削减比例越低, 这一变化越明显(图4)。模型分析发现, 在湖泊恢复、水体P浓度逐渐接近P 的过程中, T2浓度降低的速率逐渐减缓, 越过 P 之后, P浓度降T2低速率迅速加快, 直到越过 P , 逐渐接近清水稳T 1态, P浓度保持相对稳定。进一步的解释是, 在P浓度从阈值P 过渡到 P 的过程中, 维持浊水稳态的T2 T1负反馈机制逐渐变为正反馈, 湖泊迅速偏离原有状态而进入清水稳态[5,10]。模型结果还显示, 湖泊类型和形态对恢复时间有很大的影响。如果人为地改变湖泊的某些特性, 则可以加快富营养湖泊的恢复。例如常见的疏浚, 将表层P浓度较高的沉积物

移除, 使得沉积物的释放速率降低, 进而缩短恢复时间; 调水冲刷改善湖泊内水动力条件, 改善底层缺氧状况, 沉积物释放的半饱和常数升高, 同时使得流出速率增加, 恢复时间缩短; 通过生物操纵增加浮游植物的被捕食压力, 使得浮游植物对P的响应降低, 沉积物释放半饱和常数升高, 恢复时间缩

[5,7,15,32]短。这些措施不仅能缩短湖泊的恢复时间,而且由于湖泊状态得到改善, 使得湖泊发生稳态转换回到清水稳态的阈值降低, 需要的最小负荷削减比例也相应降低(图6)。湖泊恢复过程中发生稳态转换需要较长时间,是因为在生态系统存在一些负反馈机制[3,13]。负荷削减作为长期干扰输入, 使得湖泊状态逐渐变化,系统弹性逐渐减弱, 直到负荷输入低于一定阈值(T1), 变为湖泊生态系统单稳态, 湖泊状态恢复到清水稳态。然而, 实际情况下, 除负荷变化外, 来自湖泊内部和外部的变化与干扰(pulse)可能引起湖泊状态的迅速变化, 使得湖泊发生稳态转换, 例如干旱、飓风、鱼类养殖等[23,33–35]。由于这些干扰的存在, 即使在外源负荷没有控制在低于 T1的情况下, 湖泊也可能发生稳态转换。通常, 这些干扰及其效果都较难控制和预测。并且, 湖泊即使变为清水稳态, 系统弹性仍然较低, 会面临较大的再次恶化的风险[10,35–36]。所以, 在短期情况下, 由于一些干扰的存在可能会使富营养湖泊提前恢复到清水稳态。然而, 从长期系统稳定的角度来看, 负荷削减对于湖泊恢复和维持至关重要[23,26,37]。湖泊的生态系统状态变化是一个长期而复杂的过程, 本研究的重点在于以 P为例来阐释该过程及其时间尺度上的效应。为了反映长时间尺度的湖泊演化, 本文选用较简单但经典的模型, 未考虑外部环境变化对湖泊的影响。本研究的主旨在于通过稳态转换模型, 对富营养化湖泊管理中出现的一些问题(如治理时间较长、湖泊对负荷削减响应不明显等)给出可能的理论解释和计算, 可以有助于公众和决策者对湖泊的水质改善有较为合理的预期。本文对恢复时间的界定是以模型为基础的, 在复杂的湖泊边界条件下, 尚需进一步细化和完善。

3 结论

本研究构建了以 P 为核心的模型, 得到如下主要结论。

1) 对于本研究案例的湖泊, 从浊水稳态回到清水稳态的外源负荷阈值(0.825 g/(m2·a))低于富营养过程的阈值(1.35 g/(m2·a))。由于这一现象, 使得富营养湖泊的恢复时间通常较长, 甚至长达数十年。

2) 富营养化湖泊的恢复时间与负荷削减之间存在明显的非线性关系, 削减强度增加对缩短恢复时间的效果逐渐减弱。例如对于本研究案例的湖泊, 将削减后负荷从0.82降低到0.80 g/(m2·a), 恢复时间缩短25年; 将削减后负荷从0.80降低到0.10 g/(m2·a), 恢复时间仅缩短约17年。非线性关系主要由湖泊内部生态系统重建所需要的时间以及重建过程中 P浓度变化速率较缓导致。

3) 湖泊形态和状态也对恢复时间产生明显的影响; q 值越小, m 值越大, 相同削减强度下恢复时间越短, 说明寒冷地区较浅的湖泊恢复时间较短。r 值越小, 相同削减强度下恢复时间越短, 说明沉积物性质和沉积物释放对湖泊恢复具有意义。l 值越大, 恢复时间越短, 说明水动力条件对湖泊恢复有重要作用。在实践中, 通过生态修复等手段降低沉积物释放, 或通过工程措施改善水动力条件, 能够缩短富营养湖泊恢复时间。湖泊状态的改善能够降低阈值, 使得负荷削减强度要求更低, 从而降低治理的难度。 [1] Conley D J, Paerl H W, Howarth R W, et al. Controlling eutrophication: nitrogen and phosphorus. Science, 2009, 323: 1014–1015 [2] Smith V H, Tilman G D, Nekola J C. Eutrophication: impacts of excess nutrient inputs on freshwater, marine, and terrestrial ecosystems. Environmental Pollution, 1999, 100(1): 179–196 [3] Scheffer M, Carpenter S R. Catastrophic regime shifts in ecosystems: linking theory to observation. Trends in Ecology & Evolution, 2003, 18(12): 648–656 [4] Smith V H, Schindler D W. Eutrophication science: where do we go from here?. Trends in Ecology & Evolution, 2009, 24(4): 201–207 [5] Scheffer M. Ecology of shallow lakes. Berlin: Springer Science & Business Media, 2004

[6] Wetzel R G. Limnology: lake and river ecosystems. Oxford: Gulf Professional Publishing, 2001 [7] Genkai-kato M, Carpenter S R. Eutrophication due to phosphorus recycling in relation to lake morphometry, temperature, and macrophytes. Ecology, 2005, 86(1): 210–219 [8] Schindler D W, Hecky R E, Findlay D L, et al. Eutrophication of lakes cannot be controlled by reducing nitrogen input: results of a 37-year wholeecosystem experiment. PNAS, 2008, 105(32): 11254– 11258 [9] 孔繁翔, 高光. 大型浅水富营养化湖泊中蓝藻水华形成机理的思考. 生态学报, 2005, 25(3): 590–595 [10] Folke C, Carpenter S, Walker B, et al. Regime shifts, resilience, and biodiversity in ecosystem management. Annual Review of Ecology, Evolution, and Systematics, 2004, 35(1): 557–581 [11] Scheffer M, van Nes E H. Shallow lakes theory revisited: various alternative regimes driven by climate, nutrients, depth and lake size. Hydrobiologia, 2007, 584(1): 455–466 [12] Suding K N, Gross K L, Houseman G R. Alternative states and positive feedbacks in restoration ecology. Trends in Ecology & Evolution, 2004, 19(1): 46–53 [13] Nyström M, Norström A V, Blenckner T, et al. Confronting feedbacks of degraded marine ecosystems. Ecosystems, 2012, 15(5): 695–710

[14] ᖺ跃刚, Ᏽⱥ伟, ᮤⱥᮽ, ➼. 富营养化浅水湖泊稳态转换理论与生态恢复探讨. 环境科学研究, 2006, 19(1): 67–70 [15] Carpenter S R. Eutrophication of aquatic ecosystems: bistability and soil phosphorus. PNAS, 2005, 102(29): 10002–10005 [16] Andersen T, Carstensen J, Hernandez-garcia E, et al. Ecological thresholds and regime shifts: approaches to identification. Trends in Ecology & Evolution, 2009, 24(1): 49–57 [17] Janssen A B G, de Jager V C L, Janse J H, et al. Spatial identification of critical nutrient loads of large shallow lakes: implications for Lake Taihu (China). Water Research, 2017, 119: 276–287 [18] Janse J H, Domis L N D S, Scheffer M, et al. Critical phosphorus loading of different types of shallow lakes and the consequences for management estimated with the ecosystem model Pclake. Limnologica-ecology and Management of Inland Waters, 2008, 38(3): 203– 219

[19] ᮤ⋢↷ , ีỌ , 赵磊, 等. 浅水湖泊生态系统稳态转换的阈值判定方法. 生态学报, 2013, 33(11): 3280–3290 [20] 陈小华, 李小平, 王菲菲, 等. 苏南地区湖泊群的富营养化状态比较及指标阈值判定分析. 生态学报, 2014, 34(2): 390–399 [21] Scheffer M, Bascompte J, Brock W A, et al. Earlywarning signals for critical transitions. Nature, 2009, 461: 53–59

[22] ன⍞ᏹ , 张笑欣, ีᘐ玺, 等. 浅水湖泊稳态转换预警识别方法局限与展望. 生态学报, 2017, 37(11): 3619–3627 [23] Carpenter S R, Brock W A. Rising variance: a leading indicator of ecological transition. Ecology Letters, 2006, 9(3): 308–315 [24] Cottingham K L, Ewing H A, Greer M L, et al. Cyanobacteria as biological drivers of lake nitrogen and phosphorus cycling. Ecosphere, 2015, 6(1): 1–19 [25] Mooij W M, Domis L N D S, Janse J H. Linking species-and ecosystem-level impacts of climate change in lakes with a complex and a minimal model. Ecological Modelling, 2009, 220(21): 3011–3020 [26] Carpenter S R, Ludwig D, Brock W A. Management of eutrophication for lakes sugject to potentially irreversible change. Ecological Applications, 1999, 9(3): 751–771 [27] R: a language and environment for statistical computing. R Foundation for Statistical Computing [CP/OL]. Vienna: R Core Team, 2015 [20171101]. https://www.r-project.org/ [28] Soetaert K, Petzoldt T, Setzer R W. Package desolve: solving initial value differential equations in R. Journal of Statistical Software, 2010, 33(9): 1–25 [29] Soetaert K. Nonlinear root finding, equilibrium and steady-state analysis of ordinary differential equations [EB/OL]. (2015) [20171101]. https://rdrr.io/rforge/ root Solve/ [30] Grayling M J. phaser: an R package for phase plane analysis of autonomous ODE systems. R Journal, 2014, 6(2): 43–51 [31] Scheffer M, Rinaldi S, Gragnani A, et al. On the dominance of filamentous cyanobacteria in shallow, turbid lakes. Ecology, 1997, 78(1): 272–282 [32] Daskalov G M, Grishin A N, Rodionov S. et al. Trophic cascades triggered by overfishing reveal possible mechanisms of ecosystem regime shifts. PNAS, 2007, 104(25): 10518–10523 [33] Ratajczak Z, D’odorico P, Collins S L, et al. The interactive effects of press/pulse intensity and duration on regime shifts at multiple scales. Ecological Monographs, 2017, 87(2): 198–218 [34] Elmqvist T, Folke C, Nyström M, et al. Response diversity, ecosystem change, and resilience. Ecology and the Environment, 2003, 1(9): 488–494 [35] Carpenter S R, Arrow K J, Barrett S, et al. General resilience to cope with extreme events. Sustainability, 2012, 4(12): 3248–3259 [36] Sterner T, Troell M, Vincent J, et al. Quick fixes for the environment: part of the solution or part of the problem?. Environment: Science and Policy for Sustainable Development, 2006, 48(10): 20–27 [37] Carpenter S, Walker B, Anderies J M. From metaphor to measurement: resilience of what to what?. Ecosystems, 2001, 4(8): 765–781