沂沭泗流域年径流非一致性变化规律分析

Non-Stationary Variation Analysis of Annual Runoff in Yishusi Watershed

摘要: 在人类活动和气候变化的共同扰动下，沂沭泗流域的自然水文循环受到影响。因而，对流域水文序列进行非一致性分析有助于科学有效的评估流域水资源量的变化情况。本研究以沂沭泗流域王庄闸水文站上游作为研究区域，收集水文站点及泄洪闸1956~2015年的年径流总量资料及50年代以后长序列年最大洪水资料，采用Mann-Kendall检验方法进行趋势分析，利用Pettitt检验方法进行变点检验。结果表明：1) 沂沭泗流域年径流总量呈显著减少趋势；2) 流域各站点年最大洪水序列均呈下降趋势，且大部分站点通过了显著性0.05的检验。3) 流域年径流总量序列突变时间为1975年。4) 流域年最大洪水的突变时间为1975年前后。

Abstract: Under the combined disturbance of human activity and climate change, the natural hydrological cycle of Yishu River Basin has been destroyed. Therefore, the non-stationary analysis of the watershed hydrological sequence is helpful to better understand the change of water resources in the watershed. In this study, the upstream of Wangzhuang Sluice hydrological station in the Yishusi River basin was taken as the research area, and the long-term annual maximum flood data after the 1950s and the annual runoff data of each hydrological station from 1956 to 2015 were collected. The trend analysis was conducted by Mann-Kendall test method, and the change point test was conducted by Pettitt test method. The results showed that: 1) The total annual runoff in Yishu River Basin showed a significant decreasing trend; 2) The annual maximum flood series of each station in the basin showed a downward trend, and most of the stations passed the significance of 0.05; 3) The abrupt change of total annual runoff in the study area was 1975; 4) The mutation time of annual maximum flood in the study area was around 1975.

1. 引言

水资源合理开发利用是实现区域生态–经济健康持续发展的关键 [1]。流域发生洪水、干旱等极端水情会为居民生命财产安全带来严重危害，然而在城镇化、工业化等人类活动的作用下，流域自然水文循环发生了巨大变化，进而导致水文序列非一致性等水文变异情况，使得流域内极端水情发生概率提升 [2]。水文变异是指水文序列在受到外界条件干扰后发生了质的改变，变异点前后水文序列的特性表现出了明显的差异 [3]。进行各种水利工作的重要前提是水文序列具有一致性等特性。然而，水文变异致使水文序列的一致性遭到破坏，从而为水利工程建设、管理等工作带来诸多问题 [4]。因此，对流域的水文特征进行一致性分析对于科学有效的评估流域水资源量的变化情况具有重要意义。

水文序列变异诊断方法主要有Pettitt检验 [5]，Lepage检验 [6]，滑动T检验法 [7]，R/S法 [8]，Mann-Kendall (M-K)法 [9] 等。秦琳琳 [10] 等人选用M-K检验等方法分析黄泽江流域1958~2017年径流和降水的变化趋势和变异情况，结果表明，研究时段内径流呈减少趋势并且在1963、2008和2009发生突变，降雨呈不显著增加趋势；王栋 [11] 等人采用Pettitt检验等方法分析了黄河源区1961~2016年气温、径流等的趋势、突变情况，结果显示，在研究时段内气温升高，并在1997年发生突变，年径流量减少，并且突变点为1989年。

沂沭泗流域是我国重要的棉粮油生产基地和煤电能源基地，人类活动频繁，使得当地的水生态情况紧张，流域内频繁遭遇水资源短缺、洪涝等问题。张仙娥 [12] 等人研究表明沂沭泗流域1956~2010年的年降水和年地表水资源均呈下降趋势；薛丽芳 [13] 等人研究指出沂河流域1951~2002年的洪峰流量等水文资料在人类活动的影响呈现下降趋势。

以往学者对沂沭泗流域径流的研究多为趋势分析，且研究时段较早。流域的水资源情况对当地的经济发展十分重要，对流域径流资料进行变异分析，有助于更好的进行流域水资源规划管理工作。因此，本研究将以沂沭泗流域王庄闸水文站上游区域为研究区域，基于各水文站点1956~2015年年径流总量资料和长序列年最大洪水资料，采用Mann-Kendall检验法和Pettitt检验法分别对其进行趋势分析和变异分析。

2. 研究区概况与数据来源

沂沭泗流域位于东经114˚45'~120˚20'，北纬33˚30'~36˚20'之间，隶属于淮河流域，流域相对比较独立，位于淮河流域的东北部。流域横跨我国四个省份，分别是山东、江苏、安徽、河南 [14]。沂沭泗流域属于我国缺水地区，人均地表水资源量不足国内平均水平的1/6，多年平均水资源总量为210.8亿m³，其中地表水资源量为142.6 亿m³ [15]。

本文所选研究区域为沂沭泗流域王庄闸水文站上游(图1)，集水面积约1.4万km²。研究中所用数据资料为大官庄站、临沂站、新沭河泄洪闸、人民胜利堰闸1956~2015年月径流资料和20世纪50年代以来长序列年最大洪水序列资料。大官庄站是由人民胜利堰闸和新沭河闸站组合而成，是沂河东调的水利枢纽，人民胜利堰闸站控制从沂河东调西来水流，新沭河闸站控制沭河北来的水流，故本研究中大官庄站年径流量为新沭河泄洪闸和人民胜利堰闸年径流量之和。本研究中水文资料均来自淮河流域水文年鉴，DEM数据来自地理空间数据云网站。

Figure 1. Location of study area

图1. 研究区示意图

3. 研究方法

3.1. Mann-Kendall检验法

M-K检验法具有不受少数异常值的干扰等优点，常用于检测水文时间序列的趋势 [16]。假定时间序列变量为 $x_1,x_2,\cdots ,x_n$，其样本容量为n，统计量S计算公式如下：

$S={\displaystyle \underset{i=1}{\overset{n-1}{\sum }}{\displaystyle \underset{j=i+1}{\overset{n}{\sum }}\mathrm{sgn}\left(x_j-x_i\right)}}$ (1)

其中，

$\mathrm{sgn}\left(x_j-x_i\right)=\{\begin{array}{l}1if\left(x_j-x_i\right)>0\\ 0if\left(x_j-x_i\right)=0\\ -1if\left(x_j-x_i\right)<0\end{array}$ (2)

方差Var(S)的计算公式如下：

$Var\left(S\right)=\frac{n\times \left(n-1\right)\times \left(2n+5\right)}{18}$ (3)

统计量Z计算公式如下：

$Z=\{\begin{array}{l}\frac{S-1}{\sqrt{Var\left(S\right)}}S>0\\ 0S=0\\ \frac{S+1}{\sqrt{Var\left(S\right)}}S<0\end{array}$ (4)

根据Z值正负判别变化趋势，大于零为上升趋势，反之则为下降趋势。当|Z| ≥ 1.96时，即为通过了置信水平为0.05的显著性检验 [17]。

3.2. Pettitt检验法

Pettitt检验法是一种非参数变异诊断方法。该方法具有物理意义清晰等优点，常被用于水文序列的突变点检验。操作原理为假设将时间序列数据随机分割为均值相等的两部分，当前后两部分数据的分布函数不相等时，则认为分割点为突变点 [18]。计算公式如下：

$\mathrm{sgn}\left(x_j-x_i\right)=\{\begin{array}{l}1if\left(x_j-x_i\right)>0\\ 0if\left(x_j-x_i\right)=0\\ -1if\left(x_j-x_i\right)<0\end{array}$ (5)

$U_{t,T}={\displaystyle \underset{i=1}{\overset{t}{\sum }}{\displaystyle \underset{j=t+1}{\overset{T}{\sum }}\mathrm{sgn}\left(x_i-x_j\right)}}$ (6)

当序列为连续分布时，U_t,T可由下式递推计算

$U_{t,T}=U_{t-1,T}+V_{t,T}t=2,3\cdots T$ (7)

其中，

$V_{t,T}={\displaystyle \underset{j=1}{\overset{T}{\sum }}\mathrm{sgn}\left(x_t-x_j\right)}$ (8)

$U_{1,T}=V_{1,T}$ (9)

$K_t=\text{Max}\left|U_{t,T}\right|$ (10)

其对应的显著性概率计算公式为：

$P=2{\text{e}}^{\frac{-6K_t^2}{T^2+T^3}}$ (11)

当|U_t,T|取得峰值时，该点即为可能的突变点，当对应的P值小于0.05时该突变点为有效突变点。

4. 结果与讨论

4.1. 趋势性分析

采用M-K趋势检验法对沂沭泗流域各水文站点和整个研究区域的年径流量序列进行趋势分析，结果见表1。

Table 1. Trend analysis of the annual runoff

表1. 年径流量的趋势分析

Figure 2. The variation trend of annual runoff in the study area

图2. 研究区域年径流量变化趋势

由表1可知，除人民胜利堰闸站外，各站点呈现出不同的减少趋势。其中新沭河闸站径流呈不显著减少趋势；人民胜利堰闸站径流呈增加趋势，但未通过显著检验；临沂站径流呈减少趋势，并且通过了所设定的显著性检验，减少趋势显著；大官庄站呈减少趋势但未通过显著性检验，减少趋势不明显；1956~2015年间整个研究区域年径流量序列呈减少趋势，且减少趋势显著(图2)，与潘扎荣 [19] 等人得出的结论一致。流域内降水量减少等气候条件变化和城镇化等人类活动均会使得径流减少。沂沭泗流域过去几十年降水呈减少趋势，并且占研究区域大部分面积的临沂市自1978年以来城镇化速度较快，城镇建设用地面积明显增加 [20]。此外，随着经济增长，流域内各产业用水量也随之增加，这些原因均导致了研究区域年径流的减少。

Table 2. Trend analysis of the annual maximum peak discharge

表2. 年最大洪水序列的趋势分析

表2为年最大洪水序列的趋势检验结果，除大官庄站呈小幅下降趋势外，新沭河闸、人民胜利堰闸及临沂站均呈下降趋势，并且都通过了显著性检验，下降趋势显著；因此，研究区域年最大洪水呈下降趋势，且下降显著；研究区域年最大洪水与年径流总量的变化趋势相同。研究表明 [13] 研究区域多年典型暴雨–洪峰流量关系具有一致性，张爱军 [21] 等人分析沂沭泗流域1960~2011降水特性，得出结论流域极端降水特征值呈下降趋势，这与本研究中年最大洪水的变化趋势一致。由于流域气候特点、黄河改道、农业开发等原因致使流域内水患频发。为了缓解流域内水灾害等问题，自1957年以来，沂沭泗流域进行了修建水库、泄洪道、水坝等水利建筑以及增加灌溉水源等规划建设 [22]，并且成效明显。同时，由于流域内降水量减少等原因共同导致了年最大洪水的下降。

4.2. 变异分析

采用Pettitt检验对沂沭泗流域各站点及全区域年径流总量序列进行一致性检验，结果见表3。

Table 3. Abrupt change of the annual runoff

表3. 年径流总量突变检验结果表

由表3可以看出，新沭河闸站和临沂站突变时间都为1975年，人民胜利堰闸站突变时间为1969年；大官庄站未见明显突变；全区域年径流总量序列突变时间为1975年。潘扎荣 [19] 等人研究指出沂沭泗流域内临沂站径流突变时间为1976年。结合研究区域实际情况，该流域于1958年启动蓄水工程的修建，1964年前后完成了流域内第一批水库的蓄水工作，在1965~1975年之间流域内蓄水工程进一步扩大，因此，本研究中年径流总量突变时间为1975年，研究结果符合流域实际情况。

Table 4. Abrupt change analysis of the annual maximum peak discharge

表4. 年最大洪水突变检验结果表

年最大洪水序列Pettitt检验结果见表4，在显著性水平0.05下，新沭河闸站和临沂站突变显著，且突变年份分别为1976年和1974年，人民胜利堰闸站和大官庄站没有显著的突变点。薛丽芳 [13] 等人分析发现临沂站年最大洪水与同期的雨量特征变化基本一致。而研究区域20世纪50~70年代中期降水比较丰沛，70年代中期到80年代后期降水偏少，80年代后期以后为平水期。因此，研究区域年最大洪水突变时间为1975年前后。

4.3. 讨论

本研究采用M-K检验法分析研究区域水文特征序列的变化趋势，结果显示研究区域年径流和年最大洪水序列均呈下降趋势。M-K检验法是由世界气象组织推荐使用的一种趋势分析方法，该方法不受少数异常值干扰的特点，使其非常适用于水文序列的趋势分析，因而被众多学者用于水文、气象时间序列的趋势研究，并取得一定成果 [23]。但该方法仍然有未考虑序列自相关性等局限性。

同时，采用Pettitt检验方法分析发现，研究区域年径流序列在1975年发生突变，年最大洪水序列在1975年前后发生突变。Pettitt检验在水文序列变异诊断方面被广泛应用，其优点在于物理意义清晰、突变点识别准确等，帮助众多学者在不同流域的水文变异分析中取得进展 [24]。但该方法仍然有只能单点检测、易受序列长度影响等局限性。采用了M-K检验和Pettitt检验分析沂沭泗流域水文特征序列，以后的研究中要采用多种方法展开研究区域水文特征序列分析。

5. 结论

为了研究沂沭泗流域水文特征的非一致性，本文收集了沂沭泗流域王庄闸水文站上游区域4个站点1956~2015年的年径流总量和年最大洪水资料，采用Mann-Kendall检验法和Pettitt检验法对特征序列进行了趋势检验和突变检验，具体结论如下：

1) 研究区域年径流总量序列呈显著减少趋势，其中临沂站呈显著减少趋势，新沭河闸站和大官庄站呈不显著减少趋势，人民胜利堰闸站呈不显著增加趋势。

2) 各站点年最大洪水序列均呈下降趋势，并且其中临沂站、新沭河闸站、人民胜利堰闸站三个站点通过所设定的显著性检验，表现出明显的下降趋势。

3) 研究区域年径流总量突变时间为1975年，其中人民胜利堰闸站突变时间为1969年，临沂站突变时间为1975年，突变原因主要为研究区域水利工程建设等人类活动，其余站点无显著突变。

4) 研究区域年最大洪水突变时间为1975年左右，其中临沂站突变时间为1974年，新沭河闸站突变时间为1976年，突变原因为气候变化引起的降雨特征改变，其余站点无显著突变。

综上所述，沂沭泗流域的径流量近几十年以来呈现出减少趋势，故水资源量逐渐减少。此外，由于流域下垫面条件发生变化，从而改变了流域的气候条件，由此影响了流域自然的水文循环，致使流域径流的一致性也遭到了破坏。本研究揭示了研究区域内径流变化特征，以期为沂沭泗流域水资源的科学开发利用提供依据，进而为应对流域内可能出现的旱涝灾害提供科学支撑。

基金项目

国家自然科学基金青年项目(52009102)。

参考文献