1957-2016年祖厉河上游降水与水沙变化趋势

赵秀兰, 周 蕊, 张 富, 胡彦婷, 李晓雅, 包炳琛, 蒋承洋

(1.甘肃农业大学 林学院, 兰州 730070; 2.甘肃省水土保持科学研究所, 兰州 730000)

径流、输沙是河川径流重要的水文特征[1]，是衡量区域生态环境变化的重要指标。随着全球气候变化和人类活动对自然环境的影响加剧，陆地水文循环的时空格局也发生了改变。在自然状态下，降水是半干旱黄土区水资源的主要来源，对河川径流具有重要影响。径流作为水蚀的动力，其大小与输沙变化具有密切的关系。对降水与径流、输沙变化趋势及特征的研究，是做好水资源管理及流域生态建设的基础性工作。

21世纪以来，受气候变化及人类活动的影响，黄河水沙特性发生巨大的变化，对当地社会的可持续发展造成严重影响[2-3]，水沙变化原因及规律成为研究热点。赵阳等[4]对黄河流域水沙变化研究表明，在1950—2016年期间径流量和输沙量均呈显著减少趋势，径流量和输沙量的突变点分别发生在1981年、1999年，突变后减幅分别为17.93%～40.79%和46.74%。孙洋洋[5]、李二辉[6]等对黄河上中游径流量变化进行了分析，表明径流量在1919—2013年呈明显的枯(1919—1932年)—丰(1933—1985年)—枯(1986—2013年)的阶段性变化，且在1985年以后呈显著减少趋势。赵广举等[7]对黄河中游1950—2009年水沙变化进行研究，结果表明径流量和输沙量均呈减少趋势，且在1970年之后减少显著。潘彬等[8]对黄河下游径流量变化研究表明，1962—2012年降水量和径流量呈减少趋势，降水减少趋势并不显著，降水量、径流量突变时间分别为1964年和1979年；降水量和人为因素对径流影响作用分别达到32%和68%。王蕊等[9]对黄河上游湟水河小南川流域气候和土地利用变化对水沙运移关系的研究表明1963—2013年降水量逐渐增加，径流量和输沙量逐渐减少。Zhang等[10]对祖厉河流域1955—2013年降水量、径流量和输沙量变化特征进行研究，结果表明降水量、径流量和输沙量均逐年减少，降水量和径流量在1995年出现突变点，输沙量在2000年出现突变点。以上研究表明，区域尺度不同，不同时段水沙变化规律也不相同。同时研究多在较大尺度区域，而对于中小尺度区域水沙变化规律及周期变化特征研究较少。本研究基于祖厉河上游会宁水文站实测径流量、输沙量资料及雨量站降水资料，采用M-K突变检验和小波分析方法，对降水量、径流量和输沙量的变化特征、周期分布规律进行研究，以期为中小尺度区域合理开发保护流域水资源、水土流失治理及区域经济发展提供科学依据。

1 材料及方法

1.1 研究区概况

祖厉河是黄河一级支流，位于黄土高原丘陵沟壑区第五副区，总面积10 653 km2。在黄河上游面积超过1万km2的6条支流中，祖厉河年均径流量1.086亿m3，泥沙量为5 320万t，分别占入黄径流量、输沙量的1.60%，39.6%，是黄河上游泥沙的主要来源。研究区位于祖厉河上游甘肃省会宁县南部，地理位置为104°56′—105°23′E，35°22′—35°47′N，包括会宁、通渭两县的中川、新添、太平、华家岭(祖厉河的发源地)等10个乡镇，总流域面积983.86 km2(图1)。研究区丘陵起伏，沟壑纵横，地形破碎，地表被黄土覆盖，地势由南向北走向，平均海拔在2 025 m。气候属半干旱区，年均气温7.9℃，最高气温37.5℃，最低气温-26.5℃，降水较少且分布不平衡，年均降水量为406.4 mm，年平均蒸发量为1 800 mm。土壤主要是黄绵土、灰钙土和灰褐土为主，局部低洼地分布少量盐碱土。自然植被主要以冰草、针茅、蒿类等草本植物为主，粮食作物主要有小麦、马铃薯、玉米等。干旱是区域内主要的自然灾害，也伴有霜冻、冰雹、洪灾、病虫害等。土地利用类型以耕地、草地和灌木林地为主。据会宁水文站1957—2016年观测资料，研究区年均径流量为1 227.4万m3，年均输沙量为367.6万t。

图1 研究区地理位置

1.2 数据来源及方法

1.2.1 数据来源 本研究采用甘肃省水文水资源局提供的1957—2016年会宁水文站逐月径流量、输沙量及10个雨量站的降水观测数据。

1.2.2 研究方法 使用泰森多边形法对降水量观测数据进行预处理，计算平均降水量，通过Excel软件对径流量、输沙量原始资料进行整理。采用M-K趋势分析和突变检验法[11]，对研究区1957—2016年降水量、径流量和输沙量的变化趋势及突变进行研究；选取Morlet小波分析方法对降水量、径流量及输沙量时间尺度演变特性进行分析[12-14]，进行连续小波变换，得到小波系数的实部、方差信息，并对变化的主周期通过小波方差进行检验，计算平均周期；利用SPSS 20软件，通过线性回归及曲线估计方法，分析降水量、水土保持措施量与径流量、输沙量相互间的关系，建立关系式。

2 结果与分析

2.1 降水量、径流量及输沙量变化特征

2.1.1 降水量变化特征分析 根据气象资料统计，研究区60 a(1957—2016年)年均降水量为406.4 mm，其中1967年的降水量最大(607.2 mm)，1982年降水量最小(243.9 mm)，两者相差363.3 mm，最大降水量是最小降水量的2.49倍。1957—2016年降水量变化趋势分析表明，年降水量略有降低的趋势，但总体未达到显著水平(p>0.05)(图2A)。

图2 降水量变化和M-K检验

由图2B可以看出，在1957—1974年(18 a)降水量有明显增加趋势，平均降水量414.50 mm，其中1963—1968年(6 a)UFk值达到显著水平(p<0.05)，平均降水量499.20 mm；1975—1990年(16 a)平均降水量410.13 mm，在均值附近波动；1991—2016年(26 a)，降水量呈减少趋势，平均降水量398.51 mm；1975—2016年降雨量无显著变化(p>0.05)。

2.1.2 径流量、输沙量变化特征

(1) 径流量变化趋势分析。根据观测数据，研究期年均径流量为1 227.4万m3，年径流量总体出现波动下降的趋势(见图3A)。

其中1959年的径流量最大(4 549.0万m3)，2016年径流量最小(3.0万m3)，最大年径流量是最小年径流量的1 516倍。

图3 年径流量变化和M-K检验

由图3B可知，在1957—1974年(18 a)径流量呈增加趋势，平均径流量2 131万m3；在1975—1999年(26 a)径流量呈减少趋势，其中1993—1999年(7 a)达到显著水平(p<0.05)，平均径流量772万m3；在2000—2016年(17 a)，径流量呈显著减少趋势(p<0.05)，平均径流量417万m3；同时UFk与UBk两条曲线在2000年出现了交点，说明突变点出现在2000年。

(2) 输沙量变化趋势分析。年输沙量变化趋势与径流量相似，总体出现波动下降的趋势(图4A)。研究期年均输沙量为367.6万t，其中1973年的输沙量最大(1 614.5万t)，2016年输沙量最小(1.3万t)，最大年输沙量是最小年的1 238倍。年际间输沙量存在较为明显的波动。

由图4B可看出，在1957—1974年(18 a)输沙量呈增加趋势，平均输沙量619.9万t；在1975—2002年(28 a)输沙量呈减少趋势，其中1996—2002年(7 a)输沙量达到显著水平(p<0.05)，平均输沙量180.5万t；在2003—2016年(14 a)输沙量呈显著减少趋势(p<0.05)，平均输沙量123.3万t；同时UFk与UBk两条曲线的相交于2003点，说明突变点出现在2003年。

2.2 降水量、径流量和输沙量周期变化特征分析

2.2.1 降水量周期变化特征

(1) 时间尺度周期性分析。由图5A可知，年降水量存在3种尺度的周期变化：4～8 a，11～15 a和16～32 a，其中16～32 a和4～8 a周期内的“丰—枯”交替变化较明显，贯穿整个时间序列。

小尺度4～8 a的周期性变化，主要以“枯—丰”周期变换为主，出现了13次震荡；中尺度11～15 a的周期性变化在1990年以前表现较明显，出现“丰—枯”交替的4次振荡；而大尺度16～32 a周期中，降水量出现“枯—丰”交替的4次震荡，突变特性明显，且从图5A还可以看出，整个大尺度的周期变化占据了整个时间序列且状态比较稳定，具有全域性。

图4 年输沙量变化和M-K检验

注：图中黑色粗实线是小波系数零点等值线，代表突变线；细实线是小波系数正值等值线，图中浅色区域，代表降水量和径流量偏丰(输沙量多沙期)；虚线是小波系数负值等值线，图中深色区域，代表降水量和径流量偏枯(输沙量少沙期)。

图5 年降水量小波系数实部等值线、小波方差和小波周期

(2) 不同周期变化分析。图5B中存在4个峰值，分别对应5 a，8 a，22 a，27 a时间尺度，其中，27 a的周期振幅最大，为降水序列变化的第一主周期，则22 a，8 a，5 a为第2，3，4主周期。由图5C可知，降水序列在27 a时间尺度下，经历约3个波动周期，其平均变化周期约为18 a；降水量“枯—丰”的转变点在1963年、1981年和1999年。从周期变化可以预测降水量在27 a时间尺度(2017年左右)将由枯转丰。在22 a时间尺度下，约经历4个波动周期，其平均变化周期约为13 a。从变化趋势可以预测降水量在22 a时间尺度(2028年左右)将由枯转丰。降水序列在8 a和5 a时间尺度下，分别经历12个和18个周期波动，其平均变化周期分别约为5 a和3 a。从变化趋势来看，预测在8 a时间尺度下，2020年左右降水量由枯转丰，而5 a时间尺度下，2019年降水量由枯转丰。

2.2.2 径流量周期变化特征

(1) 时间尺度周期性分析。从图6A可以看出，年径流量存在3种尺度的周期变化：3～7 a，8～12 a和14～32 a，这3种尺度下径流量在整个时间序列丰、枯水期交替变化明显。小尺度3～7 a 的周期性变化主要在1957—2010年表现得较为活跃，存在“枯—丰”交替的周期性变化。8～11 a尺度上的周期性变化主要在1957—2005年变化明显，存在“枯—丰”交替的8次震荡；在1985—2005年，14～24 a时间尺度上，径流量呈现“枯—丰—枯”的状态，但在14～32 a的尺度上来说，是处于枯水期；而大尺度14～32 a来看，随着时间序列的变化，径流表现出了“枯—丰”交替的3次震荡，具有明显的突变性，具体时间为1957—1966年为枯水期，1966—1975年为丰水期，1975—1984年为枯水期，1984—1993年为丰水期，1993—2002年为枯水期，2002—2012年为丰水期，2012—2016年为枯水期，整个大尺度的周期变化占据了整个时间序列且状态比较稳定，具有全域性。

图6 年径流量小波系数实部等值线、小波方差和小波周期

(2) 不同周期变化分析。图6B中有3个峰值，分别对应4 a，9 a，27 a时间尺度，其中，4 a的周期振荡最强，为径流序列变化的第一主周期，则27 a，9 a为第2、3主周期。由图6C可知，在27 a时间尺度下，径流序列经历约3个波动周期，其平均变化周期约为18 a，径流量“枯—丰”的转变点在1965年、1982年和2003年。1969年处于27 a尺度下的偏丰极大值处。从周期变化可以预测径流量在27 a时间尺度(2021年左右)将由枯变丰；径流量在9 a和4 a时间尺度下，分别经历10个和22个周期波动，其平均变化的周期分别约为6 a和3 a。从变化周期预测在9 a时间尺度下，2020年左右的径流量由枯变为丰，而4 a时间尺度下，2018年左右的径流量由枯转丰。

2.2.3 输沙量周期变化特征

(1) 时间尺度周期性分析。年输沙量存在3种尺度的周期变化：3～7 a，8～12 a和14～32 a，这3种尺度下输沙量在整个时间序列多、少沙期交替变化明显(图7A)。在小尺度3～7 a的周期，主要在1957—2003年表现得较为活跃，存在“多—少”交替的周期性变化。8～11 a尺度上的周期变化主要在1957—1997年变化明显，存在“少—多”交替变化的7次震荡；在1972—1990年、14～24 a时间尺度上，输沙量呈现“少—多—少”的状态，但在14～32 a的尺度上来说，是处于少沙期；而大尺度14～32 a来看，随着时间序列的变化，输沙量表现出了“少—多”交替变化的3次震荡，存在明显的突变特性，具体时间为1957—1966年为少沙期，1966—1976年为多沙期，1976—1987年为少沙期，1987—1997年为多沙期，1996—2005年为少沙期，2005—2013年为多沙期，2013—2016年为少沙期，整个大尺度的周期变化占据了整个时间序列且状态比较稳定，具有全域性。

图7 年输沙量小波系数实部等值线、小波方差和小波周期

(2) 不同周期变化分析。图7B中有4个峰值，分别对应4 a，6 a，8 a，27 a时间尺度，其中，6 a左右的周期振荡最强，为输沙序列变化的第一主周期，则4 a，8 a，27 a为第2、3、4主周期。由图7C可知，在27 a时间尺度下，输沙序列经历约3个波动周期，其平均变化周期约为18 a；输沙量“少—多”的转变点在1968年、1987年和2005年。从周期变化可以预测输沙量在27 a时间尺度(2021年左右)将由少变多；在8 a的时间尺度下，输沙量经历12个周期波动；输沙量在6 a和4 a尺度下，分别经历15个和22个周期波动，其平均变化周期分别约为4 a和3 a。从变化周期来看，6 a和3 a时间尺度下，预测输沙量在6 a时间尺度下，2020年左右由少变多，而3 a时间尺度下，2018年左右输沙量由少变多。

对年降水量、径流量和输沙量比较可知，尽管径流量、输沙量和降水量的周期变化并不完全一致，但很相似，存在包含和部分包含的关系，如降水量的16～32 a的时间尺度包含于径流量和输沙量的大尺度14～32 a中，但径流量和输沙量的3～7 a和降水的4～8 a小尺度与只是部分包含。而径流量和输沙量的变化周期尺度一致，具有同步性。同时，从图中还可以发现，无论是降水量还是径流量、输沙量，都有一个大尺度下的丰水期(多沙期)或者枯水期(少沙期)，存在小尺度下的丰—枯水期(多—少沙期)嵌套的现象，小尺度下的降水量和输沙量、径流量转变点要多于大尺度，且不同尺度下的转变点时间及个数都不相同。

2.3 影响水沙量变化的主要因素

2.3.1 降水量与径流量、输沙量的关系 影响径流量、输沙量的主要因素是降雨和以水土保持措施为主的人类活动[15]。按照时间周期将降水量与径流量、输沙量的时间序列划分为4个时间段，第Ⅰ时段1957—1974年，第Ⅱ时段1975—1992年，第Ⅲ时段1993—2010年和第Ⅳ时段2011—2016年。不同时段降水量、径流量和输沙量平均值见表1。由回归分析可知(表2)，在第Ⅰ时段径流量、输沙量对降水量变化响应强烈，趋势协同性强，相关系数较高，降水量和径流量达到显著水平(p<0.05)，降水量对径流量的决定系数R2为0.482，降水量和输沙量未达到显著水平(p>0.05)，降水量对输沙量的决定系数R2为0.385，平均径流模数和输沙模数分别为2.17万m3/(km2·a)和0.63万t/(km2·a)；第Ⅱ时段径流量、输沙量对降水量变化响应减弱，趋势协同性弱化，相关系数降低，降水量和径流量达到显著水平(p<0.05)，降水量对径流量的决定系数R2为0.398，降水量和输沙量未达到显著水平(p>0.05)，降水量对输沙量的决定系数R2为0.293，平均径流模数和输沙模数分别为1.29万m3/(km2·a)和0.39万t/(km2·a)；第Ⅲ时段径流量、输沙量对降水变化的响应迅速衰减，趋势协同异化，降水量与径流量、输沙量均未达到显著水平(p>0.05)，降水量对径流量、输沙量的决定系数R2分别达到为0.114，0.012，平均径流模数和输沙模数分别为0.65万m3/(km2·a)和0.21万t/(km2·a)；说明在1957—2010年随时间延长，降雨对径流量和输沙量的影响逐渐降低；与其相比，径流量与输沙量相关关系趋势协同性强，略有减弱，但不明显，径流量对输沙量的决定系数R2分别为0.920，0.874，0.860，均达到极显著水平(p<0.01)。

表1 不同时段年降水量、径流量和输沙量平均值

注：第Ⅰ，Ⅱ和Ⅲ时段为18 a的完整周期，第Ⅳ时段是6 a的不完整周期。

表2 不同时段年降水量与径流量、输沙量相关关系

注：x1为降水量；y1为径流量；y2为输沙量；*为p<0.05；**为p<0.01。R为相关系数；R2为贡献系数；F为方差检验值，Sig.为显著性检验值。

2.3.2 水土保持措施与径流量、输沙量的关系 根据祖厉河流域水利水土保持措施对入黄水沙变化的影响及发展趋势研究(1957—1989年)[15]、甘肃省水土保持综合治理措施效益研究(1979—2012年)[16]及甘肃省水土保持年报资料(2013—2016年)会宁县历年水土保持措施保存面积资料，采用按比例分摊法，按照研究区面积占会宁县的比例，计算得到历年研究区水土保持措施面积。研究区按照水土流失治理速率，大体上可分为1957—1984年、1985—2000年和2001—2016年3个阶段(图8A)。1957—1984年为分散零星治理阶段，治理速度小、措施比较单一，主要措施以造林和梯田为主体，至1984年底水土保持措施面积达到1.09万hm2，其中梯田面积增加到0.53万hm2，造林面积增加到0.30万hm2，种草面积增加到0.14万hm2；1985—2000年以小流域为单元综合治理阶段，这个时期实施的国家级、省级梯田建设工程、小流域治理项目较多，治理速度加快、措施质量提高，至2000年底水土保持措施面积达到4.01万hm2，其中梯田面积由0.53万hm2增加到1.23万hm2，造林面积由0.30万hm2增加到1.16万hm2，种草面积由0.14万hm2增加到1.25万hm2；2001—2016年，为稳定提高时期，此阶段除实施梯田、小流域综合治理项目外，开展了大规模退耕还林(草)工程，至2016年底水土保持措施面积达到5.93万hm2，其中梯田面积由1.23万hm2增加到2.33万hm2，造林由1.16万hm2增加到1.51万hm2，种草面积由1.25万hm2增加到1.51万hm2，水土流失治理程度达到60.31%。

图8 水土保持措施量(A)和径流量、输沙量变化(B)

由图8B可以看出，随着水土保持措施量的逐渐增加，径流量和输沙量呈现波动减少趋势。由回归分析可知(表3)，在1957—1984年，径流量、输沙量对水土保持措施量变化响应不强烈，都未达到显著水平(p>0.05)，水土保持措施量对径流量的决定系数R2为0.076，对输沙量的决定系数R2为0.027。在1985—2000年，径流量对水土保持措施量变化响应强烈，达到显著水平(p<0.05)，水土保持措施量对径流量的决定系数R2为0.309，对输沙量的决定系数R2为0.115。在2001—2016年，径流量、输沙量对水土保持措施量变化响应增强，分别达到极显著水平(p<0.01)和显著水平(p<0.05)，水土保持措施量对径流量的决定系数R2为0.390，对输沙量的决定系数R2为0.352。平均径流模数和输沙模数分别由1957—1984年1.90万m3/(km2·a)和0.56万t/(km2·a)，减少为2001—2016年0.42万m3/(km2·a)和0.13万t/(km2·a)。说明水土保持措施对径流量和输沙量的影响逐渐增大，随着水土保持治理措施增加，水土流失量迅速下降。

表3 不同时段水土保持措施量和径流量、输沙量相关关系

注：x2为水土保持措施量。

3 讨 论

研究区与祖厉河流域相比，其上游(研究区)降水量、径流量和输沙量变化趋势，与韩通[17]、张富[11]等对祖厉河流域1957—2006年及1955—2013年的研究结果一致，但祖厉河流域降水量呈显著减少趋势，在上游地区这种趋势没有达到显著水平；而径流量、输沙量均显示了相同的显著减少的变化趋势，但祖厉河流域径流量、输沙量的突变点分别出现在1995年和2000年但上游分别出现在2000年和2003年，分别滞后了5 a和3 a。黄维东[18]在对祖厉河流域水沙变化研究中发现，降水对水沙减少量影响是30%，人类活动影响占70%。刘淑燕[19]、李宜坪[20]等对黄土丘陵沟壑区水土保持措施对径流输沙的影响进行研究，结果表明增加水土保持措施可有效的减少流域内径流输沙量，水沙关系受降水量、降水强度和土地利用变化的显著影响，并说明坡改梯是引起水沙变化的主要原因，坡面工程措施增加了地表形态。本研究中，降水对径流的决定系数R2由0.482减为0.114，水土保持措施量对径流量的决定系数R2为由0.076增为0.390；降水量对输沙量的决定系数R2由0.385减为0.012，水土保持措施量对输沙量的决定系数R2由0.027增为0.352，表明随着时间推移，水土保持治理措施的增加和治理程度的提高，对水沙变化的影响由降水量变为以水土保持治理措施为主的人类活动。

4 结 论

(1) 研究区年降水量呈现波动的不明显减少趋势，年际存在明显的“丰—枯”周期性变化。最大年降水量是最小年降水量的2.5倍；年降水量存在3种尺度的周期变化，分别是4～8 a，11～15 a和16～32 a，有5 a，8 a，22 a，27 a4个变化主周期，存在3 a，5 a，13 a，18 a平均变化周期。在18 a平均周期下，预测降水量将由枯转丰。

(2) 径流量、输沙量呈现相似的波动减少趋势，年际存在明显的“丰—枯(多—少)”周期性变化。最大年径流量是最小年径流量的1 516倍，最大年输沙量是最小年输沙量的1 238倍；同时径流量和输沙量分别在2000年和2003年出现突变点。年径流量和输沙存在3种相同尺度的周期变化，为3～7 a，8～12 a和14～32 a，径流量有4 a，9 a，27 a变化主周期，存在3 a，6 a，18 a平均变化周期；而输沙量有4 a，6 a，8 a，27 a变化主周期，存在3 a，4 a，5 a，18 a平均变化周期。在18 a平均周期下，预测径流在2021年左右将由枯变丰，输沙在2021年将由少变多。在平均周期为18 a的3个周期下，径流、输沙量对降水量变化响应的趋势协同性由强逐渐减弱，降水量对径流量的相关系数R由0.69减小为0.34，水土保持措施量对径流量的相关系数R为由0.28增为0.62；降水量对输沙量的相关系数R由0.62减小为0.11，水土保持措施量对输沙量的相关系数R由0.16增为0.59，说明对水沙变化的主要影响因素已由降水转变为人为活动。本文对降水量、径流量和输沙量都是以年为时间尺度进行研究，对于年内的变化和水土保持措施以外的人类活动对水沙变化的影响还需进一步研究。