本文主要收集了瀑河流域水文气象站平泉站(图1)1981— 2020年月均径流深的实测数据作为水文数据, 实测数据较为完整, 且研究时间序列较长, 可以满足研究的需要。

此外, 本文还收集了瀑河流域水文气象站平泉站以及雨量站大庙站及沙坨子站(图1)1981— 2020年月均降水量数据以及平泉站1984— 2013年月均蒸散发数据。研究区的面降水量用3个站点的算术平均值表示, 3个站点均满足研究的需要, 且各站点雨量数据无明显的突变点和随机变化, 数据变化较为一致, 数据可靠, 可以代表区域气候状况。遥感数据均为在地理空间数据云平台(https://www.gscloud.cn)和空天院遥感数据服务系统(http://eds.ceode.ac.cn/nuds/freedataquery)网站下载的Landsat系列卫星影像数据。

2.2.1 趋势分析方法

M-K法是一种非参数检验方法, 因其不受样本值和样本分布类型的影响, 通常用于水文、气象和序列的其他非正态分布。本文利用M-K法分析瀑河上游年均径流深序列的突变点及变化趋势^[6]。

重标度极差分析法(Rescaled Range Analysis, R/S法)^[13]是由英国水文学家Hurst在研究尼罗河水坝工程时提出的方法, 后来被用在各种时间序列分析中, 主要用来预测时间序列的变化趋势。Hurst指数的取值范围为0~1, 当H∈ [0, 0.5)时, 表明时间序列的变化具有反持续性, 即未来变化趋势与过去相反, 且H值越接近于0, 反持续性就越强; 当H=0.5时, 序列表现为独立同分布的随机序列; 当H∈ (0.5, 1]时, 表明时间序列具有长程相关性, 即序列在未来的变化趋势与过去相同, 且H值越接近1, 持续性就越强。

2.2.2 基于Budyko假设的气候弹性系数法

根据全球水和能源平衡, Budyko指出蒸散发受降水和蒸发能力的控制, 中国气象学家傅抱璞^[14]提出了实际蒸散发的计算公式

E=E01+rE0-1+rE0ω1/ω 。(1)

式中: E为实际蒸散发, mm; E₀为潜在蒸散发(最大可能蒸散发), mm; r为降水量, mm; ω 为下垫面特性参数, 主要反映下垫面条件, 如地形、植被等条件。

根据傅抱璞^[14]的计算, 山地和高原地区下垫面特性参数ω 计算式为

ω=1+0.2928r1.2101-YY0.393 。(2)

式中: Y为径流系数; 其余符号含义同式(1)。

通过计算, 本研究区ω =1.29。

Li等^{[16, 17]}在式(1)的基础上进行了改进, 计算公式为

EP=1+ω(E0/P)1+ω(E0/P)+(E0/P)-1 。(3)

式中: P为多年平均降水量, mm; 其余符号含义同式(1)。

对于某个封闭区域, 其多年平均降水量的表达式为

P=E+Q+ΔS 。(4)

式中: Q为年均径流深, mm; Δ S为流域蓄水变量, m³, 在平均规模的年份中, 可以忽略^[15]; 其余符号含义同式(3)。

根据Li等^{[16, 17]}提出的方法, 气候变化和人类活动引起的径流变化的计算公式为

ΔQP=∂Q∂PΔP ,(5)

ΔQE0=∂Q∂E0ΔE0 ,(6)

ΔQh=ΔQ-ΔQc 。(7)

式中: Δ Q_P为年均降水量的变化引起年均径流深的变化量, mm; Δ QE0为年均潜在蒸散发的变化引起年均径流深的变化量, mm; Δ Q_h为人类活动引起年均径流深的变化量, mm; Δ P和Δ E₀分别为年均降水量和年均潜在蒸散发的变化, mm; Δ Q为年均径流深的变化量, mm; Δ Q_c为气候变化引起年均径流深的变化量, mm; ∂Q∂P和 ∂Q∂E0分别为降水弹性系数和潜在蒸散发弹性系数^[18]。

∂Q∂P=1+2x+3ωx(1+x+ωx2)2 ,(8)

∂Q∂E0=-1+2ωx(1+x+ωx2)2 。(9)

式中: x为干旱指数, 为年均潜在蒸散发与年均降水量的比值; 其余符号含义同上。

气候变化和人类活动对年均径流深变化的贡献率计算公式为

ηc=|ΔQc||ΔQc|+|ΔQh|×100% ,(10)

ηh=1-ηc 。(11)

式中: η _c为气候变化对年均径流深变化贡献率; η _h为人类活动对年均径流深变化贡献率; 其余符号含义同上。

据M-K法, UF指的是正序列趋势统计量, UF> 1.96表明时间序列有着显著的增加趋势, UF< 1.96表明时间序列有着显著的减少趋势, UF> 0对应年均径流深呈增加趋势, UF< 0对应年均径流深呈减少的趋势。由图2可知, 年均径流深在1984— 1999年间有着显著的上升趋势, 2001年后均处于下降趋势。UB指的是逆序列趋势统计量, UF与UB的交点为突变点, 发生在2006年(图2)。依据M-K法检验结果, 将1981— 2006年定为基准期, 2007— 2020年定位变化期。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 瀑河上游1981— 2020年年均径流深M-K检验图Fig.2 M-K trend test of annual average runoff depth in the upper reaches of the Bao River from 1981 to 2020

通过点绘Log(R/S)散点图(图3), 得知瀑河上游的H值为0.4026(接近0.5), 根据2.2.1中Hurst指数取值范围含义, 表明未来年均径流深与1981— 2020年间年均径流深的变化趋势相反但会出现随机分布情况。1981— 2020年均径流深呈下降趋势(图4), 即未来年均径流深呈相对平稳的增加趋势。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 R/S法散点图Fig.3 Scatter diagram of R/S analysis method

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 瀑河上游1981— 2020年均降水量、年均潜在蒸散发、年均径流深的趋势性检验Fig.4 Trend test of annual average precipitation, potential evapotran-spiration and runoff depth in the upper reaches of the Bao River from 1981 to 2020

气候变化对瀑河上游径流的影响不容小觑, 主要包括降水量和潜在蒸散发两方面因素。气候变化中降水的弹性系数最大, 是影响径流变化的主要因素。瀑河上游1981— 2020年平均降水量为515.6 mm, 最大降水发生在1990年, 为797.8 mm, 最小降水量发生在2000年, 为333 mm, 年际变化较大(图4)。2007— 2020年年平均降水量较1981— 2006年减少了38.71 mm, 降水改变导致径流深降低了33.13 mm, 贡献率为71.31 %, 降水条件的变化直接影响径流的变化趋势。

另外, 有关潜在蒸散发, 现今也有较多相关研究: 王庆明等 ^[19]在海河流域山区植被覆盖度与水热因子关系的研究中提出, 植被覆盖度与潜在蒸散发呈正相关性, 说明植被覆盖度的降低会减少潜在蒸散发; 赵勇等^[20]通过Choudhury-Yang公式计算得出降水量的降低与潜在蒸散发减少呈正相关性。本研究区内同时存在着降水量减少与植被覆盖度增加两种因素, 在这两种因素共同作用下导致了潜在蒸散发的减少。瀑河上游1981— 2020年平均潜在蒸散发为380.29 mm, 最大年均潜在蒸散发发生在2002年, 为426.97 mm, 最小年均潜在蒸散发发生在2019年, 为313.98 mm, 年际变化较大(图4)。2007— 2020年年均蒸散发较1981— 2006年减少了31.16 mm, 潜在蒸散发的减少使径流深增加了 14.41 mm, 贡献率为-31.01%, 对径流产生积极的影响。

本研究中人类活动主要指土地利用类型的变化。遥感解译选取瀑河上游1985— 2020年8期土地利用类型数据, 结果统计可知(图5): 研究区乔木林地面积呈现上升趋势, 2020年较1985年增长约67.61%, 2020年乔木林地面积占流域总面积的51.26%; 草地面积呈现下降趋势, 2020年较1985年下降约52.16%, 2020年草地面积占流域总面积的13.92%; 耕地面积呈现下降趋势, 2020年较1985年下降约36.94%, 2020年耕地面积占流域总面积的28.51%; 建设用地、灌木林地及水域面积较小, 共占流域总面积的6.21%, 对径流的影响较小。整体来看, 乔木林地面积和草地、耕地面积存在较为明显的负相关关系, 说明耕地、草地面积与乔木林地面积存在相互转化的关系。总体来看, 研究区植被覆盖率高, 生态稳定性良好。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 瀑河上游8期遥感解译土地利用类型柱状图Fig.5 Remote sensing interpretation of the land use type histogram for the eight phases of the upper reaches of the Bao River

由瀑河上游1985— 2020年土地利用类型变化转移矩阵(表1)可知, 流域内土地利用类型变化较为活跃, 共有30.28%的土地发生了相互转化。由图6可知, 乔木林地增长的67.61%中, 有48.63% 都是由草地转化而来, 剩余18.98%由耕地转化。根据王勇等^[21]研究, 乔木较草地有更好的水源涵养性能。研究区水源涵养能力较差的耕地及草地大量转换为水源涵养能力较强的乔木林地, 使得乔木林地面积占到了总研究区面积的一半, 证明研究区在1985— 2019年间的水源涵养能力有所增加, 水分利用效益的提高会导致径流深的减少。

表1

Tab.1

表1(Tab.1)

表1 瀑河上游1985— 2020年土地利用类型变化转移矩阵 Tab.1 Land use type transformation matrix of the upper reaches of the Bao River from 1985 to 2020 地表覆盖层利用类型 2020年各土地利用类型面积/km2 耕地 乔木林地 灌木林地 草地 建设用地 裸地 水域 总计 1985年各土地利用 类型面积/km2 耕地 92.80 20.27 0.00 8.91 14.21 0.00 0.34 136.53 乔木林地 0.00 106.84 0.00 0.00 0.00 0.00 0.00 106.84 灌木林地 0.00 0.00 0.29 0.00 0.00 0.00 0.00 0.29 草地 6.66 51.96 0.05 39.72 3.08 0.05 0.14 101.66 建设用地 0.00 0.00 0.00 0.00 3.45 0.00 0.00 3.45 裸地 0.00 0.00 0.00 0.003 0.001 0.00 0.00 0.004 水域 0.13 0.003 0.00 0.003 0.01 0.00 0.43 0.576 总计 99.59 179.073 0.34 48.636 20.751 0.05 0.91 349.35

表1 瀑河上游1985— 2020年土地利用类型变化转移矩阵 Tab.1 Land use type transformation matrix of the upper reaches of the Bao River from 1985 to 2020

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 瀑河上游1985— 2020年土地利用转化情况桑基图Fig.6 Sankey diagram of the land use transformation in the upper reaches of the Bao River from 1985 to 2020

通常降水径流双累积曲线应呈现一条直线的趋势, 只有当存在其他因素影响时, 径流深的累积量才会发生偏移。由图7中可知, 与基准期相比, 变化期累积径流量呈减少状态, 说明除气象因素外, 人类活动也在发生着较大的变化。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 瀑河上游降水-径流双累积曲线图Fig.7 Double accumulation curve of precipitation and runoff for the upper reaches of the Bao River

通过第2节公式, 计算得到干旱指数x等(表2)。由表2可知, 与基准期相比, 变化期年均潜在蒸散发减少了31.16 mm, 年均降水量减少了 38.71 mm, 年均径流深减少了46.47 mm, 约占基准期径流深总量的77.10%, 干旱指数减少了0.01。

表2

Tab.2

表2(Tab.2)

表2 瀑河上游1981-2020年间不同时间周期的参数信息 Tab.2 Parameters for different time periods of the upper reaches of the Bao River from 1981 to 2020 时间 年均潜在蒸 散发/mm 年均降水量/ mm 年均径流深/ mm 干旱 指数 基准期 410.33 529.10 60.27 0.78 变化期 379.17 490.39 13.80 0.77

表2 瀑河上游1981-2020年间不同时间周期的参数信息 Tab.2 Parameters for different time periods of the upper reaches of the Bao River from 1981 to 2020

通过式(8)和式(9), 分别得到了降水和潜在蒸散发的弹性系数, 通过式(10)和式(11), 得到气候变化和人类活动对年均径流深变化贡献率(表3)。如表3所示: 年均降水量的减少使年均径流深减少了33.13 mm, 贡献率为71.31%, 但年均潜在蒸散发的减少使年均径流深增加了14.41 mm, 贡献率为-31.01%, 综合来看, 气候变化对径流影响贡献率为40.30%; 人类活动使年均径流深减少27.74 mm, 贡献率为59.71%。在陈鑫^[12]对滦河流域近60 a径流变化归因分析中, 人类活动贡献率为59%, 人类活动为主导因素。本文计算结果与之相对应, 人类活动是径流深减少的主要因素。

表3

Tab.3

表3(Tab.3)

表3 基于Budyko假设的瀑河上游1981— 2020年间气候变化和人类活动对径流的贡献 Tab.3 Contribution of climate change and human activities to runoff depth in the upper reaches of the Bao River from 1981 to 2020 based on the Budyko hypothesis 项目 ∂ Q/∂ P ∂ Q/∂ E0 Δ QP/mm η P/% Δ QE0/mm η E0/% Δ Qh/mm η h/% 敏感系数及贡献率 0.86 -0.46 -33.13 71.31 14.41 -31.01 -27.74 59.71 注: η P表示年均降水量对年均径流深的贡献率; η E0表示年均蒸散发对年均径流深的贡献率; 其余符号含义同式(1)— (11)。

表3 基于Budyko假设的瀑河上游1981— 2020年间气候变化和人类活动对径流的贡献 Tab.3 Contribution of climate change and human activities to runoff depth in the upper reaches of the Bao River from 1981 to 2020 based on the Budyko hypothesis

基于Budyko假设的气候弹性系数法分析了气候变化和人类活动对径流深的影响, 人类活动是径流深减少的主要原因。研究区人类活动主要体现在土地利用方式的改变, 即对下垫面的改变。由8期遥感解译结果得知, 1985— 2020年, 瀑河上游区域土地利用类型的主要变化为耕地、草地转换为乔木林地, 主要受退耕还林、矿山修复等水土保持措施影响。由于耕地、草地的植被顶冠较小, 根系发育程度也不及乔木林地发达, 持水能力较差, 因此土地利用类型发生改变后, 大大增强了研究区的水源涵养能力, 更多的水资源得以留存下来导致地表径流深减少。近年来, 很多研究也表明, 受退耕还林等政策影响, 乔木林地的增加是我国诸多流域径流深减小的主要原因^{[22, 23, 24]}。在后期的工作中, 可通过水文模型的建立, 细化造成的影响是积极的还是消极的, 量化水源涵养功能的变化趋势, 为京津冀水源涵养功能区的保护修复提供数据支撑, 同时为地方政府在水资源配置等工作中提供数据支撑。