基于Budyko假设的洮河与大夏河径流变化归因识别

杨 林,赵广举,2,*,穆兴民,2,田 鹏,高 鹏,2,孙文义,2,张丽梅

1 西北农林科技大学水土保持研究所,黄土高原土壤侵蚀与旱地农业国家重点实验室, 杨凌 712100 2 中国科学院水利部水土保持研究所, 杨凌 712100 3 西北农林科技大学资源环境学院, 杨凌 712100

日趋剧烈的人类活动深刻改变流域下垫面,进而影响流域的产汇流机制,导致河流的径流量等水文要素变化,从而影响水循环过程和水量平衡[1- 2]。近几十年来,我国在黄河上中游地区实施了大规模的水土保持措施,建设了众多水利工程。截止2011年,累积建设梯田超过4.5万km2,淤地坝5.5万座,植被覆盖度由20世纪80年代的28%增长到近年的63%[3-4]。这些措施极大的改变了黄河干支流水文过程,由此导致径流输沙急剧锐减[5- 8]。黄河水沙变化规律及其归因识别已引起相关行业部门的极大关注。定量分析其变化的原因,对深入理解流域水文演变机理和规律具有重要的科学意义,对流域水资源科学管理具有现实意义。

针对河川径流变化影响因素的评估方法已有较多研究,如广泛应用的水文模型法、水文统计法等[9- 10]。基于Budyko假设的水热耦合平衡方程可以定量评估气候变化和人类活动对径流的敏感性和贡献量[11- 13],目前,已有许多研究学者利用黄河流域的气象资料和水文资料拟合水热关系,证实了Budyko水热耦合平衡方程理论的适用性[14- 16]。

洮河、大夏河作为黄河上游区重要的一级支流,受气候变化和人类活动双重作用,河川径流发生显著的变化[17-18]。径流的丰枯变化必定影响区域的水源供水格局,进而牵动区域内的社会经济的发展[19-20]。目前,针对洮河、大夏河流域的水文研究多集中在水文气象要素变化特征方面,而对其径流变化归因分析研究较少,且多集中在2015年之前[20-23]。鉴于此,本文以洮河、大夏河流域为研究区域,采用Budyko水热耦合平衡方程理论量化降雨、潜在蒸散发、下垫面特征系数对径流变化的影响程度,识别引起流域径流变化的原因,对深刻理解水循环过程和改善流域水资源管理措施具有理论和实践意义。

1 研究区概况和数据来源

1.1 研究区概况

洮河流域位于甘肃省西部,介于101°36′— 104°20′E、34°06′—36°01′N(图1),发源于西倾山北麓,由西向东依次经过碌曲、卓尼、岷县,岷县受阻后北拐穿过九甸峡和海甸峡,经东乡、临洮等地,在永靖注入黄河。流域总面积24973 km2,干流全长673.2 km,流域以白石山、太子山一线为界,分为青藏高原和黄土高原两大地貌单元,地貌单元间地理—生态格局天然分异,由地势平缓的草原区、高山耸立的森林草原区向植被稀疏、水土流失严重的黄土区过渡。红旗水文站是洮河流域的把口站,多年平均径流量为46.8亿m3。

大夏河发源于甘南高原,河源由夏河和咯河(又名：德乌鲁河)两大支流组成(图1)。大夏河流经夏河县、合作市、临夏、东乡县注入刘家峡水库,流域面积6808 km2,主河道长201.9 km,流域内多年平均降水量630 mm,平均蒸发量705.8 mm。折桥水文站是大夏河流域的把口站,多年平均径流量为4.479亿 m3。

图1 研究区水文、气象站点空间分布Fig.1 Locations of hydro-climatic stations in the study area

1.2 数据来源与处理

本文选择洮河流域的红旗水文站和大夏河流域的折桥水文站(图1),通过黄河水文年鉴、《黄河泥沙公报》及黄土高原数据中心(http://loess.geodata.cn/hyd/inedx.jsp)获取两个水文站的月径流数据。通过中国气象数据共享网(http://data.cma.cn)下载获得流域范围内12个气象站点1961—2017年的基础气象数据,包括最高温度、最低温度、平均温度、日照时数、平均相对湿度、日降雨量、平均风速。采用世界粮农组织(FAO)修正的Penman-Monteith公式获取研究区的潜在蒸散发,该方法充分考虑了作物地表特征和空气动力学的参数的变化[24],广泛应用于蒸散发计算中[25-26],计算公式如下：

(1)

式中,ET0为潜在蒸散量(mm/d);Δ为饱和水汽压曲线斜率(kPa/℃);Rn为太阳净辐射(MJ m-2d-1);G为土壤热通量(MJ m-2d-1);γ为干湿常数(kPa/℃);U2为2 m高处的风速(m/s);T为平均气温(℃);es为平均饱和水汽压(kPa),ea为实际水汽压(kPa),(es-ea)为饱和水汽压差(kPa)。

2 研究方法

2.1 Mann-Kendall 趋势检验和突变检验

Mann-Kendall(MK)趋势检验和Mann-Kendall突变检验分析是水文时间序列常用的方法,能够很好地反映河川径流变化发展规律,具有精确度高、不受异常值干扰、适用于非正态分布数据的优点。计算方法详见参考文献[8, 27]。

2.2 径流变化归因识别

基于Budyko水热耦合平衡理论量化气候变异和人类干扰对径流变化的敏感性和贡献率,能够准确、快速识别引导径流变化的主导因素,在研究气候变化和下垫面变化对径流影响程度分析中广泛应用[15, 28]。

2.2.1敏感性分析

基于Budyko水热耦合平衡方程计算长时间尺度下的实际蒸散发[29],其方程基本形式为：

(2)

式中：ET表示流域年实际蒸散发(mm),P表示年降雨量(mm),ET0表示年潜在蒸散发(mm),参数ω代表流域的下垫面参数(无量纲),与土壤性质、地形和植被相关[30-31]。

在多年尺度上,闭合流域满足水量平衡关系,即：

R=P-ET

(3)

式中,R为多年平均径流深(mm)；P为多年平均降雨量(mm)；ET为多年平均实际蒸散量(mm)。将方程(2)代入方程(3)可以得到：

(4)

假设P、ET0和ω是独立的变量,公式(3)可以改写为R=f(P,ET0,ω),径流深可以表达为：

(5)

(6)

进一步推导,各弹性系数表达式为(Ø=ET0/P)：

(7)

(8)

(9)

降雨、潜在蒸散发、下垫面的径流弹性系数由(7)—(9)式求得,弹性系数的正负表示径流与该因子的正负相关关系,弹性系数的绝对值大小表示该因子对径流变化的影响程度。

2.2.2径流变化的贡献率分析

如果径流发生显著性变化,基于突变点,研究时期分为基准期和若干人类活动期,基准期和人类活动期的径流深变化为ΔR,理论上径流深变化可以归因于降雨、潜在蒸散发、下垫面状况改变三方面,由公式(10)计算单个因素引起的径流深变化表达为：

(10)

式中i代表降雨P、潜在蒸散发ET0、下垫面系数ω。

通过计算可以分别获得人类活动期降雨、潜在蒸散发、下垫面系数相对基准期引起的径流深变化为ΔRP、ΔRET0、ΔRω,从而计算单个因素对径流变化贡献率为：

ηi=ΔRi/ΔR×100%

(11)

其中ηi表示降雨、潜在蒸散发、下垫面系数对径流变化的贡献率。

3 结果与分析

3.1 水文气象要素特征

洮河和大夏河流域的径流特征如表1所示,结果表明：洮河流域的多年平均径流深为177.71 mm,是大夏河121.49 mm的1.41倍。大夏河年际径流深极值比为22.66,是洮河的5—6倍,表明大夏河年际径流深差异明显,这主要与大夏河南北地域差异性和季风性气候引起的降雨时空分布有关。采用非参数MK趋势检验流域径流深的变化趋势可知：洮河、大夏河流域的径流深均呈现减少趋势,显著性超过99%。图2为流域年径流深的线性拟合,由此可知,洮河、大夏河1961—2017年径流深减少速率分别为-1.85 mm/a、-1.36 mm/a。Mann-Kendall突变点检验发现,洮河、大夏河径流突变分别发生在1987、1985年,突变之后人类活动期径流相对变化率超过30%。

表1 流域水文气候变化特征

图2 径流变化趋势和突变检验Fig.2 Runoff changing trends and abrupt changesUFk：MK正序检验统计量 MK normal series test statistics； UBk：MK逆序检验统计量 MK inverse series test statistics

流域降雨、蒸发等气象因素直接影响着河川径流变化。采用MK趋势检验分析研究区降雨、潜在蒸散发变化特征发现(表1),洮河、大夏河年降雨量无显著的变化趋势,年降雨量变化速率分别为-1.08、0.005 mm/a(图3)；而洮河和大夏河流域年潜在蒸散发显著增长(P≥0.01),增长速率分别为0.6、1.46 mm/a,2000年之后年潜在蒸散发均超过850 mm,这与黄河流域整体增温变暖一致。

图3 研究区气象变化趋势Fig.3 Climate changing trends in the study area

3.2 径流变化归因分析

3.2.1径流对气候和下垫面变化的敏感性分析

基于Budyko水热耦合平衡方程理论计算洮河流域和大夏河流域各时段水文特征以及径流对降雨量、潜在蒸散发、下垫面参数的弹性系数(表2)。分析发现：洮河、大夏河流域在人类活动期降雨量、径流深相对基准期均减少,年潜在蒸散发增加,径流系数减少,降雨转换为蒸散发的比例增加。

表2 气象水文变量特征

总体而言,洮河、大夏河流域径流变化与ET0、ω呈负相关,与P呈正相关。对比突变点前后的两个时段,洮河流域的降雨弹性系数由1961—1987年的1.67,增加到1988—2017年的1.92,说明同样的降水增加10%时,1987年以前会使洮河径流增加16.7%,则会使1987年之后增加19.2%,表明降水对径流的影响程度进一步加强；同理,潜在蒸散发或者下垫面参数ω增加10%时,1987年以前分别会使洮河径流减少6.7%或者10.8%,1987年之后分别会使径流减少9.2%或者13.1%。同理,大夏河流域的降雨、潜在蒸散发、下垫面特征参数分别增加10%,1985年以前分别会使大夏河径流增加19.2%、减少9.2%、减少11%,1985年之后径流分别增加23.1%、减少13.1%、减少16%。通过对比降雨、潜在蒸散发、下垫面参数ω弹性系数绝对值大小,降雨弹性系数最大,下垫面次之,潜在蒸散发最小,揭示了洮河、大夏河对降雨补给最敏感,降水量的增减,直接决定流域水量的多少；下垫面对流域降水拦蓄截留,改变产流过程的填洼、下渗和土壤蓄水能力,汇流过程中的地表糙率、河流水文路径等水文过程。

图4反映了洮河、大夏河流域1961—2017年干旱指数(ET0/P)和降雨、潜在蒸散发、下垫面弹性系数的年际变化趋势,干旱指数呈现上升趋势,流域内暖干化现象显著,降雨、潜在蒸散发弹性系数绝对值呈现显著的上升趋势(P<0.01),反映了流域径流变化对气象要素的敏感性增强。

图4 干旱指数、弹性系数年际变化Fig.4 Interannual changes of drought index and elasticity coefficients

3.3 径流变化的归因识别

基于Budyko水热耦合平衡方程理论量化气候变化和下垫面条件对径流过程的影响程度(表3)。由表可知,理论计算的径流深变化(dR)和实际径流深变化(dR＇)数值接近(二者误差的绝对值δ<5),表明本研究评估径流变化归因的方法是可行的。洮河流域和大夏河流域人类活动期相对于基准期,降雨、潜在蒸散量和下垫面状况对径流的贡献有较大差异。与基准期相比,洮河流域的径流减少31%,下垫面变化对径流减少的贡献占比最高(63.68%),其次为降雨,贡献33.99%,潜在蒸散发贡献最小。大夏河流域人类活动期相比于基准期,下垫面变化对径流贡献为76.87%,降雨和潜在蒸散发贡献接近,分别贡献10.50%和12.63%。

综上,下垫面条件是引起洮河、大夏河流域径流减少的主要因素,降雨次之,潜在蒸散发影响程度最小。流域2000年之后降水颇丰,但径流深基本保持在150 mm左右,可见下垫面状况的改变对流域的产汇流过程产生了深刻的影响。

表3 径流变化归因识别

4 讨论

4.1 流域特征参数ω的含义

Budyko方程的唯一参数ω,是一个区域下垫面特征综合体现,一般认为主要受地形、土壤、植被因素的影响。已有学者针对人类活动影响较弱的地区建立了区域相对入渗能力、相对土壤水储量、平均坡度三者的函数求解ω的关系[14],研究发现流域在较短时间内地形和土壤性质的变异性相对较弱,参数ω主要反映下垫面植被的变化。20世纪末,黄土高原实施退耕还林(草)、荒山造林等一系列大规模生态修复工程,植被状况显著改善,覆盖度由上世纪80年代的28%增加至近年的63%[21]。同时,学者针对洮河流域内1980—2010年归一化植被指数(NDVI)研究分析发现,春、秋季植被生长旺盛,流域范围内NDVI显著增加[32]。洮河、大夏河流域范围内梯田的修建、植被覆盖的增加深刻地改变了流域下垫面的状况,人类活动期的下垫面参数ω相比于基准期分别增长了19.05%、25.6%,表明流域中降水转化为蒸散发的比例增加,这与流域中植被覆盖度变化密切相关,植被覆盖度的增加导致冠层截留增加,同时改变了产流过程中的填洼、下渗和土壤蓄水能力和汇流过程中地表糙率、地表蓄水量、河道汇水路径等水文过程,延长产汇流时间,增加了蒸散发。1990年之后径流减少30%以上,植被变化是控制流域的下垫面状况的主导因素,但是流域中植被变化驱动的水文过程变化机理研究相对较少,后期需要进一步加强。

4.2 河川径流变化主要驱动因素

采用Budyko水热耦合平衡方程理论对径流影响因素的弹性系数进行分析发现,气候因子的弹性系数绝对值最大,表明气候变化直接影响着流域水文循环过程。洮河、大夏河流域人类活动期的径流较基准期均呈显著减少,尽管下垫面的弹性系数介于1—1.3之间,但下垫面状况的改变对径流深变化的贡献率超过了60%,与已有的研究结果一致[23]。流域内降雨年际变化不均匀,平水年和枯水年交替出现,洮河、大夏河流域在2003年之后,降水量相对颇丰,但径流深分别减少26.08%、28.61%,表明流域内径流的改变与降雨变化不同步,因此,流域下垫面变化是河川径流锐减的主导因素。统计研究区的治理措施发现,20世纪70年代末,洮河、大夏河流域均实施小流域综合治理,修筑梯田、淤地坝等水土保持措施,并于1999年开展退耕还林草工程,截止2015年,洮河流域修筑梯田2116 km2,退耕还林面积达561 km2,修建水利水电工程及水库62座,中型淤地坝6座,灌溉面积达到6.03万hm2[23, 33]；大夏河流域至2010年底,兴修梯田2.11万hm2,营造林草面积达2.41万hm2,修建谷坊380座、淤地坝7座,水土流失治理率初步达到43.5%[34],这些举措深刻地改变了流域下垫面条件,影响了产汇流机制,由此可以看出,研究区下垫面的变化深刻地改变了流域的水文过程。

5 结论

本文以黄河上游支流洮河、大夏河流域为研究对象,采用Mann-Kendall趋势检验、突变等方法分析流域1961—2017年径流序列变化,基于Budyko水热耦合方程理论量化气候变化和人类活动对径流变化的贡献,主要结论如下：

洮河流域、大夏河流域径流深均显著减少,变化速率分别为-1.85 mm/a、-1.36mm/a。两个流域的年降雨量未发生显著变化,但年潜在蒸散发量显著增加(P<0.01)。采用Mann-Kendall突变检验识别洮河、大夏河年径流阶段变化特征,人类活动影响较小的时期分别为1961—1987、1961—1985年,人类活动影响剧烈期为1988—2017年、1986—2017年。

依据Budyko水热耦合平衡方程识别径流变化的主导因素,研究发现人类活动引起的洮河、大夏河流域下垫面变化是径流变化的主要因素,分别贡献63.68%和76.87%。洮河、大夏河流域的径流减少主要受下垫面变化的控制,降雨量次之,潜在蒸散发影响最弱。