不同时间尺度岔巴沟流域水沙关系变化

康靖羚, 丁文峰,2, 韩昊宇, 郭宜薇, 王一然

(1.长江水利委员会长江科学院, 武汉 430010; 2.水利部 山洪地质灾害防治工程技术研究中心, 武汉 430010)

黄河是我国的第二大河，华夏民族在此繁衍生息，孕育出了古老灿烂的中华文明。长期以来，强烈的人为活动以及全球气候变化等因素导致出现水少沙多、水沙异源等问题，水沙关系极不协调[1]，黄土高原已成为我国土壤侵蚀最严重的地区[2]。鉴于水土流失的严重性，自20世纪50年代以来，我国政府将黄河流域治理作为全国水土保持工作的重点，开展了一系列生态治理工程，合理开发利用土地资源，有效控制水土流失。20世纪末期，经过治理的黄土高原小流域示范区水土流失治理度由46.1%增加至80.2%，综合治理成效显著[3]。2000年，全国陆续启动了“退耕还林”等生态工程，植被覆盖率迅速增加[4]，截至2015年，黄河中游的水土流失治理度已达60%[5]。潼关控制站实测年均输沙量也从20世纪前半叶的15.92亿t减少到2000—2018年的2.4亿t，其中2015年实测沙量仅有0.55亿t[1,6]。

河流输沙量作为水土流失的集中反映，已成为全球水沙关系研究的重要内容，也是间接研究流域土壤侵蚀产沙的主要手段。由于大理河地处黄河中游粗泥沙集中来源区[7]，作为大理河的一级支流，岔巴沟流域成为研究黄土丘陵沟壑区土壤侵蚀以及水沙关系的理想区域[8]，此地区长期的观测数据也为研究提供了良好基础，了解岔巴沟流域水沙关系及其变化趋势和影响因素对治理黄土丘陵沟壑区乃至黄土高原的水土流失具有重要意义。对岔巴沟流域的相关研究从20世纪80年代就已经开始，目前为止许多学者在此流域对输沙量变化展开研究，并在流域面积与水沙关系的空间变异性[9-11]、不同空间尺度下的泥沙量或泥沙浓度[12-14]、各类流域产沙模型[15-19]以及泥沙动力学等[20-21]方面取得许多成果。但大多数研究侧重于不同空间尺度下的水沙关系及模型建立，对不同时间尺度下的流域水沙关系研究较少。本文通过分析岔巴沟流域1970—2016年的径流泥沙变化规律，探讨年时间尺度以及次洪水时间尺度下流域水沙关系变化趋势和影响因素，以期为流域水沙模型建立提供依据。

1 研究区概况

岔巴沟流域位于陕西省子洲县，属于黄土高原丘陵沟壑区第一副区，地处东经109°47′，北纬37°31′。流域面积205 km2，沟壑密度1.05 km/km2，主沟比降7.57‰，海拔高度900～1 100 m[22]。属半干旱区，为温带大陆性季风气候，干燥少雨，降雨年内分配不均，大部分集中于6—9月，且多是高强度短历时暴雨[23]。

岔巴沟流域土壤由黄土发育而成，结构疏松，易受侵蚀[21]。早在1959年，岔巴沟流域的水土流失治理工作就已开始，70年代后，流域内大量水库、淤地坝建成，并逐渐发挥效益，至1977年末，所修建库坝数量达444座，总库容2 651万m3。至1980年，岔巴沟流域治理面积已达30.7 km2，占总面积的15%[23]。截至1993年，水土保持治理程度已达到70%以上[19]。90年代后期，岔巴沟流域响应国家号召，开展了大规模的水土保持工作，生态环境明显好转[24]。

2 数据来源及处理

2.1 数据来源

本文所用数据资料出自《黄河流域水文年鉴》，包括雨量站降水数据以及水文站径流和泥沙数据(1970—2016年)。本次研究选择4个测站，分别为和民墕站、刘家坬站、姬家硷站和曹坪站，具体信息见表1。

表1 测站信息

2.2 数据处理

2.2.1 次洪水数据处理 一次完整的洪水事件是指从径流量增加至基流以上，到仅基流构成径流量为止[25]，故采用数字滤波法[26-27]去除基流，得到新的洪峰流量和径流量。研究变量包括次洪水输沙量(ESY)，采用式(1)计算；洪峰流量(Qw，m3/s)，径流过程的最大值；径流深(H，mm)，即径流量与流域面积之比；次降雨量(P，mm)，即一次降雨总量；降雨历时(D，h)，次降雨持续时长；平均降雨强度(I，mm/h)，即降雨量与降雨历时之比以及径流系数(Rc)，即径流深和降雨量之比。

ESY=∑SqΔt

(1)

式中:ESY为次洪水输沙量(万t);Δt为水文站水文观测的时间间隔(s);q为Δt时间间隔内的流量(m3/s);S为Δt时间间隔内的径流含沙量(t/m3)。

从所有筛选出的数据中，选取洪峰流量大于1 m3/s、有对应降雨记录的洪水过程进行分析，共计106场洪水事件。

2.2.2 年尺度洪水数据处理 为研究年时间尺度下降雨、径流和泥沙变量之间的关系，采用年平均输沙量(ESYa，万t)、年最大洪峰流量(Qwa，m3/s)，即一年中全部洪水事件中的最大洪峰流量值、年平均径流深(Ha，mm)、年平均降雨量(Pa，mm)、年平均径流系数(Rca)以及年最大24 h雨强(I24 h，mm/h)等研究变量进行分析。采用Mann-Kendall检验法[28-30]对年尺度研究变量进行趋势分析和突变点检测。数据处理均采用Excel 2016以及SPSS 25，图像绘制采用ArcGIS 10.2以及Origin 2019。

3 结果与分析

3.1 年尺度研究变量分析

表2显示了岔巴沟流域降雨量、径流量和输沙量的年际变化，可以看出降雨量呈增加趋势，径流量与输沙量呈减小趋势，降雨量、径流量和输沙量的最大值分别大于平均值13%，23%和117%，最小值分别比平均值小13%，20%和80%，年际变化较大。

表2 岔巴沟流域年降雨量、径流量和输沙量统计

3.1.1 年降雨量趋势分析 根据实测降雨资料，在ArcGIS 10.2中采用泰森多边形法计算流域面降雨量，并绘制年降雨量变化过程线，见图1。M-K检验曲线见图2。

图1 年降雨量变化过程

图2 年降雨量M-K检验曲线

统计量Z=2.75，|Z|>2.32，表明年降雨量在99%置信水平上呈显著增加趋势。显著水平临界线内出现的交点分别位于1972—1973年与1986—1987年，但1972—1973年无前期研究数据，无法确定在此时间段内发生突变，所以将第2个交点作为突变点，即降雨量于1986年出现突变。

3.1.2 年径流量变化分析 图3即流域年径流量变化过程线，M-K检验曲线见图4，统计量Z=-2.18，|Z|>1.96，表明年径流量在95%置信水平上呈显著减少趋势。显著水平临界线内出现的两个交点分别位于1991—1992年和2001—2002年，即年径流量突变点为1991年、2001年。

图3 年径流量变化过程

图4 年径流量M-K检验曲线

3.1.3 年输沙量变化分析 流域年输沙量变化过程线见图5，M-K检验曲线见图6，统计量Z=-3.13，|Z|>2.32，表明年输沙量在99%置信水平上呈显著减少趋势。显著水平临界线内出现一个交点，位于1982—1983年，故1982年是年输沙量的突变点。

图5 年输沙量变化过程

图6 年输沙量M-K检验曲线

年降雨量的突变不受人为因素控制，但年输沙量在1982年发生突变的原因是70年代大量修建的淤地坝在这个时期产生效用，从而大幅减少水土流失。通过调查岔巴沟流域水土保持措施得知[31]，1991年后，岔巴沟流域再无新建淤地坝，同时岔巴沟流域Landsat MSS/TM遥感资料[32]显示，近年来植被恢复措施也在大力推进，植被覆盖率逐渐由中低覆盖度向中高覆盖度转化，至2007年，中覆盖度、中高覆盖度及高覆盖度分别为51.33%，26.95%和20.65%。这就是年径流量产生突变在1991年后短暂上升后又下降的原因。通过M-K检验得知，年降雨量呈增加趋势，而径流量和输沙量却始终呈现减少趋势，这也印证了淤地坝建设以及植被恢复等水土保持措施对控制流域水土流失所起的积极作用。

3.2 次洪水研究变量趋势分析

采用年输沙量突变点1982年将研究时段划分为两个时期，第一时期是1970—1981年，共60场洪水，考虑到剩余时期时间序列过长，将其划分为第二时期和第三时期，第二时期是1982—2000年，共28场洪水，第三时期是2001—2016年，共18场洪水。3个时期研究变量的平均值见表3。

表3 研究变量平均值

由表3可知，输沙量和径流系数呈减小趋势，且减小幅度较大，与1970—1981年相比，2001—2016年输沙量和径流系数分别减少了约67%和83%；洪峰流量、径流深、降雨量、降雨历时在1982—2000年降低，后在2001—2016年又增加，但洪峰流量仍低于1970—1981年均值，约减少13%，径流深、降雨量、降雨历时均大于1970—1981年均值，分别增加了约40%，69%和55%。

比较不同时间尺度下研究变量的变化趋势，发现降雨量和输沙量在年尺度和次洪水尺度下表现一致，但径流深却在次洪水尺度下呈增加趋势，究其原因，虽近年来径流总量减少，但洪水发生次数也在减少，2001—2016年径流总量约占70年代径流总量的3/4，而发生的洪水次数仅为70年代的1/3，其减小幅度远大于径流总量的减小幅度。

3.3 不同时间尺度水沙关系回归分析

3.3.1 次洪水时间尺度回归分析 构建岔巴沟流域各研究时期次洪水降雨、径流和泥沙变量的皮尔逊相关系数矩阵，1970—1981年的计算结果详见表4。

表4 1970-1981年次洪水研究变量的相关系数矩阵

1970—1981年ESY与Qw，H和P具有极显著的相关性，与D呈显著相关性，表明Qw，H，P以及D在此研究时期是ESY的主要影响因素。1982—2000年ESY与Qw，H和P具有极显著的相关性，2001—2016年ESY仅与Qw和H具有极显著的相关性，与P呈显著相关性，选择上述与输沙量具有显著相关性的研究变量进行回归分析，回归结果见表5。

表5 研究变量与次洪水输沙量的回归方程

ESY在3个时期均与Qw，H以及P具有显著的相关性，但整体来看，ESY与Qw和H的相关关系最为显著和稳定，虽然决定系数有所浮动，但始终在p<0.01水平上呈显著相关性。从而建立ESY与Qw和H的回归方程，回归模型的决定系数均大于0.9，说明ESY始终与Qw和H存在良好的幂函数关系。H和Qw可以很好的反映输沙能力，这在许多学者的研究成果中已经体现，但与Zhang[11]、Zheng等[15]分析单个因素与输沙量呈线性关系或比例关系不同，本文研究得出洪峰流量和径流深与输沙量呈幂函数关系，这与肖学年等人[10]的分析一致。考虑到本文采用H作为影响因素分析，与P，D同为次降雨指标，它们之间具有相关性，故未将P和D作为回归因子进行分析。

回归方程中，回归系数代表研究变量对输沙量的影响程度，回归系数在1982—2000年最大，在2001—2016年最小，表明上述研究变量在3个时期对输沙量的影响整体呈先增大后减小的趋势；变量所对应的指数越大，表明它对输沙量的影响越大，前两个时期径流深对输沙量的影响大于洪峰流量，然而2001—2016年洪峰流量的指数明显增大，径流深则无显著变化，表明洪峰流量对输沙量的影响更大。同时，回归模型的决定系数逐时期下降，说明岔巴沟流域的水沙关系逐渐变的复杂，而在上述影响输沙量的相关变量的研究中，并未将水土保持因子考虑在内，尤其是随着时间推移，淤地坝与植被恢复工程发挥的水土保持效益应呈现相反态势，这也使得2001—2016年回归模型的决定系数低于前两个时期。

根据Li等[33]研究，暴雨是造成黄土高原水土流失的主要原因，一场严重暴雨造成的土壤流失可占年土壤流失的60%～90%。采用极端降雨阈值[34]评估岔巴沟流域3个时期的降雨量，具体而言，就是将研究时段内降雨量大于0.1 mm的降雨事件升序排列，将第90个百分位数值认定为该站点研究时段内极端降雨的阈值，当次降水量超过此阈值时，便认为此次降雨为极端降雨事件。极端降雨事件在3个研究时期的发生频率分别为10%，10.7%和11.1%，相关特征值见表6。

表6 极端降雨事件特征值

对比不同时期相似的极端降雨条件下的径流量和输沙量，如降雨事件2009-07-19与1987-08-26，以及2006-08-29与1978-08-07，发现2001—2016年极端降雨事件所产生的径流量和输沙量明显减少，表明水土保持措施对减少径流泥沙起着重要作用，但输沙量占总输沙量的比例却明显增多，分别为7%，18%和43%。在发生频率大致相同的情况下，输沙量占总输沙量的比例远大于前两个时期，这表明极端降雨事件在2001—2016年内对流域水沙关系的影响明显变大，导致水沙关系发生了显著变化，极大的增加了水沙关系预估的不确定性，这与Zhao等[34]研究所得结论一致。

3.3.2 年时间尺度回归分析 构建岔巴沟流域各研究时期年降雨、径流和泥沙变量的皮尔逊相关系数矩阵，1970—1981年的计算结果详见表7。

表7 1970-1981年研究变量的相关系数矩阵

1970—1981年输沙量与H和Qwa具有极显著的相关性，与Rca呈显著相关性，表明H，Qwa和Rca在此研究时期是输沙的主要影响因素。1982—2000年输沙量与H呈极显著的相关性，与Rca呈显著相关性，2001—2016年输沙量与Qwa呈极显著的相关性，与I24 h呈显著相关性。对3个时期的年尺度水沙关系分别建立回归函数模型，结果见表8。

表8 研究变量与年平均输沙量的回归方程

考虑到Pa，Rca与H同为次降雨指标，故未将Pa和Rca作为回归因子进行分析。回归方程中决定系数均大于0.96，说明ESYa始终与Qwa和Ha存在良好的多元幂函数关系。与目前Zheng等[15]在年时间尺度上通过单因素分析得出的线性关系相比，研究变量为径流深和最大洪峰流量时回归方程的决定系数更大，回归效果更好，同时回归模型也与次洪水尺度下相同，均为幂函数。然而回归系数整体呈降低趋势，表明研究变量对输沙量的影响越来越小；在前两个时期中径流深的指数大于最大洪峰流量的指数，表明径流深对输沙量的影响更大，而2001—2016年回归关系发生改变，最大洪峰流量的指数远大于径流深的指数，此时最大洪峰流量对输沙量的影响更大。

上述回归分析表明，不同时间尺度下水沙关系变化趋势不同。在年尺度下回归方程的决定系数呈增加态势，而次洪水尺度下呈降低趋势。这种差异主要是由于研究数据和变量的时间尺度不同而导致的，年尺度回归分析中采用的研究事件除已选取的次洪水事件之外，还包含其他因洪峰流量小、降雨历时短以及缺乏降雨资料等原因未定义为一次洪水事件的径流过程。统计所选取的次洪水事件，发现3个研究时期中次洪水事件的总输沙量占年平均输沙量的比例分别为76%，71%和59%。由此可见，年尺度回归分析同次洪水回归分析相比，关于流域输沙过程研究的全面性和普遍性逐时期递增，所以年尺度回归分析的决定系数逐渐增加并大于次洪水回归分析的决定系数。然而年尺度回归分析无法细致的表征洪水事件的相关特征变量，这从回归系数数值较低且逐时期减小的趋势也可以看出，表明年尺度分析对流域水沙关系的研究无法深入，故目前的研究多集中于次洪水尺度下的水沙关系分析。

4 结 论

(1) 岔巴沟流域1970—2016年降雨量呈显著增加趋势，年径流量和年输沙量呈显著下降趋势，年输沙量、年降雨量和年径流量分别于1982年、1986年、1991年、2001年发生突变，106场次洪水事件分析结果表明降雨量和输沙量的变化趋势与年尺度上表现一致，但径流深呈增加趋势，主要是洪水次数大幅减少所致。

(2) 在次洪水时间尺度上，输沙量与洪峰流量和径流深始终呈现显著相关关系，其中洪峰流量对输沙量的影响逐渐增大，径流深则无显著变化。

(3) 在年时间尺度上，输沙量整体与最大洪峰流量以及径流深关系显著，前期径流深对输沙量的影响较大，后期最大洪峰流量对输沙量的影响较大。