水源涵养是评估生态系统功能的重要内容之一^[1], 是生态系统在特定时空条件下对水分的调蓄和保存能力, 具有拦蓄洪水、削减洪峰、调节径流和净化水质等功能^[2]。水源涵养能力对区域气候、土壤、水文和植被状况产生直接影响^[3], 在维持生态系统健康和稳定方面发挥着重要作用。水源涵养能力通常采用水源涵养量表征, 反映了一段时间内生态系统在水文循环中存储降水的能力^[4]。探明水源涵养能力的变化规律及影响因素对提高区域水资源供给能力、保持生态系统稳定具有重要意义^[5]。

水源涵养能力的定量研究是生态系统水源涵养功能研究的重点和热点^[4]。水源涵养能力的定量评估传统方法主要有水量平衡法、土壤蓄水量法、林冠截留剩余法、综合蓄水能力法和年径流法等。根据研究区特点以及不同研究目的构建的模型逐步成为水源涵养能力评价的主流方法^[6], 其中, InVEST模型, 即“生态系统服务功能与权衡交易综合评价模型”, 是目前应用较为广泛的量化多种生态系统服务功能的评估模型^[7]。通过对模型参数的修正, InVEST模型在中国不同地区生态系统服务功能评估方面有较好适应性^[8], 目前已在国内三江源地区^[9]、黑河流域^[10]、黄土高原^[11]、大别山区^[12]、北京^[13—14]等地区成功应用, 在估算水源涵养能力, 模拟生态系统水源涵养功能变化方面起到重要作用。

气候变化和人类活动是影响水源涵养能力的重要因素^[15]。气候变化主要通过降水量和蒸散量的改变对区域水源涵养能力产生影响^[16—17], 人类活动通过改变土地利用和覆被状况等地表特征直接影响着生态系统的水源涵养能力^[18]。龚诗涵等^[19]通过相关分析发现降水、温度、蒸散发与水源涵养能力存在显著正相关关系。邵全琴等^[20]对三江源生态保护和建设一期工程的生态成效进行了评估, 结果表明, 生态保护工程区水源涵养能力恢复程度优于非工程区。宁亚洲等^[21]通过量化不同土地利用类型的水源涵养量, 表明林地、草地水源涵养能力高于其它土地利用类型。情景模拟法可通过假设理想情景, 模拟单因素或多因素变化的影响。王盛等^[22]假设降水不变而潜在蒸散发生定量改变和潜在蒸散不变而降水发生定量改变两种情景, 对张承地区水源涵养功能进行评估, 结果表明, 降水对水源涵养服务的影响远大于潜在蒸散。张利利^[23]假设土地利用情况未发生变化, 模拟石羊河流域产水深度变化并分别计算土地利用对产水量的影响量, 结果表明, 土地利用与气候变化对产水深度分别起反向抑制作用和正向促进作用。

三江源是长江、黄河和澜沧江的发源地, 是国家重点水源涵养功能保护区^[24]。受全球变暖和过度放牧等人类活动的影响, 三江源地区生态状况持续退化, 水源涵养能力受到影响^[25]。2005年国家启动三江源生态保护与建设工程以来, 生态退化现象得到有效控制^[21], 水源涵养能力有所恢复^[9]。三江源生态恢复过程中水源涵养能力的变化及其驱动因素研究尚不充分。因此, 本研究基于InVEST模型, 定量评估1990—2020年三江源生态系统水源涵养能力的时空变化, 对比研究生态治理前后水源涵养能力变化及其对气候变化和人类活动的响应, 对揭示三江源区水源涵养能力变化规律和驱动机制具有较好参考价值, 为三江源地区生态保护与建设工作提供重要的参考依据。

1 数据来源与方法 1.1 研究区概况

三江源区(89°45'—102°23'E、31°39'—36°12'N)为长江、黄河和澜沧江的源头汇水区^[9], 总面积约36.3万km²(图 1)。以山地地貌为主, 海拔为2554—6826 m。气候为典型高原大陆性气候, 雨热同期, 年均温-5.4—7.5℃^[26], 气温日较差大、年较差小；年均降水量300—500 mm, 由东南向西北递减, 主要集中于6—8月^[27]。三江源区径流资源丰富, 雪山、冰川与湿地分布广泛, 素有“中华水塔”之称, 是重要的水源涵养地, 也是我国和东南亚国家生态安全和可持续发展的重要生态屏障, 在我国生态文明建设中占有重要地位。党的十八大以来, 三江源水源涵养区着力推进一系列生态保护与建设工程, 保证了生态系统质量和稳定性不断提升。

1.2 数据来源

降水量、潜在蒸散量、土地利用/覆被、高程、土壤属性和其它数据来源及相应处理如表 1所示。

1.3 研究方法 1.3.1 产水量

InVEST模型产水模块是基于水量平衡的估算方法, 将每个栅格上的降水量与实际蒸散发量之差定义为该栅格的产水量, 并求取子流域上每个栅格单元的产水量之和。计算方法如下：

(1)

(2)

(3)

(4)

(5)

(6)

式中, Y_xj表示j类土地类型栅格单元x的产水量, mm; AET_xj表示j类土地类型栅格单元x的年实际蒸散发量, mm；P_x表示栅格单元x的年降雨量, mm；R_xj为潜在蒸散发与降水量的比值；ω_x为表达气候-土壤属性关系的非物理参数；AWC_x表示植物有效利用含水量, mm；MaxsoilDepth为最大土壤深度；RootDepth为根系深度；PWAC_x为植被可利用水。sand为土壤砂粒含量, %；silt为土壤粉粒含量, %；clay为土壤粘粒含量, %；OM为土壤有机质含量, %。Z即Zhang系数, Zhang系数是表示研究区降水特征的常数, 又称“季节常数”, 能够代表区域降水分布及水文地质, 特征取值范围1—30。ET_0x表示单元x的潜在蒸散发量, mm。K_xj表示某植被类型的蒸散系数。

1.3.2 水源涵养量

计算得到产水量的基础上, 结合研究区地形指数、地表流速系数和土壤饱和导水率计算栅格尺度的水源涵养量。计算方法如下：

(7)

(8)

式中, Rentention_xj为j类土地类型栅格单元x的水源涵养量, mm；Vel_xj是j类土地类型栅格单元x的流速系数, 由模型参数表得到。K_satx为栅格单元x的土壤饱和导水率, mm/d；TI_x为栅格单元x的地形指数；drainagearea为集水区栅格数量；soil_depthx为栅格单元x的土壤深度, mm；slope_x为栅格单元x的百分比坡度。根据表 2中土地类型与生态系统类型的对应关系, 可计算得到相应生态系统类型的水源涵养量。

1.3.3 Mann-Kendall趋势分析

Mann-Kendall(MK)检验是一种非参数检验方法, 在气象、水文研究中有广泛应用^[33]。本文借助MATLAB R2018B软件实现栅格尺度水源涵养能力变化的Mann-Kendall趋势分析。公式如下：

(9)

x_k, x_i为连续的数据变量, n为数据资料时间长度。

(10)

(11)

若Z大于0, 数据序列呈上升趋势；Z小于0, 数据序列呈下降趋势。当|Z| ≥1.96时, 表明趋势变化显著；当|Z| ≥2.58时, 表明趋势变化极显著。

1.3.4 驱动因素情景分析

在假定情景下, 认为一种或多种变量保持不变, 在模型中输入该种或多种变量的固定值及其它变量的实际值进行运算, 核算出的结果即为该模拟情景的水源涵养量^[34]。InVEST模型产水量模块参数中, 土地利用/覆被参数代表相应年份土地利用/覆被状况, 降水量和潜在蒸散量参数代表相应年份气候状况。为评估气候变化和土地利用/覆被变化单一因素对水源涵养量的影响, 本研究设计了两种变化情景：一是分析气候变化对水源涵养的驱动效应, 设定气候为驱动因素, 土地利用/覆被不发生变化；二是分析土地利用/覆被变化对水源涵养的驱动效应, 设定土地利用/覆被变化为驱动因素, 气候因素不发生变化。情景一条件下, 输入的土地利用/覆被参数为1990年初始值, 输入的气候数据为相应年份的降水量、潜在蒸散量；情景二条件下, 输入的土地利用/覆被参数为1990—2020年每5年的动态数据, 降水量、潜在蒸散量设定为1990年初始值。基于以上两种情景, 分别模拟产水量并计算水源涵养量。

2 结果与分析 2.1 三江源生态系统分布及其水源涵养量

三江源区生态系统主要类型为草地生态系统、灌丛生态系统、荒漠生态系统、湿地生态系统(图 2)。草地生态系统分布广泛, 面积占约26.07万km², 占三江源总面积的72.09%。灌丛生态系统面积约为1.08万km², 主要分布在长江源西南部及东南部、黄河源东部及东南部、澜沧江南部地区。荒漠生态系统面积约为6.04万km², 主要分布在长江源西北部、黄河源西北部及中部地区。湿地生态系统面积约为2.55万km², 主要分布在长江源西北部及西南部、黄河源北部区域。

草地生态系统为三江源水源涵养功能的主体(图 3), 水源涵养量为120.04亿m³, 占总量的70.59%；荒漠生态系统水源涵养量为35.61亿m³, 占总量的20.95%；灌丛生态系统水源涵养量为7.05亿m³, 占总量的4.14%；湿地生态系统水源涵养量为5.03亿m³, 占总量的2.96%。以单位面积的水源涵养量表征水源涵养能力, 森林、灌丛、草地生态系统水源涵养能力显著高于其它类型。三江源森林生态系统水源涵养能力最高, 为70.05 mm；灌丛生态系统次之, 为64.98 mm；草地生态系统水源涵养能力为46.04 mm。

2.2 三江源水源涵养量的年际变化特征

1990—2020年三江源区及长江源、黄河源、澜沧江源三个子流域水源涵养量均呈增长趋势(图 4)。三江源多年平均水源涵养量为163.84亿m³；水源涵养量显著增加, 速率为1.80亿m³/a(P < 0.05), 从1990年的119.49亿m³增加至2020年的169.19亿m³。三个子流域中, 长江流域和黄河流域水源涵养量呈显著增长趋势, 变化速率分别为0.86亿m³/a(P < 0.05)和0.81亿m³/a(P < 0.05)；澜沧江流域增长速率最小, 且趋势不显著, 为0.13亿m³/a。

以2005年生态治理前后为节点, 将研究期分为生态工程实施前(1990—2005年)、实施后(2005—2020年)两个时段。三江源整体在两个时段内水源涵养量均呈增长趋势, 2005—2020年变化速率(0.16亿m³/a)小于1990—2005年变化速率(3.34亿m³/a)。三个子流域中, 黄河流域和澜沧江流域水源涵养量均呈增长趋势, 生态工程实施前水源涵养量变化速率大于生态工程实施后；长江流域生态工程实施前水源涵养量增加速率为1.82亿m³/a(P < 0.05)；生态工程实施后水源涵养量呈下降趋势, 但变化不显著。

2.3 三江源水源涵养能力的空间分布及变化特征

三江源地区1990—2020年平均水源涵养能力表现出东南高、西北低的空间分布特点(图 5)。长江流域水源涵养能力整体较低, 平均水源涵养能力为39.52 mm, 大部分区域低于40 mm；黄河流域水源涵养能力空间差异显著, 南部区域部分高于100 mm, 北部区域部分低于20 mm, 平均水源涵养能力为49.15 mm；澜沧江流域水源涵养能力整体较高, 大部分区域水源涵养能力大于60 mm, 平均水源涵养能力为62.98 mm。

三江源地区1990—2020年水源涵养能力显著增长面积为22.07万km², 占全区总面积的60.79%。三个子流域中, 黄河流域显著增长区域范围最广, 占该流域面积的77.82%；长江流域中东部区域水源涵养能力显著增长, 占该流域面积的58.64%；澜沧江流域23.05%的区域水源涵养能力显著增长。

2.4 气候变化和土地利用对水源涵养能力的影响

1990—2020年三江源降水量、潜在蒸散量、实际蒸散量、实际蒸散比变化趋势如图 6所示。降水量、潜在蒸散量和实际蒸散量呈显著增加趋势, 变化速率分别为2.67 mm/a(P < 0.05)、1.84 mm/a(P < 0.05)和0.95 mm/a(P < 0.05)。实际蒸散比(实际蒸散量占降水量的比率)呈显著下降趋势(P < 0.05)。

Pearson相关分析表明, 降水、实际蒸散与水源涵养能力存在显著正相关关系, 实际蒸散比与水源涵养能力存在显著负相关关系(表 3)。三江源降水、实际蒸散与水源涵养能力相关系数分别为0.991(P < 0.01)、0.876(P < 0.01)。实际蒸散比与水源涵养能力的负相关系数为-0.953(P < 0.01), 表明相同降雨条件下, 实际蒸散量越大, 水源涵养能力越低。

气候变化情景下, 受降雨量增加影响, 1990—2005年三江源水源涵养能力整体明显增加(图 7)。其中, 长江流域平均增加21.08 mm, 西南部地区呈极显著增加趋势(P < 0.01)；黄河流域平均增加32.25 mm, 东北部地区呈显著增加趋势(P < 0.05)；澜沧江流域平均增加17.35 mm, 西北部地区呈显著增加趋势(P < 0.05)。2005—2020年三江源水源涵养能力整体降低, 除三江源东部小部分区域水源涵养能力增加(0—20 mm)外, 其它地区水源涵养能力普遍下降, 在P=0.05水平上变化趋势不显著。

土地利用/覆被变化情景下, 生态治理前, 1990—2005年三江源水源涵养能力整体降低(图 7), 但变化幅度较小, 长江流域、黄河流域和澜沧江流域水源涵养能力分别变化-0.03 mm、-0.47 mm和-0.01 mm, 在P=0.05水平上变化趋势不显著。生态治理后, 2005—2020年三江源水源涵养能力整体增加, 长江流域西北部及南部地区、黄河流域东部及南部地区、澜沧江大部分地区增长明显, 水源涵养能力增加大于30 mm；长江流域、黄河流域和澜沧江流域水源涵养能力平均增长17.00 mm、19.60 mm和33.76 mm, 呈显著变化趋势(P < 0.05)。

3 讨论

水源涵养能力是生态系统维持健康稳定状况的重要指示器^[3]。三江源是国家重点水源涵养功能保护区^[24], 探明该地区水源涵养能力变化趋势及其影响因素对提高区域水资源供给能力、保障生态系统健康和维持功能稳定具有重要意义^[5]。本研究表明, 1990—2020年三江源水源涵养能力呈上升趋势, 其中黄河流域东部区域水源涵养能力增加显著, 与吕乐婷等^[35]、潘韬等^[9]和张媛媛^[36]研究结果相似。三江源区20世纪九十年代以来降雨不断增多, 温度持续上升^[37], 气候趋向暖湿, 对水源涵养能力提高具有一定正向作用。然而, 本世纪以来, 降水波动幅度大, 强降雨增加, 气温上升趋势显著高于降水, 气候条件逐渐向暖干化发展, 对植被生长有抑制作用, 尤其导致高寒草甸退化^[38], 抑制了水源涵养能力提高。三江源由北至南、由西至东植被覆盖度逐渐增加；长江流域未利用地广布, 西部地区有湿地、冰川分布, 东部地区有草地分布；黄河流域以草地为主, 东部分布有耕地和林地；澜沧江流域植被条件较好, 林、草地覆盖广泛^[39]。国家于2005年启动三江源生态保护与建设工程以来, 退耕还林、休牧育草等措施的实施使生态退化现象得到一定遏制^[40], 林、草植被增加, 水热条件较好的东部地区植被恢复情况有明显好转, 西北部地区植被恢复较为缓慢^[41]。黄河流域东部草地覆盖度显著提高, 部分草地向森林转化；长江流域东部至黄河流域西部草地面积显著增长；澜沧江流域变化较小^[42]。此外, 在气候变化影响下, 21世纪以来长江流域湿地面积呈下降趋势, 人类活动同时干扰湿地景观破碎化程度加剧, 湿地水源涵养能力显著降低^[43], 对长江流域水源涵养能力影响较大。气候变化与人类活动共同影响着三江源水源涵养能力的空间分布格局。

气候变化通过改变降水量和蒸散量来影响流域的水源涵养能力^[16—17], 土地利用/覆被变化引起流域下垫面植被格局的变化, 从而改变了区域水源涵养能力的空间格局。本研究表明, 生态保护与修复工程实施前, 假设土地利用/覆被状况不发生变化, 气候变化对水源涵养能力变化呈正影响, 假设气候状况不发生变化, 土地利用/覆被变化对水源涵养能力变化呈负影响；生态保护与建设工程实施后, 假设土地利用/覆被状况不发生变化, 气候变化对水源涵养能力变化呈负影响, 假设气候状况不发生改变, 土地利用/覆被变化对水源涵养能力变化呈正影响。龚诗涵等^[19]全国水源涵养能力驱动因子的研究表明水源涵养能力与降水、蒸散、温度等气候因素呈显著正相关关系, 这与本研究降水、实际蒸散与水源涵养能力存在显著正相关关系结论基本一致。生态治理前, 三江源土地利用类型变化不显著^[24], 水源涵养能力主要受到降水量增加、实际蒸散量增加和实际蒸散比减少等气候条件变化的影响, 呈增长趋势。生态治理后, 三江源林、草地面积显著增加^[42], 人为干预下, 气候变化对水源涵养能力的负影响得到抑制, 由于森林和草地生态系统水源涵养能力显著高于其它类型生态系统, 退耕还林和休牧育草等措施促进了水源涵养能力提高。

本文模型模拟的中间参量和模拟结果与现有研究进行了大量的对比验证, 实际蒸散量结果与吕乐婷等^[35]、尹云鹤等^[44]的研究结果进行了比较；水源供给量的计算结果与潘韬等^[9]、乔飞等^[45]的研究进行了比较, 并与1990—2020年青海省水资源公报相关数据相比较；水源涵养量的估算结果与张媛媛^[36]、刘敏超等^[46]的结果进行了对比, 研究结果均具有较高的一致性。在模型参数合理条件下, 本文基于情景模拟法, 在生态治理前后两个时段分析了三江源气候与土地利用/覆被变化分别对水源涵养能力的影响。但三江源水源涵养能力的模拟研究还有待进一步提升, 本研究未能充分考虑土地利用/覆被变化对三江源区土壤物理和生态水文过程的影响, 模型中输入的土壤参数未考虑其时间变化特征。此外, 受全球气候变暖影响, 冰川融水和冻土退化作为三江源河流径流量的重要补给来源, 是影响三江源区水源涵养能力的一个重要因素, 这在InVEST模型中未能具体量化, 在未来研究中, 将继续深入。

4 结论

本文基于InVEST模型, 对三江源区1990—2020年水源涵养能力的变化趋势、时空分布及其影响因素进行分析, 取得如下结果：

(1) 草地生态系统为三江源水源涵养功能主体, 多年平均水源涵养量为120.04亿m³。

(2) 1990—2020年三江源区水源涵养量年际变化呈显著上升趋势, 1990—2005年水源涵养量增长速率为3.34亿m³/a, 大于2005—2020年增长速率(0.16亿m³/a)。

(3) 三江源区水源涵养能力空间分布上表现为东南高西北低的特征。1990—2020年三江源水源涵养能力显著增长区域面积22.07万km²。三个子流域中, 黄河流域水源涵养能力增长面积最广。

(4) 对比生态治理前后两时期, 1990—2005年三江源降水量增长和实际蒸散比下降等气候变化对水源涵养能力增长起到主要正向促进作用；2005—2020年三江源生态保护与修复工程实施后, 林草植被覆盖面积增加, 土地利用/覆被变化对水源涵养能力增长起到主要正向促进作用。

林草植被与湿地在生态系统水源涵养功能中起到重要作用, 为扎实推进三江源地区生态保护与修复工程, 应结合地区气候与土壤条件, 科学规划退耕还林、休牧育草以及湿地保护工作, 制定更合理的生态保护政策。