不同突变检验方法在萧濉新河流域径流序列中的研究

Study on Runoff Sequence of Xiaosuixin River Basin by Different Mutation Test Methods

摘要: 径流序列的突变性研究对于探讨气候变化与人类活动下的水文循环的演变规律具有重要的意义。本文以安徽省萧濉新河流域1973~2022年径流量为例，采用T方法、有序聚类法、滑动t检验、Mann-Kendall检验法、Yamamoto法、Pettitt法等6种方法对比分析水文序列突变特性，结果表明：1973~2022年流域年径流序列的突变点为1987年、1995年、2010年，T方法、有序聚类法、滑动t检验、Mann-Kendall检验法、Pettitt法结论基本一致，其中Mann-Kendall检验法检测突变点最精确，包括了其他检验方法的所有结果，此方法也反映出径流序列的趋势性。

Abstract: The study of the abrupt change of runoff series is of great significance to the study of the evolution of hydrological cycle under climate change and human activities. Taking the runoff of Xiaosuixin River basin in Anhui Province from 1973 to 2022 as an example, six methods, including T method, ordinal clustering, sliding t-test, Mann-Kendall test, Yamamoto method and Pettitt method, were used to analyze the characteristics of hydrological sequence mutation, the results showed that the abrupt change points of annual runoff series from 1973 to 2022 were 1987, 1995 and 2010. The conclusions of T method, ordered Clustering Method, sliding t test, Mann-Kendall test and Pettitt method were basically consistent, the Mann-Kendall test is the most accurate method to detect the mutation point, which includes all the results of other test methods. This method also reflects the trend of runoff series.

1. 引言

随着全球气候变化和人类活动影响加剧，水文循环物理机制和水文循环过程均发生了显著变化，水文序列呈现出非一致性特征 ‎[1] ，尤其表现在水文要素上发生了突变，可当作是跳跃的一种特殊情况，是一种水文序列中少数值偏离原序列的现象。从水文统计角度看，是均值、Cs值、Cv值等发生了改变，使得水文序列不满足一致性的特征。目前，国内外学者对降水序列、径流序列的突变研究较广泛，而研究方法主要有滑动t检验法、Mann-Kendall检验法 ‎[2] ‎[3] 、T方法 ‎[4] 、有序聚类法 ‎[5] 、Yamamoto法 ‎[6] 、TFPW-MK检验法 ‎[7] 、Pettitt法 ‎[8] 、启发式分割算法 ‎[9] 等。

萧濉新河流域在皖北地区已有60多年的煤炭开采历史 ‎[10] ，先后建有大大小小矿区几十个，面积约0.2万km²，其中探明煤炭储量约40亿t，现年产原煤约5800万t。煤矿矿区的大规模开采，下垫面条件发生了显著的变化，从上世纪50年代开始产生了煤陷区，但当时煤炭开采量较小，沉陷区面积较小。自上世纪80年代开始，煤炭开采逐步增大，沉陷区面积增长迅速，近几年已积水成湖，相互连通，占据了城市和流域较大的面积。目前对皖北地区煤陷区的研究集中于湿地建设、水资源利用、生态修复等方面 ‎[11] ‎[12] ‎[13] ，对降雨–径流相关关系、流域径流变化特征研究较少。由于突变物理机制研究尚不十分成熟，在确定水文序列发生突变现象时选择多种方法进行比较较好。为了系统了解气候变化及人类活动下萧濉新河流域径流突变特点，本文采用中位数、有序聚类法、滑动t检验、Mann-Kendall检验法、Yamamoto法、Pettitt法等方法对比分析水文序列突变检验结果的异同，确定萧濉新河流域径流序列的突变点，为流域水文序列演变规律研究和区域社会经济发展提供支撑。

2. 流域概况

萧濉新河位于淮河中游，东经116˚25'~117˚10'，北纬33˚44'~34˚26'，东西约80 km，南北约110 km，流域坡度17.4%，平均高程37.2 m，流域不对称系数0.66。东北与江苏省接壤，西与河南省毗邻，涉及淮北、宿州、徐州、商丘4个市级行政区，流域面积2626 km²。经淮北市城区西南部绕城而过，由濉河入新汴河。河道从瓦子口到新汴河段长62.1 km，多年平均径流量为1.9亿m³。利用萧濉新河流域下游符离集站观测的1973~2016年径流序列，分析流域水文特征的突变性。

萧濉新河流域多年平均降水量为756.4 mm，降水年际变化较大，年内分配不均匀，最大年降水是最小年降水的1.6倍，多年平均汛期(6~9月)降水量占年降水总量比例超过60%；多年平均蒸发量为890.4 mm。多年平均气温14.6℃，年际间变化不大。

3. 研究方法

1) T方法

Tukey提出了基于“中位数”的方法进行分析和剔除突变项的方法，先产生一个曲线的平滑估计，再把它从数据中去掉，若差值大于预定值，则剔除或用内插值替代。具体说明如下：设原始序列为{x_i}(i = 1,2,…,n)，取{x_i}(i = 1,2,…,5)的中位数 ${x^{\prime }}_3$ ，取{x_i}(i = 2,3,…,6)的中位数 ${x^{\prime }}_4$ ，依次类推，得到新序列{ ${x^{\prime }}_i$ }(i = 3,4,…,n−2)；以类似的方法在新的序列相邻的三个数据中选取中位数构成序列{ ${x^{″}}_i$ }(i = 4,5, …,n−3)；最后，按以下方式构成新序列{ ${x^{‴}}_i$ }(i = 5,6,…,n−4) ‎[4] ；

${x^{‴}}_i=\frac{1}{4}{x^{″}}_{i-1}+\frac{1}{2}{x^{″}}_i+\frac{1}{4}{x^{″}}_{i+1}$ (1)

分析序列 ${\left\{x_i-x^{‴}\right\}}_i\left(i=5,6,\cdots ,n-4\right)$ ，绝对值大于预定值δ处可能有突变项，剔除或以内差值替代。

2) 有序聚类分析法

有序聚类分析法是一种统计的估计方法，通过统计分析得出水文时间序列可能的突变点，再结合实际情况进行具体分析。其实质是使同类之间的离差平方和较小而类与类之间的离差平方和较大，从而求最优分割点。设分割点为t，计算方法为 ‎[5] ：

$S_t={{\displaystyle {\sum }_{i=1}^t\left(y_i-\overline{y_t}\right)}}^2+{{\displaystyle {\sum }_{j=t+1}^N\left(y_j-\overline{y_{N-t}}\right)}}^2$ (2)

式中，S_t为序列离差平方和与突变点后的离差平方和之和；N为序列总数， $\overline{y_t}$ 为分割点t之前样本的均值； $\overline{y_{N-t}}$ 为分割点t之后的样本均值，最优分割点为使得序列离差平方和 $\overline{y_t}$ 最小的点。

3) Mann-Kendall法

Mann-Kendall检验法是由Mann (1945年)和Kendall (1975年)提出的非参数检验方法，此法不需要样本遵从一定的分布，不受少数异常值的干扰，计算简便，目前被很多学者广泛应用于水文的时间序列趋势分析上。基本原理是：原假设 $H_0$ 为时间序列 $\left(Y_1,Y_2,\cdots ,Y_N\right)$ 是N个独立的随机变量同分布的样本；备择假设 $H_1$ 是双边检验，对于所有的 $i,j\le N$ ，且 $i\ne j$ ， $y_i$ 和 $y_j$ 的分布是不相同的，检验的统计变量S按下式计算 ‎[3] ：

$S={\displaystyle \underset{i=1}{\overset{N-1}{\sum }}{\displaystyle \underset{j=i+1}{\overset{N}{\sum }}\mathrm{sgn}\left(Y_j-Y_i\right)}}$ (3)

其中 $\mathrm{sgn}\left(Y_j-Y_i\right)=\{\begin{array}{l}1\left(Y_j-Y_i>0\right)\hfill \\ 0\left(Y_j-Y_i=0\right)\hfill \\ -1\left(Y_j-Y_i<0\right)\hfill \end{array}$ 。

统计变量S服从正态分布，统计变量S的均值与方差为：

$E(S)=0$ $E\left(S\right)=0，vas\left(S\right)=\frac{N\left(N-1\right)\left(2N+5\right)}{72}$ (4)

标准统计变量： $Z=\{\begin{array}{l}\frac{S-1}{\sqrt{vas\left(S\right)}},S>0\hfill \\ 0,S=0\hfill \\ \frac{S+1}{\sqrt{vas\left(S\right)}},S<0\hfill \end{array}$ 。

在给定的α置信水平上，如果 $\left|Z\right|\ge Z_{1-\alpha /2}$ ，则原假设是不可接受的，即在α置信水平上，时间序列数据呈上升或下降趋势，且Z大于0时呈上升趋势，Z小于0时呈下降趋势。当 $\left|Z\right|\ge 1.96$ 时，表示通过了95%的显著性检验。

4) 滑动t检验

滑动t检验是通过序列中两段子序列的均值有无显著差异来验证2个总体均值的差异性，如果差异性超过一定的显著水平，即表示水文序列发生突变。对于序列{x_i}(i = 1,2,…,n)，将某一时刻设置为基准点，该序列前后两段子序列x₁、x₂，对应的样本个数、平均值、方差分别为n₁和n₂、 $\overline{x_1}$ 和 $\overline{x_2}$ 、 $s_1^2$ 和 $s_2^2$ 。滑动统计量t的计算公式 ‎[3] ：

$t=\frac{\overline{x_1}-\overline{x_2}}{s\times \sqrt{\frac{1}{n_1}+\frac{1}{n_2}}}$ (5)

$s=\sqrt{\frac{n_1{s}_1^2+n_2{s}_2^2}{n_1+n_2-2}}$ (6)

遵循自由度 $v=n_1+n_2-2$ 的t分布。

5) Yamamoto法

从气候信息与气候噪声两部分来讨论突变问题，文中用以分析水文要素的突变性。为了检验突变转折年是否达到气候突变的标准，对其计算信噪比R_SN ‎[14] ：

$R_{SN}=\frac{\left|\overline{X_1}-\overline{X_2}\right|}{S_1+S_2}$ (7)

式中： $\overline{X_1}$ 、 $\overline{X_2}$ 和 $S_1$ 、 $S_2$ 分别为转折年前后两阶段的平均值和标准差，R_SN为信噪比，当R_SN > 1.0则认为有突变发生，当R_SN > 2.0时，则认为有强突变发生，最大信噪比所对应的年份为突变出现年。这种检验突变的方法称为Yamamoto法。

6) Pettitt法

Pettitt法也是一种非参数检验方法，样本容量为n的序列{x_i}(i = 1,2,…,n)，对应的秩序列为{r_i}(i = 1,2,…,n)，构成统计量：

$S_k=2{{\displaystyle {\sum }_{i=1}^{k}r}}_i-k\left(n+1\right),\text{k}=\text{1,2,}\cdots ,n$ (8)

S_k为检验结果，若在某年出现突变，则 $S_E=\max \left|S_k\right|$ ，E点为突变点，计算统计量 $P=2\exp \left\{-6S_E{}^2/\left(n^3+n^2\right)\right\}$ ，若P ≤ 0.5，则认为突变点在统计意义显著 ‎[8] 。

4. 研究结果

4.1. 符离集水文站年径流序列变化特征

从图1符离集站1973~2022年水文站年径流序列可以看出，年降水和年径流的趋势基本一致，具有较好的相关性，整体呈现下降趋势 ‎[15] 。分析原因，主要由于人类活动的影响，流域内显著的人类活动为煤矿开采，尤其是80年代以来，流域内煤矿大量开采，出现了大小不一的塌陷区，对流域的降水径流关系产生一定的影响。其中历年最大径流量为5.01亿m³，出现在1982年，最小径流量为2002年的0.02亿m³，最大径流量为多年平均径流量的3.7倍。在2007年和2008年出现2个高峰值，但从2009年之后持续下降。从符离集站年径流实测序列趋势特征来看，径流序列突变点发生在20世纪80年代、21世纪10年代。

Figure 1. The relation curve between annual precipitation and annual runoff at the departure station

图1. 符离集站年降水量–年径流关系曲线

4.2. 突变检验分析研究

Figure 2. The result of annual runoff clustering method in the runoff collection station

图2. 符离集站年径流有序聚类法检验结果

由T方法，认为存在3个可能突变的时间，分别是1987年、1995年和2001年。符离集站1973~2022年年径流量有序聚类法检验结果见图2，从图中可以看出，从21世纪10年代初开始年径流序列离差平方和开始出现低于往年数值，虽后有小幅上升，但年径流量序列离差平方和在2009年开始呈现窄谷底的变化过程，2010年数值最低，最小点较显著，t = 38时，S_t最小，这与年径流演变特征也相吻合，2010年可能是对应的变异点，说明径流突变不是一个突然发生的跳跃过程，而且持续的现象。

符离集站年径流Mann-Kendall检验法结果见图3，设置显著性水平α = 0.05，U_α = ±1.96，可以看出，自1975年开始至1983年，U_Fk低于0，之后有3年小幅波动大于0，从1987年开始直至2006年，U_Fk一直低于0，2007~2012年有小幅大于0，2013~2022年又小于0，说明整体呈现减少趋势。按照突变点的特征，年径流量的几乎每次变化都会产生潜在的突变点，分别为1975年左右由增加到减少的突变点、1982年为减少到增加的突变点、1987年为增加到减少的突变点、1995年为减少到增加的突变点、2001年为增加到减少的突变点、2010~2011年为增加到减少的突变点，这一突变结果在t检验法和有序聚类法中均有体现。

Figure 3. Mann-Kendall test of annual runoff in Fuli Sluice station

图3. 符离集闸站年径流量的Mann-Kendall检验

Figure 4. The Pettitt statistic curve of annual runoff from the collecting station

图4. 符离集站年径流的Pettitt统计量曲线

采用滑动t检验法对符离集站年平均径流量序列进行检验，按照公式计算统计量，1973~2022年符离集站年径流量滑动t统计量有4处超过了α = 0.01显著水平的临界值2.878，分别是1987年、1995年、2005年和2010年。说明符离集站1973~2022年年径流量经过了4次较明显的突变，最显著的是2010年的突变，滑动统计量超过了α = 0.01显著水平，而其他3次突变统计量略微超过α = 0.01显著水平，说明在2010年前后径流量变化趋势发生了较大的变化，这与有序聚类法结论一致；其次是1987年前后年径流量变化趋势也发生了很大的变化。此处结论与上述方法表述相一致。

而由Yamamoto法，RSN = 0.10 < 1.0，此法认为没有突变，这与以上4种方法计算结果差距较大，说明此法不适用萧濉新河实测径流突变性研究。

Pettitt法突变检验结果为2010年突变点，突变点统计量P值为0.49，说明此突变点具有显著的统计意义，图4为符离集站S_k随时间变化曲线图。

6种检验方法对符离集站1973~2022年径流系列可能存在的突变点检验结果见表1。

Table 1. Statistical table of mutation test

表1. 突变性检验统计表

根据上述分析，结合萧濉新河流域社会经济发展实际情况，分析得出符离集年径流系列的突变点分别为1987年、1995年、2010年。可见Mann-Kendall检验法用于此流域突变性检验分析得到的结果最精确，其次是滑动t检验、T方法、Pettitt法，精度最差的是Yamamoto法。几种方法检测的突变点的意义具有差异性，T方法的优点是不受极端值影响，以水文序列与中位数差异的大小作为区分突变的判定条件。有序聚类法检测的突变点反映的是年径流序列的聚类特性，以离差平方和的值将结果分为前后两类。滑动t检验法检测突变性反映的是年径流序列的显著差异性，需要进行t分布检验，本次研究中检测出3个突变点，整体与其他检测方法结论相符。Pettitt法检验出的突变点也需要进行显著性检验，但缺点是不适用与检验突变点较多的水文序列。Mann-Kendall检验法检测的突变点反应的是径流序列的综合突变性质，也反映出年径流序列的趋势性，在本次研究中分析结果结合了有序聚类法、滑动t检验、T方法、Pettitt法的所有结果。

水文序列分析结果表明，萧濉新河流域年径流量整体呈现下降趋势，结合降雨分析来看，结论一致。从20世纪60~70年代开始，经济快速发展，对水资源需求迅速增加，修建了大量的水利工程，下垫面条件改变了流域的水文循环特性，20世纪90年代之后人类活动对区域水资源影响更大。因此，符离集站1973~2022年年径流量的突变时间的确定结果基本符合区域水文情势发展的特征。

5. 结论

研究径流序列的突变性对于探讨变化环境下水文循环的演变规律具有重要的意义和价值。在以往的水文序列突变检验中，往往使用单一的检验方法判定突变点，但由于每种突变方法均有优缺点，不能完全反映出序列的突变特性，因此可采用多种突变方法进行检验，以确定结果的准确性。本文以萧濉新河流域1973~2022年年径流量为例，采用T方法、有序聚类法、滑动t检验、Mann-Kendall检验法、Yamamoto法、Pettitt法等方法对比分析水文序列突变检验结果的异同，同时结合趋势分析，得出年径流序列的突变点为1987年、1995年、2010年，结论更符合实际情况。

参考文献