三江源区地处青藏高原腹地, 是长江、黄河、澜沧江的发源地, 是中国以及东南亚的重要产水区, 同时也是亚洲乃至北半球气候变化的“感应器”和“敏感区”。三江源区特殊的地理位置、丰富的自然资源、重要的生态功能使其成为青藏高原生态屏障的重要组成部分。然而, 由于该地区海拔高、自然条件恶劣, 其生态环境十分脆弱, 生态系统质量的变化非常敏感^[1—3], 尤其是20世纪70年代以来, 在气候变暖和人为活动不断扩大的双重影响下, 冰川消退、江河断流、草地退化、水土流失等问题相继涌现, 三江源区生态环境质量呈恶化趋势^[4—6]。党和国家高度重视三江源区生态保护工作, 通过建立自然保护区、成立国家公园等体制机制的逐步完善和一系列生态保护与恢复工程的实施, 三江源区生态系统退化趋势初步得到遏制^[7—8]。

气候变化和人类活动的双重作用对三江源区水资源供给和径流调节功能造成了重要影响, 不少学者针对三江源区的水资源时空分布、演化规律和径流调节能力等进行了研究分析, 如：罗玉等^[9]运用滑动平均、累积距平曲线、集中期、集中度、滑动t检验、小波分析等方法统计分析了长江源区径流量变化特征；刘希胜等^[10]采用过程线法、集中度和集中期等方法分析了黄河流域降水、径流的演变特征, 研究了径流对降水的响应特征；张士锋等^[11]借助气候驱动模型分析了三江源区径流变化的主要气候驱动因素；蒋冲等^[12]通过双累积曲线、相关分析和贡献率分析等方法探讨了影响三江源区年际和年内径流量变化的气候气象因素；张小咏^[13]基于春汛期、夏汛期、枯水期三个不同时段水文参数, 分析了气候变化和生态系统退化或转好对生态系统径流调节功能的影响。

径流调节反映了生态系统通过对降水截留、过滤、吸收等, 对其时空再分配的过程, 在一定程度上起到削峰补枯、缓和地表径流、增加地下径流等的作用, 通常采用水量平衡法、降水贮存量法、综合蓄水能力法等方法评估^[14—16]。当前针对径流调节等反映生态质量、生态功能的评估多是基于可观测值或相关可观测值的评估, 虽然能够直观给出某一特定时期的变化情况, 但由于其缺乏对应时期评价指标的参考值, 难以量化其与该区域最优状态之间的差距, 量化生态系统的恢复潜力^[17], 在一定程度上制约了生态保护措施的开展实施。生态系统参照系是指某地区在长期气候、地形、海拔等自然因素的影响下, 经自然演替形成的最优生态系统及其所具有的生理生态特征, 具有生物群落发展到顶级植被的生态系统特征, 它反映了一种自然原始、未受或极少受人为活动干扰的条件状况。在开展生态系统恢复潜力或生态系统保护成效分析中, 通过选取生态系统参照系, 量化参照系条件下关键指标的数值, 并把它作为长期的评估标准, 就可以科学地量化生态系统的现状与理想状态或最优状态之间的差距。通过建立客观的生态系统评估参照系, 定量评估参照系和实际条件下的径流变化过程, 有利于厘清自然因素和人为因素的影响, 量化生态系统的恢复进程, 实现不同时间、不同区域生态系统径流调节功能的可比较。

对此, 本研究以科学可比较为评估目标, 以径流量和径流系数作为径流调节功能表征指标, 建立了基于历史参照系的三江源区径流调节功能评估模型, 分析了2000—2017年三江源区径流调节功能变化特征；采用参照条件下的径流系数与实际条件下的径流系数的比值作为径流系数质量指数, 计算了基于参照系下的三江源区径流系数和地表径流系数质量指数得分变化, 分析了三江源区径流调节功能现状和变化情况。

1 研究区概况

三江源区位于我国青海省南部, 是长江、黄河、澜沧江的发源地, 被誉为“中华水塔”。多年平均水资源量约为500亿m³, 分别占长江流域总水量的20%, 黄河流域总水量的40%和澜沧江流域总水量的15%^[18]。其行政范围包括玉树、果洛、海南、黄南藏族自治州的21个县和格尔木市的唐古拉山镇, 区域总面积39.5万km²。三江源区是典型的高原大陆性气候, 冷热两季交替、干湿两季分明。年平均气温由东南向西北随着纬度和海拔的升高而逐渐下降, 区内温度梯度差异明显, 降水量分布总体上也呈现出由东南向西北逐渐递减的趋势^[18—19]。三江源区土壤以高山草甸土为主, 且土壤空间分布具有明显的垂直地带性规律, 区域内以高寒草甸和高寒草原等草地生态系统为主, 其次是湿地生态系统, 此外还少量分布着森林、灌丛等植被类型^[20—21]。植被随降水变化呈地带性分布, 降水量最低的长江源区, 其西部最干旱区以荒漠草原为主, 逐渐向东部过渡为以针茅为主的高寒草原, 东部主要是以高寒嵩草为主的高寒草甸。降水量最高的澜沧江源区, 其植被类型以高寒草甸为主, 东部边界谷地分布着落叶阔叶灌丛和针叶林。黄河源区内植被整体仍以高寒草甸为主, 其北部以温带草原和高寒草原为主, 中部偏南部有落叶灌丛分布。

2 数据与方法 2.1 数据

径流数据来源于中国水利统计年鉴, 包括长江流域直门达站、黄河流域吉迈站和唐乃亥站、澜沧江流域香达站的逐月径流量。气象数据来源于中国气象共享网, 包括三江源区内气象站点的日降水量、日平均气温、日最高气温、日最低气温、日照时数、日平均相对湿度和日平均风速数据, 及辐射站点的日太阳辐射数据(站点分布图见图 1)。植被类型数据来源于中科院植物所的1:100万中国植被类型图^[22]。生态系统类型数据来源于中国科学院资源环境科学数据中心, 空间分辨率为30 m, 共包括2000, 2005, 2010和2015年4期生态系统类型图, 本研究以各期生态系统类型图代表最近5年的生态系统类型情况(2000—2002年采用2000年数据, 2003—2007年采用2005年数据, 2008—2012年采用2010年数据, 2013—2017年采用2015年数据)。地形高程数据来源于地理空间数据云, 空间分辨率为90 m；土壤类型图来源于国家青藏高原科学数据中心“基于世界土壤数据库的中国土壤数据集(v1.1)(2009)”。

2.2 研究方法

径流系数是表征流域降雨产流能力的重要参数, 它是指某一时段内汇水面积上的产流量与相应时段内的降雨量的比值, 反映了流域内自然地理因素对降水形成径流过程的影响。径流系数主要受降水总量、降水强度和下垫面影响^[23—25], 具有明显的区域特征。本研究利用SWAT (Soil & Water Assessment Tool)模型模拟区域径流量, 选取径流系数和地表径流系数来反映三江源区生态系统的径流调节功能, 以参照系径流系数与实际径流系数的比值作为径流系数质量指数, 评价实际生态系统相较于优良生态系统径流调节功能的偏离情况, 分析三江源区生态系统径流调节功能参照值、现状值及变化量的时空变化规律。

2.2.1 参照系确定

生态系统参照系的选取应尽量体现未经过自然与人为因素扰动后的生态系统特征^[17]。研究表明, 三江源区草地退化格局在20世纪70年代中后期初步形成并快速发展^[5]。20世纪70年代中期—2000年, 达日县发生退化的草地面积约占全县面积的29.39%, 玛多县约有70%的天然草场面积发生退化, 其中大部分为重度退化；黄河源头20世纪80—90年代草地退化速度相比20世纪70年代增加了1倍^[26—30]。对此，本研究选取成图在上世纪70年代的1:100万植被类型图作为生态系统参照图层。结合植被类型图, 将三江源区内生态系统类型重分类, 划分为常绿针叶林、落叶阔叶林、落叶阔叶灌丛、常绿阔叶灌丛、常绿针叶灌丛、高山草地、草原、草甸、农田、湿地、冰川积雪、水体、裸露石山、荒漠14个二级类别(表 1), 结合2000—2015年生态系统类型图, 进一步归为森林、灌丛、草地、农田、冰川积雪、湿地、水体和裸地8个一级类别。

2.2.2 径流模拟与验证

SWAT模型是美国农业部农业研究局开发的流域尺度分布式水文模型。SWAT模型对流域水文循环过程的模拟主要包括坡面产流过程和河道汇流演算过程, 主要模拟的水分分量包括地表径流、蒸散发、壤中流、地下水四个部分, 其中地表径流采用SCS(Soil Conservation Service)曲线法进行估算, 蒸散发考虑了层截留蒸发、植物蒸腾和土壤水分蒸发, 壤中流采用动力贮水方法计算, 地下径流主要考虑河道地下水和蓄水层补给量、基流退水系数等计算。

由于1970年以前气象数据缺测较多, 本研究基于1970年以后的各水文站逐月径流数据进行参数率定与检验。具体设置为：1970年作为模型缓冲期, 1971—2000年作为模型率定期优选模型参数, 对于直门达、唐乃亥站以2000—2017年为模型检验期；香达和吉迈站以2000—2012为模型检验期。并采用纳什效率系数(Nash-Sutcliffe efficiency coefficient, NSE)、相对误差(Re)和决定系数(R²)3个指标评价模拟精度^[31—32]。

2.2.3 径流系数质量指数

以参照径流系数与实际径流系数的比值作为径流系数质量指数, 评价实际生态系统相较于优良生态系统径流调节功能的偏离情况, 具体计算公式如下：

式中, RA_rcⁱ为第i个像元径流系数质量指数, R_i为第i个像元的实际径流系数, R_mi为第i个像元的参照径流系数。当R_i等于0时或RA_rcⁱ>1, RA_rcⁱ取为1。

S_rcⁱ计算公式如下:

式中, S_rcⁱ为第i个像元径流系数质量指数得分, S_R为实际径流系数质量指数偏离参照生态系统的比例系数, 取值为10。当S_rcⁱ大于或等于100时按100计, 径流系数质量指数得分越高, 说明实际生态条件下的径流系数越接近参照条件下的径流系数, 两者偏离越小。

县域尺度径流系数质量指数得分为该县域内各像元径流系数质量指数得分的均值：

式中, S_rc^c为县域径流系数质量指数得分, S_rcⁱ为第i个像元径流系数质量指数得分, n为县域像元总个数。

3 结果与分析 3.1 参数率定与验证

结合已有文献, 并利用SWAT-CUP软件对模型参数进行敏感性分析, 发现三江源区各水文站控制范围内的地表水文过程参数CN₂(径流曲线系数)、ESCO(土壤蒸发补偿系数)、SOL_AWC(土壤可利用水量)较为敏感；地下水文参数ALPHA_BF(基流α系数)、GW_DELAY(地下水滞后系数)较为敏感。此外长江流域和澜沧江流域CH_K2(主河道有效导水率)不全为0, 3个流域的主要模型参数情况见表 2。表 3中各水文站点径流量验证结果显示, 率定期和检验期各站点径流模拟结果均具有较好的精度^[33], SWAT水文模型在三江源区具有较好的适用性。

3.2 生态系统类型变化

结果显示(图 2), 20世纪70年代, 草地和湿地为三江源区主导生态系统类型, 草地面积占比高达80.58%, 湿地面积占比为7.37%。70年代至2000年期间, 三江源区草地生态系统退化严重, 大量草地向裸地转化, 草地面积减少比例为16.77%, 主要集中在长江源西部和北部区域, 黄河源区北部区域；湿地生态系统面积也大幅减少, 减少比例为51.63%, 主要集中在长江源西南部区域和黄河源西部源头区。2000—2015年, 草地和裸地为主导生态系统类型, 草地、裸地面积占比分别在66.86%—67.07%之间, 21.22%—21.36%之间。各生态系统面积整体变化不大, 其中, 2000—2005年, 森林、灌丛、草地和农田面积共减少685.37 km², 裸地增加584.99 km²。2005—2015年, 草地和裸地面积分别减少153.88 km²、324.17 km², 森林、灌丛和农田面积基本稳定。

3.3 径流调节功能时空动态变化 3.3.1 径流量时空变化

2000—2017年三江源区实际多年平均径流量为495.15亿m³, 比参照系高26.78亿m³, 各年实际径流量均高于参照系径流量, 2012年两者差值最大, 约为48.83亿m³, 2015年差值最小, 约为11.47亿m³, 参照系径流和实际径流的历年差别率在4%—9%之间(图 3)；2000—2017年实际多年平均地表径流量为96.64亿m³, 比参照系高28.04亿m³, 各年实际地表径流量均高于参照系, 参照系地表径流和实际地表径流的历年差别率在22%—58%之间(图 4), 表明三江源区生态系统的退化增加了总径流, 尤其是地表径流显著增加, 而壤中流和地下径流则明显减少。

2000—2017年三江源区多年平均径流空间分布图(图 5)显示, 参照系和实际生态系统条件下的径流空间分布特征基本一致, 总体呈现东部和南部高, 西部和北部低的分布特征, 降水空间差异是影响径流空间差异的主要因素。2000—2017年, 三江源区多年平均径流为125.35 mm, 其中43.71%的区域低于100 mm, 37.81%的区域在100—200 mm之间。从图 5可看出, 三江源区大部分地区实际径流高于参照径流, 其中黄河源区曲麻莱东部和玛多县实际径流与参照径流差别最大；实际径流与参照径流差别较小的地区主要分布在湿地面积明显增加的唐古拉山镇(约增加了3810 km²)、草地面积增加明显的玉树县(约增加了1663 km², 增长约15.09%)以及部分草地转化为森林生态系统的囊谦县中东北部、兴海县南部等地区。

2000—2017年三江源区多年平均地表径流空间分布图(图 6)显示, 参照系和实际生态系统条件下的地表径流空间分布特征基本一致, 但和总径流空间分布有较大差异, 总体呈现西部高, 东部低的分布特征。2000—2017年, 三江源区多年平均地表径流为24.43 mm, 其中61.72%的区域低于20 mm, 35.80%的区域在20—100 mm之间。从图 6可看出, 三江源区大部分地区实际地表径流高于参照地表径流, 两者差值大于20 mm的高值区主要在治多县北部的裸地生态系统分布区、称多县东部和海南州东部的农田生态系统分布区、杂多县的湿地大面积退化区。降水量小且湿地面积增幅明显的唐古拉山镇、植被覆盖度增加的杂多县东南部、玉树县、治多县东南部, 共和县、贵德县和贵南县实际地表径流与参照值差别较小。

3.3.2 径流系数时空变化

2017年三江源区径流系数为0.26, 地表径流系数为0.05(图 7), 参照系条件的径流系数为0.23, 地表径流系数为0.04。2000—2017年, 三江源区实际条件下的径流系数和地表径流系数均高于参照系条件下的数值, 且两种条件下的时间动态变化较为一致, 年际变化均不明显。从年均值来看, 2000—2017三江源区多年平均径流系数为0.27, 地表径流系数为0.06, 均比参照系条件下的数值高出0.02。

三江源区多年平均径流系数与径流量的空间分布特征基本一致(图 8), 地表径流系数与地表径流量空间分布特征基本一致(图 9)。实际条件下的径流系数和参照系条件的径流系数差别不大, 但前者明显高于后者。相较于参照生态系统, 2000—2017年间, 三江源区产流能力尤其是地表产流能力明显增加。

3.4 基于参照系的径流系数质量指数变化

2000—2017年, 三江源区径流系数质量指数得分均值为98.63, 历年数值变化在98.29—99.07分之间；地表径流系数质量指数得分均值相较径流系数质量指数低了1.68分, 历年数值变化在96.43—97.71分之间。三江源区径流系数质量指数得分呈先降低后基本平稳趋势, 地表径流系数质量指数得分呈先降低后上升趋势, 表明区域内生态系统质量呈先退化后逐渐恢复的变化特征。从得分直方图可知(图 10), 各县(镇)径流系数得分≥98.86的数量占比(59.09%)最多, ≤97.65的数量占比(13.64%)最少；地表径流系数质量指数得分≥98.33的数量占比(40.91%)最多, ≤96.39的数量占比(22.73%)最少。

各县(镇)径流系数质量指数得分和地表径流系数质量指数得分有较大差别(表 4)。2017年, 三江源区径流系数质量指数得分高值区为班玛、甘德、久治县, 低值区为唐古拉山镇、治多和杂多县。地表径流系数质量指数得分高值区为久治、河南、甘德县, 低值区为玛多、囊谦和曲麻莱县。2000—2017年三江源区各县(镇)径流系数质量指数得分变化介于-1.94—0.69之间, 得分变化率介于-1.95%—0.71%之间；各县(镇)地表径流系数质量指数得分变化介于-0.34—1.84之间, 得分变化率介于-0.35%—1.90%之间。

4 结论与讨论 4.1 讨论

降水数据是SWAT模型径流模拟最重要的气象输入数据, SWAT模型假设一个HRU内的下垫面具有相同的水文响应属性, 在同一子流域内部的所有HRU采用同一雨量站数据, 即降水量的空间精度为一个子流域。由于三江源区域内及周边国家基本气象站较少, 本研究根据子流域划分情况, 采用临近子流域内已有气象数据的加权平均对无气象站点的子流域进行气象数据插补, 在一定程度上减少了气象站点不足对径流模拟结果的影响。此外, 光照也是植被所必需的气候因子, 光照强度直接影响植被光合作用的强弱。尤其是在水分充沛的地区, 植被活动受到光照的影响最为明显^[34—35]。但三江源区域内及周边太阳辐射站不足10个, 站点插值也难以弥补, 对径流空间模拟结果有一定影响。

基于植被类型图和生态系统类型面积变化的分析显示, 20世纪70年代至2000年期间, 三江源区草地生态系统退化严重, 草地面积减少了约5.25万km², 减少比例为16.77%。草地面积的大幅减少使得生态系统对降水径流的蓄纳能力降低, 导致实际生态系统的径流量和地表径流量均高于参照系。通过参照系和实际径流的差值分析发现, 参照系径流和实际径流的差别率、参照系地表径流和实际地表径的差别率均在2003—2007年最大, 2008年后有所减少。主要原因可能如下：三江源区生态保护与修复工程开始于2005年, 2005年之前草地退化严重, 植被生态质量较差^[36], 此后, 受生态保护与修复工程的推进和气候暖湿化的影响, 区域内草地退化问题基本得到遏制, 生态质量显著好转^[37—39]。此外, 受数据限制, 未获取到历年生态系统类型图, 本研究分别以2000, 2005, 2010和2015年生态系统类型图作为近5年基础数据, 也导致2003—2007年的评估结果较为接近。

基于参照系的评估体系不仅能够有效避免气候波动等造成的生态系统功能年际变化的不确定性, 也可以更科学准确地反映出生态系统功能现状及其变化。在参照系条件下, 径流调节功能的年际空间变化主要由降水等气象因素变化引起, 而在实际生态系统条件下还受生态系统变化影响。对于特定年份, 降水量、生态系统类型、土壤类型和地形空间分布是影响径流调节功能空间分布的主要原因, 其中降水量对径流系数贡献量最大, 呈现径流系数随降水量增加而增加的空间分布特征, 但对径流系数质量指数的影响则有所区别, 生态系统类型变化是影响其空间分布的主要原因；对于不同年份, 生态系统类型的变化是影响径流调节功能空间分布变化的主要原因。单纯基于实际径流的变化不能代表其现状值和变化的真实情况, 结合参照系和实际径流差别率的变化可以更加准确地反映出三江源径流调节功能的年际变化, 根据径流系数质量指数可以反映出实际径流调节功能与参照系径流调节功能的差距。

4.2 结论

(1) 2000—2017年, 三江源区多年平均径流量为495.15亿m³, 地表径流量为96.64亿m³, 参照系条件下数值分别为468.37亿m³和68.60亿m³。整体而言, 三江源区年径流量和地表径流量均高于参照系。2000—2017年, 三江源区参照生态系统和实际生态系统的多年平均径流量整体呈东南高西北低的空间分布特征, 地表径流量则呈西部高东部低的空间分布特征。

(2) 2000—2017年三江源区多年平均径流系数为0.27, 地表径流系数为0.06, 均比参照系条件下的数值高0.02。2000—2017年, 三江源区历年实际条件下的径流系数和地表径流系数均高于参照系条件下的数值, 且两种条件下的时间动态变化较为一致, 年际变化均不明显。三江源区径流系数与径流量空间分布特征基本一致, 地表径流系数与地表径流量空间分布特征基本一致。

(3) 2000—2017年, 三江源区径流系数质量指数和地表径流系数质量指数得分均值分别为98.63和96.86。区域径流系数质量指数得分呈先降低后基本平稳趋势, 地表径流系数质量指数得分呈先降低后上升趋势, 地表径流系数质量指数相较径流系数质量指数变化更为明显。