ABSTRACT

Quantitative analysis of the causes of runoff changes in the basin can make us better use of water and soil resources. The Mann Kendall, Pettitt mutation point test and cumulative anomaly curve methods were used to process the meteorological data in the Kuye River basin from 1971 to 2016. Based on Budyko water and heat balance theory, the attribution of runoff change in Kuye River was identified, and the contribution rate of meteorological factors and active factors to runoff change was obtained. The results show that: rainfall and potential evaporation are not significantly increased, runoff is significantly decreased, and the contribution rate of human factors to runoff reduction is 98%. In human activities, through the change of land use type map, it can be seen that a large number of coal mining, returning farmland to forest and grassland, and urban and rural residents’ water use led to a sharp decrease in runoff.

Keywords:Kuye River Basin, Budyko Theory, Water and Heat Balance, Attribution Analysis, Land Use, Climate Change

基于Budyko水热平衡理论的窟野河径流变化归因识别

郑嘉昊，蒋晓辉^*

西北大学城市与环境学院，陕西 西安

收稿日期：2020年9月24日；录用日期：2020年10月6日；发布日期：2020年10月13日

摘 要

黄河流域水沙变化剧烈，为更好地利用和开发水资源，定量的分析出引起河川径流变化的原因。采用Mann-Kendall、Pettitt突变点检验和累积距平曲线法，对1971~2016年窟野河流域内气象数据进行处理，基于budyko水热平衡理论对窟野河径流变化归因识别，得出气象因素和认为活动因素对径流变化的贡献率。结果表明：降雨和潜在蒸发均为非显著性增加，径流显著性下降，人为因素的下垫面变化对径流减少的贡献率达到98%是最主要的影响因素。人为活动中，通过土地利用类型图的变化，看出大量的煤矿开采，退耕还林还草，以及城乡居民用水导致了径流的急剧减少。

关键词 :窟野河流域，Budyko理论，水热平衡，归因分析，土地利用，气候变化

Copyright © 2020 by author(s) and Wuhan University.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 概述

由于气候和认为的因素导致近50年来黄河流域径流发生着剧烈的变化，引起人们的广泛关注，黄河流的气候条件却没有出现显著性的变化，人类活动就起到了很关键的影响因素。窟野河流域有大量的煤矿和村庄，所以探究各种因素对流域径流变化的贡献率就非常重要，可以有效地对水资源进行合理可持续的利用。

流域径流的变化在年尺度上受到通过气候变化和人类因素的影响，Bukyko假设构建了流域水量与能量耦合平衡方程 [1]。近年来很多研究人员基于这一假设提出了一系列反应下垫面参数的经验公式并验证了假设合理性 [2]。孙福宝等 [3] 通过对黄河中游多个子流域气候变化和下垫面变化进行分析验证了Budyko水热耦合平衡方程在黄土高原流域的应用合理性。郭巧玲等 [4] 认为窟野河流域1997~2010年，该阶段人类活动对径流减少的影响所占比例高达93.62%。刘二佳等 [5] 认为窟野河流域径流变化在1980年~2005人类活动占78.25%其中煤炭开发的影响大于水保措施导致的径流减少。白乐等 [6] 采用SCRAQ法估算认为在1997~2010年的变化期人类活动对降水的影响为87.99%。本文在Budyko水热耦合理论的基础上计算降水、潜在蒸散发人类活动的弹性系数，最终定量化的得出三者对径流变化的贡献率，为窟野河水资源利用提供理论支持。

2. 研究区概况与数据获取

2.1. 研究区概况

窟野河流域，黄河一级支流，坐标范围为108˚28'E~110˚45'E，38˚22'N~39˚50'N (地理位置见图1)，地处黄河中游支流部分，发源于内蒙古自治区巴定沟，经伊金霍洛旗和陕西省府谷县境内，于神木县注入黄河。两条主要支流乌兰木伦河和悖牛川河自神木境内汇合，由温家川站流出，全长242 km，流域面积达8706 km²，窟野河流域因属于干旱半干旱的大陆性气候，径流量随季节降雨影响大，因此径流量具有明显的季节性变化的特点。降水时空分布极不均匀，从时间上看，汛期6~9月降水量占年降水量的75%~81%；7、8两个月的降水量占年降水量的50%~60%，暴雨为主；从空间上来看，降水量南多北少、东多西少。窟野河流域多年平降雨量约为399.13 mm (1971年~2016年)，多年平均潜在蒸散发量为1273.66 mm (1971年~2016年)。

图1. 窟野河流域地理概况图

2.2. 数据来源与处理

本文选取窟野河流域的把口水文站温家川站(见表1)，年径流量数据来源于控制面积最大的温家川水文站(1971~2016年)的日流量数据。本文考虑到站点空间分布以及气象站资料的时间序列和完整性选择了窟野河流域周边的神木、兴县、河曲、东胜4个国家基本气象站的1971~2016年逐日气象资料。

表1. 水文站基本情况

3. 方法

3.1. Mann-Kendall趋势检验

Mann-Kendall趋势检验法，作为世界气象组织推荐并已经广泛应用的统计方法，能够有效的确定某一自然过程的变化趋势的显著性 [7]。

3.2. Pettitt突变点检验和累积距平曲线法

气候的变化和人类活动会影响到水文序列在某个时间点前后发生显著突变，从而导致水文数据出现显著的阶段性特征。利用Pettitt突变点检验和累积距平曲线法，在一定的显著水平或置信区间内，通过统计量的计算，判断趋势或突变变化，比较确定显著变化的年份 [8]。

3.3. 径流变化归因识别

根据Budyko水热平衡理论 [9]，建立了能够反映气候和人类活动造成地表变化对径流影响的分析方法。经多流域和长序列验证能够应用与窟野河流域，是一种准确适合的方法。

3.3.1. Budyko水热平衡理论

在一定的气候和植被条件下，流域长期的水文气候特征服从水分和能量平衡原理 [10]，也就是基于Budyko假设推导出的Choudhury-Yang公式的流域水热耦合平衡方程 [13]，表达式如下：

$ET=\frac{P\times E{T}_0}{{\left(P^n+E{T}_0^n\right)}^{1/n}}$ (1)

式中：ET₀为多年平均潜在蒸散发量(mm)；n为下垫面特征参数。

潜在蒸散发量ET₀采用FAO修正的Penman-Mnoteith公式 [11]，计算公式：

$E{T}_0=\frac{0.408\Delta \left(R_n-G\right)+\gamma \frac{900}{T+273}{U}_2\left(e_s-e_a\right)}{\Delta +\gamma \left(1+0.34U_2\right)}$ (2)

式中：ET₀潜在蒸散发量(mm/d)； $\Delta$ 为饱和水气压曲线斜率(kPa/℃)；R_n太阳净辐射(MJ m⁻²d⁻¹)；G为土壤热通量(MJ m⁻²d⁻¹)； $\gamma$ 为干湿常数(kPa/℃)；U₂为2 m高处的风(m/s)；T为平均气温(℃)；e_s为平均饱和水汽压(kPa)，e_a为实际饱和水汽压差(kPa)。

下垫面特征参数(n)是影响流域水平衡的一个重要直接因素，是Budyko曲线整体形状的特征参数参数 [12]。本文研究的流域中下垫面特征主要是指流域内的植被、地形、土壤等。

在长序列的水文过程中，我们可以认为 $P=ET+R$，故结合(2)式可以得出水量平衡公式以及他的全微分形式：

$R=P-\frac{P\times E{T}_0}{{\left(P^n+E{T}_0^n\right)}^{1/n}}$ (3)

$\text{d}R=\frac{\partial f}{\partial P}\text{d}P+\frac{\partial f}{\partial E{T}_0}\text{d}E{T}_0+\frac{\partial f}{\partial n}\text{d}n$ (4)

3.3.2. 气候因素和下垫面的弹性分析

在水文过程中，某气象因素的弹性指单位的该因素变化导致流域径流量的变化程度 [14] 如降水量增加1%导致年径流量相对于多年平均值的变化百分比。同样下垫面的单位变化导致径流的变化。在这里下垫面变化也就是通过下垫面特征参数n来反映出由于人类活动、自然原因导致的下垫面特征发生变化。由式(4)可知，P、ET₀、n都是对径流的特定独立变量。那么径流对独立变量x的弹性系数可以表达为 [15]：

${\epsilon }_{x_i}=\frac{\partial R}{\partial x_i}\times \frac{x_i}{R}$ (5)

式中 ${\epsilon }_{x_i}$ 是径流对特定独立变量 $x_i$ 的弹性系数， $x_i$ 可以表示为P、ET₀、n。

根据式(5)的定义，用式(4)除多年径流深R，可以得到：

$\frac{\text{d}R}{R}={\epsilon }_P\frac{d}{P}+{\epsilon }_{E{T}_0}\frac{\text{d}E{T}_0}{E{T}_0}+{\epsilon }_n\frac{\text{d}n}{n}$ (6)

令 $\varphi =\frac{E{T}_0}{P}$ 最终可得出各变量弹性系数计算如下：

${\epsilon }_P=\frac{{\left(1+{\varphi }^n\right)}^{1/n+1}-{\varphi }^{n+1}}{{\left(1+{\varphi }^n\right)}^{1/n}\left[{\left(1+{\varphi }^n\right)}^{1/n}-\varphi \right]}$ (7)

${\epsilon }_{E{T}_0}=\frac{1}{\left(1+{\varphi }^n\right)\left[1-{\left(1+{\varphi }^{-n}\right)}^{1/n}\right]}$ (8)

${\epsilon }_n=\frac{\ln \left(1+{\varphi }^n\right)+{\varphi }^n\ln \left(1+{\varphi }^{-n}\right)}{n\left(1+{\varphi }^n\right)\left[1-{\left(1+{\varphi }^{-n}\right)}^{1/n}\right]}$ (9)

根据式(7)~(9)可以求得气象因素和下垫面因素对径流的弹性系数。

3.3.3. 径流变化归因分析

根据突变点将研究时段分为两个子时段，时段1记为基准期，时段2记为径流变化期，时段1的径流深记为R¹，时段2的径流深记为R²，从时段1到时段2的多年平均径流深变化可表示为：

$\Delta R=R^2-R^{\text{1}}$ (10)

同理降水量变化、潜在蒸散发变化和下垫面变化可表示为

$\Delta P=P^2-P^1$ (11)

$\Delta E{T}_0=E{T}_0^2-E{T}_0^1$ (12)

$\Delta n=n^2-n^1$ (13)

由因子引起的径流变化可以用因子在两个时期的变化与它和径流偏导数的乘积来计算，因此，每个因子对径流的变化贡献可以通过公式表示为：

$\Delta R^{\prime }=\frac{\partial R}{\partial P}\Delta P+\frac{\partial R}{\partial E{T}_0}\Delta E{T}_0+\frac{\partial R}{\partial n}\Delta n$ (14)

$\Delta R^{\prime }$ 表示通过计算得到的径流深变化，也可以得到每个因素对径流变化的相对贡献dR_p、dR_ET0、dR_n。

每个因素对径流的变化相对贡献率可以计算如下：

$C_{X_i}=\frac{\text{d}{R}_{x_i}}{\text{d}{R}^{\prime }}\times 100\%$ (15)

式中： $x_i$ 表示P、ET₀、n， $C_{x_i}$ 表示各因子对径流变化的贡献率。

4. 结果与分析

4.1. 窟野河流域径流、降水及潜在蒸发变化趋势分析

图2为窟野河流域径流、降水、潜在蒸散发、下垫面系数1971~2016年的逐年变化情况从简单的线性拟合趋势，可以看出降水量和潜在蒸散发都有所增加但是并不明显，下垫面系数有明显的增加趋势，径流深有明显的下降趋势。

图2. 窟野河流域径流、降水、潜在蒸散发、下垫面系数1971~2016年的逐年变化

根据图3和图4可以看出在采用Mann-Kendall趋势检验时(显著性水平α = 0.05、临界线值为±1.96)，窟野河流域径流在研究时间序列内基本处于减少状态，在1972年~1976年以及1997年之后均处于显著减少(α = 0.05，UF统计量 < −1.96)，其余时间为波动性减少(统计量为负值)。而降水在1972年~1973年显著性减少(UF统计量 < −1.96)，其余时间的减少没有达到显著性水平，部分年份稍有增加。潜在蒸发变化趋势则是在1971年~1999年减少，其中1975年~1987年显著性减少(UF统计量 < −1.96)，而从1999年变为增加，分析原因可能是由于1997年后国家开始退耕还林还草政策的执行，植被覆盖度的增加导致了植物截留耗水，进一步导致潜在蒸散发的增加。

图3. 窟野河流域多年径流深Mann-Kendall趋势检验

图4. 窟野河流域多年降水、潜在蒸散发Mann-Kendall趋势检验

4.2. 窟野河径流量突变特征

窟野河流域径流量变化特征采用累计平距法和Pettitt检验共同确定，气候变化和人类活动可引起水文序列在某个时间点前后的统计规律发生显著突变，从图4中发现窟野河流域径流量在1997年发生了明显的变化，所以我们将1997年以前作为基准期。将研期分为两段1971~1996年(26年)、1997~2016(20年)。各时期水文变量特征值见表2。

Pettitt突变点检验的结果见图5左显示：窟野河流域径流深在1971年~2016年发生了两次突变，1980年为一级突变，1997年为二级突变，统计值U_T分别为249、488分别达到了置信度为95%、99%的显著性检验。累积距平的突变分析见图右：在研究时段内1979年上升趋势变缓，进入波动增加的阶段，1996年则是第二个突变点之后开始进入下降趋势。

图5. 窟野河流域径流深的Pettitt突变点检验和累积距平突变点识别

表2. 窟野河基准期和径流变化期水文要素特征值

注： ${\epsilon }_P$ 降水量弹性系数， ${\epsilon }_{E{T}_{\text{0}}}$ 潜在蒸散发的弹性系数， ${\epsilon }_n$ 下垫面特征参数弹性系数。

4.3. 窟野河径流量变化的归因分析

4.3.1. 径流对气候因素和下垫面变化的敏感性分析

上表2中列出窟野河流域各时期水文气候数据以及下垫面特征参数值及径流对气候因素和下垫面参数的弹性系数。流域潜在蒸散发量较基准期增加，降水略有增加，径流深较基准期明显减少；

干旱指数(ET₀/P)较基准期增大但并不显著，可见气候在两个时段内的变化并不明显，径流系数(R/P)较基准期减小超过一倍，降水的产流机制在人类活动期受到了很大的影响，雨水下渗或被直接拦蓄而不是产生径流。径流与ET₀、n呈负相关，但与P呈正相关。ET₀的弹性系数范围为−1.950~−1.606、P为2.606~2.95，ω为−3.222~−2.846。这些范围表明，ET₀增加1%将导致1.606%~1.950%的下降，P增加1%将导致2.606%~2.95%的增加以及n增加1%将导致2.846%~3.222%的径流减少。降水量偏少，径流系数小，蒸散发量与降水量的比值指数达到半干旱带标准，多年平均径流深约为46.0 mm，符合黄河中游的气候条件。

图6为1971~2016年窟野河径流对气候因素和下垫面因素的弹性系数趋势，反映了降水、潜在蒸散发量和下垫面特征参数对径流的影响，由图可知窟野河流域降水的弹性系数处于波动性增加但未达到显著性水平和潜在蒸散发的弹性系数的绝对值略有增加，而下垫面特征系数的弹性系数在1997开始绝对值有显著增加(P < 0.05)这也表明了窟野河流域径流深人类活动影响期对下垫面情况的变化更加敏感，而对气候变化敏感性相对较低。

图6. 1971~2016年窟野河径流对气候和下垫面因素的弹性系数趋势

4.3.2. 窟野河径流变化归因识别

降水、潜在蒸散发和下垫面特征条件的变化对径流变化的贡献率见下表3，其中 $\delta$ 值很小，表示在计算各要素对径流变化贡献率的方法正确有效。

降水的减少，潜在蒸散发的增加或是下垫面特征参数的增加会引起径流量的减少，反之则会引起径流量的增加。由表3的结果可以看出与其他研究不同的是再今年的降雨数据加入后人类活动期对比基准期降雨开始稍有增加，潜在蒸散发和下垫面参数都处于增加的状态；窟野河流域的径流变化主导因素是人类活动导致下垫面条件的变化对径流的贡献率达到了98%，气候因素方面潜在蒸散发对径流的影响大于降水所带来的径流变化。

综上，本文研究区流域径流量减少的主要因素就是下垫面的变化，其次是潜在蒸散发，降水量变化较小。

表3. 窟野河径流变化归因识别

dR_P：降水量P引起的径流变化；dR_ET0：潜在蒸散发量ET₀引起的径流变化；dR_n：下垫面特征参数n引起的径流变化；dR：人类活动期径流深与基准期之差；dR'：通过计算求得的径流深变化； $\delta$ ：dR'和dR的差；C_p：降水变化对流域径流变化的贡献率；C_ET0潜在蒸散发对径流变化的贡献率；C_n下垫面变化对径流变化的贡献率。

5. 讨论

窟野河流域径流深的急剧减少受到气候因素和人为地形变化因素影响。通过分析流域内4个国家气象站得出气象方面降水和潜在蒸散发均呈非显著性的变化，表明窟野河流域的气候在近45年没有显著变化，对窟野河流域的径流变化的贡献占比仅为−1.64%~2.76%。

而在下垫面因素方面，从表4、图7可以看出：由于国家开展退耕还林还草的政策影响，窟野河流域内耕地面积不断减少，林地草地的面积不断增加，林地草地对窟野河流域内的水土流失影响显著。通过植被的固土作用防止水土流失，减少径流中的含沙量，但同时也导致窟野河就流量减少，增加了流域内的蒸散发量，并且值得注意的是窟野河流域内的煤矿用地在1996年之后激增约200倍，且大部分分布流域中北部的河流沿岸和水域附近，在开采煤矿的时候会造成地下形成采空区，导致地表径流向地下，填补采空区使径流量急剧减少。并且根据Budyko假设认为长时间序列流域储水量 $\Delta S$ 为0，但是由于1980年以来窟野河流域内修建闸道、淤地坝等，对径流起到了拦蓄的作用，1996年之后的径流向煤矿采空区下渗，都会影响到 $\Delta S$ 并不能完全忽略，在未来的计算中应该考虑引入参数来减少或避免这类水保工程以及径流损失所带来的误差影响。

表4. 窟野河流域土地利用类型及面积(km²)

图7. 1980年和2015年窟野河流域土地利用类型

6. 结论

本文采用Mann-Kendall趋势检验、Pettitt突变点检验和累积距平曲线突变点检验法对窟野河流域内气象数据处理，应用budyko水热平衡理论对窟野河流域径流变化归因识别，得出了人类活动及气候因素对窟野河径流量变化的贡献，结论如下：

1) 在45年内，窟野河流域的降水无显著变化，稍有增加。潜在蒸散发呈非显著性增加，径流则在1971~1996呈非显著性减少，97年后呈显著性减少。他们增加(或减少)的平均速率为1.22 mm/a、0.77 mm/a、−1.42 mm/a。

2) 采用两种突变点检验中得出1996年后是人类活动影响的时期，将研究时段分为了1971~1996年的基准期和1997~2016年的人类活动期，并分别计算得出了人类活动期降水、潜在蒸散发和下垫面特征值对窟野河径流变化的贡献值分别为：−1.64%、2.76%、98.88%，可以看出下垫面的变化对窟野河径流减少起到非常显著的影响，接下来是潜在蒸散发和降雨带来的影响。

基金项目

国家自然科学基金51779209；陕西水利科技计划项目(2019slkj-20)。

文章引用

郑嘉昊,蒋晓辉. 基于Budyko水热平衡理论的窟野河径流变化归因识别
Attribution Identification of Runoff Variation in Kuye River Based on Budyko’s Theory of Water and Heat Balance[J]. 水资源研究, 2020, 09(05): 471-482. https://doi.org/10.12677/JWRR.2020.95050

参考文献