黄河上中游水土保持减沙效果研究

2016-10-09李敏

中国水土保持 2016年9期

关键词：输沙量回归方程降水量

李　敏

(黄河上中游管理局，陕西西安 710021)

黄河上中游水土保持减沙效果研究

李敏

(黄河上中游管理局，陕西西安 710021)

年输沙量；水土保持；减沙效果；黄河上中游

采用逐步回归方法，对1954—2012年黄河年输沙量、降水量和水土流失治理度数据分3个时段进行了计算，分析认为：水土流失治理度较好地“解释”了20世纪80年代以来黄河年输沙量锐减的现象，回归方程的显著性水平均达到α=0.001；在1967—1987年时段，相同降雨条件下，水土流失治理度的提高不仅显著减少了黄河年输沙量，而且使年输沙量对降雨的响应趋于“迟缓”，在1988—2012年时段，仅由水土流失治理度一个因子就能够很好地“解释”黄河年输沙量的变化；在水土流失治理度超过55%时，计算年输沙量仅1.2亿t，以80%的置信度，黄河年输沙量将小于5亿t。

1　概述

1.1黄河上中游水土保持发展进程

根据各省区年度统计上报，从20世纪80年代初期开始，黄河上中游每年治理水土流失面积1.2万km2左右，到2015年累计治理水土流失面积超过37万km2，连同20世纪80年代以前的治理面积，合计统计上报水土流失治理面积超过45万km2。受人为和自然因素影响，统计上报的水土流失治理面积数据中存在一定的“水分”。1985年以来，通过对典型小流域和重点治理区进行实地测量调查，不断核实实际保存的水土流失治理面积。经过大量的核实工作，黄河中游一般梯田保存率为90%～95%、水保林保存率为70%～80%、种草保存率约为40%，故采用60%左右的综合保存率对统计上报数据进行核实。根据黄河中游水土保持委员会第十二次会议报告推算，截至2015年底，黄河上中游共有水土流失综合治理面积约26万km2(表1)。从表1可以看出，黄河上中游地区自1980年以来水土流失综合治理面积基本呈线性增加趋势。

1.2黄河年输沙量变化

表1　黄河上中游不同时期的水土流失综合治理面积

1954—2012年的59年间黄河年输沙量波动变化大体分为3个阶段：①1954—1966年输沙量持续高位，年输沙量超过25亿t的有5个年份，时段平均年输沙量为15.66亿t；年际间变化剧烈，1964年输沙量超过26亿t，1965年下降到仅4.26亿t，1966年又超过26亿t[1-2]。②1967—1987年为输沙量骤降阶段，年输沙量从1967年的26亿t以上波动减少到1987年的4亿多t，期间出现1977年的特殊高输沙量年份。③1988—2012年为年输沙量持续递减阶段，时段平均年输沙量5.7亿t，进入21世纪后平均值仅为2.6亿t(见图1)。

图1　1954—2012年黄河年输沙量变化

2　研究内容

由上述分析知，黄河年输沙量变化过程和黄河中游水土保持发展进程的变化趋势具有相似性，而黄河泥沙又基本上来自流经黄土高原地区的渭河、北洛河、汾河和河口镇(头道拐)到龙门区间(以下简称“河龙区间”)的各支流，其中主要产沙区域为河龙区间。因此，本研究选择渭河(华县站)、北洛河(氵状头站)、汾河(河津站)控制区域和河龙区间(以下简称“4站区域”)为研究区，探讨黄河中游水土流失综合治理措施的减沙效果。4站区域总面积30.05万km2。黄河自1919年设站观测，距今不到100年，选取资料较完整的1954—2012年作为研究时段，共计59年。

3　研究方法

根据水文学和水土保持学基本原理，河流年输沙量的变化遵循统计学规律，水土流失综合治理措施的减沙效益亦遵循统计学规律。因此，本研究将研究区域黄河年输沙量的产生过程视为一个灰箱或者黑箱，采用逐步回归分析方法，筛选具有物理成因概念和含义的解释因子(自变量)，对不同时段的资料进行计算处理，获得符合水文学和水土保持学基本原理、仿真模拟效果好的计算结果，分析水土保持对黄河年输沙量的影响。理论和研究成果表明，输沙量与降雨和水土流失治理度等因子之间表现为非线性相关关系。为此，本研究对年输沙量取对数形式，考虑到降雨对输沙量影响的复杂性，对部分自变量也增加对数形式，以获得更好的结果。

3.1降雨因子选择

由于已知黄河泥沙主要源自河龙区间和汛期降雨(暴雨)产沙，因此对于降雨因子，初步选取4站区域年降水量及河龙区间年降水量、7—8月降水量、最大3日降水量、最大1日降水量等。河龙区间是黄河中游三大暴雨区之一，是黄河泥沙的主要来源区。为了精细地仿真模拟黄河年输沙量的变化，将河龙区间划分为5个亚区(图2)[3]，其中又将河龙区间的陕西北片和晋西北片合并为“河吴区间”片，采用这些亚区7—8月降水量、最大3日降水量和最大1日降水量等作为“解释”黄河年输沙量的备选自变量。以上共初步选择22个降雨因子。考虑到59年系列降雨与输沙量关系的复杂性，增加了降雨因子的自然对数(ln)，使参与计算筛选的降雨变量达到44个。

图2　河龙区间降雨分区示意

3.2水土保持因子选择与预处理

为了使坡面治理措施与淤地坝相互协调，既不重复计算，又不产生矛盾，计算中对坡面措施面积、大型坝(骨干坝)坝控面积进行“耦合”与综合处理。首先，将核实的坡面治理措施面积加上大型淤地坝(骨干坝)的坝控面积，并扣除二者的重叠面积，形成“水土流失综合治理面积”(图3)，其中大型淤地坝(骨干坝)按20年设计淤积年限动态计算面积[4-5]。

图3　水土流失综合治理面积示意

为了使今后对不同时间和空间尺度的计算结果具有可比性，本研究将水土流失治理度作为参选自变量。水土流失治理度的计算公式为

式中：D为水土流失治理度，%；A1为水土流失治理面积，万km2；A2为水土流失面积，万km2。

4　计算结果与分析

4.11954—1966年计算结果

1954—1966年研究区域黄河年输沙量较大、年际间波动剧烈，水土保持工程规模较小，而且毁林开荒、毁林炼钢，人为水土流失剧烈。对这一时段数据的计算结果表明，仅河吴区间7—8月降水量一个因子就可以“解释”该时段年输沙量近80%(调整的R2=0.776)的变化，河吴区间7—8月降水量与年输沙量回归计算结果见表2。

表2　1954—1966年计算结果

比较回归计算值与实测值(图4)可以看出，回归计算年输沙量与实测年输沙量拟合较好，散点图显示计算年输沙量为实测年输沙量的95.58%。其中，1966年标准化残差接近2、计算值明显小于实测值。初步分析，造成这一现象的原因是该年黄河流域发生较大的局部暴雨，冲毁了部分水土保持工程措施，造成以往多年拦淤减少的泥沙在当年被冲入黄河，使“河吴区间7—8月降水量”不能够“解释”这些输沙量，因而计算年输沙量明显小于实测年输沙量[6-8]。

图4　1954—1966年实测年输沙量与计算年输沙量比较

4.21967—1987年计算结果

4.2.1总体分析

对1967—1987年数据采用逐步回归方法，通过4次计算，先后筛选出由河吴区间7—8月降水量的自然对数、水土流失治理度、河龙区间最大1日降水量和区域年降水量的自然对数等4个自变量组成的回归方程，其相关系数R=0.959 0(表3)，显著性水平α=0.001。

在逐步回归计算中，河吴区间7—8月降水量的自然对数第一个进入回归方程(表3)，相关系数R达到0.832 0；第二次计算增加了水土流失治理度，使回归方程的相关系数R超过0.9，调整的R2达到0.789 9，说明仅由河吴区间7—8月降水量的自然对数和水土流失治理度两个因子就“解释”了该时段黄河年输沙量近79%的变化；第三次计算又增加了河龙区间最大1日降水量，使回归方程的相关系数R达到0.937 6，调整的R2达到0.857 7；第四次计算增加了区域年降水量的自然对数，回归方程的相关系数R达到0.959 0，调整的R2达到0.899 5。

表3　1967—1987年计算结果

筛选出的水土流失治理度自变量为负值，符合“水土保持是黄河泥沙治本措施”的基本原理；筛选出的河吴区间7—8月降水量的自然对数、河龙区间最大1日降水量和区域年降水量的自然对数均为正值，符合降雨产沙的基本规律[9-10]。

根据回归方程，得出计算年输沙量，将其与实测年输沙量比较可以看出(图5)，计算年输沙量与实测年输沙量拟合很好，散点图显示计算年输沙量为实测年输沙量的95.76%。

图5　1967—1987年实测年输沙量与计算年输沙量比较

4.2.2标准化回归系数分析

对第四次计算结果的4个自变量系数进行标准化处理(表4)，去掉不同量纲和不同数量级的影响，分析各自变量在回归方程中的相对重要性。水土流失治理度的标准化回归系数为-0.435 5，河吴区间7—8月降水量的自然对数、河龙区间最大1日降水量和区域年降水量的自然对数的标准化回归系数分别为0.315 3、0.324 9和0.250 1。水土流失治理度标准化回归系数的绝对值最大，这一现象说明在由这4个自变量组成的回归方程中，水土流失治理度对年输沙量的影响最大，是该时段黄河年输沙量骤降的主要因素。

表4　1967—1987年回归系数

4.2.3三维空间图形分析

由于水土流失治理度与年输沙量之间是非线性关系，因此不同水土流失治理度对年输沙量的影响程度不同，影响的过程(趋势)也不同。根据回归方程做出在相同降雨条件下不同水土流失治理度不同视角的三维曲面图形，可直观对比分析不同水土流失治理度对年输沙量的影响程度和影响趋势。为了做出三维图形，将贡献率最小的区域年降水量的自然对数固定为时段的平均值，选择时段起始年(1967年)和时段末尾年(1987年)，绘制水土流失治理度对年输沙量影响的曲面进行比较分析(图6)。

由于是对时段内河吴区间7—8月降水量的自然对数和河龙区间最大1日降水量两个降雨因子全部组合的计算，因此在1967年回归曲面中，两个降雨因子同时取最大值时计算的年输沙量超过35亿t。从图6可以看出，1967年和1987年的水土流失治理度对年输沙量影响显著不同。在曲面的最高点，即两个降雨因子均处于最大值时，1987年输沙量响应曲面比1967年输沙量响应曲面降低约20亿t；在曲面的最低点，即两个降雨因子均处于最小值时，两个年输沙量响应曲面相差幅度约为5亿t。对比分析说明，随着水土流失治理度的增加，一方面总体上年输沙量显著减少，另一方面年输沙量响应曲面在高降雨条件下的减少幅度大于低降雨条件。

比较图6中1967年和1987年两个年份的输沙量响应曲面可以看出，1967年治理水平年输沙量响应曲面陡峭，计算年输沙量随降雨的变化而急剧变化，计算年输沙量最高与最低值相差约30亿t；1987年治理水平年输沙量响应曲面不仅较低，而且显著平缓，计算年输沙量最高与最低值相差约15亿t。这一现象说明水土流失治理不仅大幅度减少了河龙区间黄河年输沙量，而且使年输沙量对降雨变化响应的敏感度降低[11]。

图6　年输沙量对不同水土流失治理度的响应曲面

4.31988—2012年计算结果

4.3.1回归计算结果

回归分析表明，1988—2012年时段仅水土流失治理度一个因子就能够“解释”黄河年输沙量77.65%的变化(表5)，回归方程的相关系数R达到0.886 5，显著性水平α=0.001。水土流失治理工程的大规模开展是该时段黄河年输沙量持续处于低位的根本原因。20世纪80年代初期实施的水土保持工程措施持续发挥水土保持效益，同时国家又在重点水土流失地区实施了多个治理项目，水土保持治理工程规模进一步扩大，水土保持对黄河年输沙量的影响更加显著(图7)。

表5　1988—2012年计算结果

4.3.2减沙分析

截取1999—2012年的计算结果，分析计算年输沙量变化趋势(图8)。可以看出：1999年水土流失治理度达到40%，计算年输沙量小于5亿t，其80%置信度上限不到15亿t； 2010年水土保持治理度达到55%，计算年输沙量仅1.2亿t，其80%置信度上限小于4.5亿t。根据1988—2012年数据分析，在当前水土保持治理程度，除非出现稀遇降雨事件，以80%的置信度，黄河年输沙量将小于5亿t。