非一致性条件下嘉陵江上游水文频率分析

非一致性条件下嘉陵江上游水文频率分析
郭明;杨志勇;周政辉;袁丽婷
【摘要】受全球气候变化及高强度人类活动的影响,流域暴雨特征、产汇流特性发生了显著改变,破坏了水文频率分析的基础---水文系列统计的一致性.以嘉陵江流域上游的略阳站为例,应用Hurst系数法、Lee-Heghinian检验法、有序聚类法、滑动T检验法等多种方法,对略阳站1940-2016年的年最大日流量序列统计的一致性进行了检验,识别了其主要的变异点,计算了对应不同时期气候及下垫面特征的水文系列统计参数.结果表明:略阳站1940-2016年的年最大日径流序列呈不显著递减趋势,递减速率为-68.25 m3/(s·10a);综合多种非一致性检验结果,1994年为可信度最高的变异年份;经一致性修正,对应近期气候与下垫面条件下不同频率洪水设计值较历史值有显著变化,百年一遇、五十年一遇及二十年一遇洪水的设计值均有所降低,变幅分别为-2.94%、-6.51%及-12.49%.
【期刊名称】《华北水利水电学院学报》
【年(卷),期】2018(039)001
【总页数】6页(P42-47)
【关键词】非一致性;变异识别;水文频率;嘉陵江;略阳站
【作者】郭明;杨志勇;周政辉;袁丽婷
【作者单位】中国水利水电科学研究院,北京100038;流域水循环模拟与调控国家重点实验室,北京100038;中国水利水电科学研究院,北京100038;流域水循环模拟与调控国家重点实验室,北京100038;中国水利水电科学研究院,北京100038;河海
大学,江苏南京210098;中国水利水电科学研究院,北京100038;流域水循环模拟与调控国家重点实验室,北京100038
【正文语种】中文
【中图分类】TV12
气候变化和频繁的人类活动,破坏了水文序列原本的一致性特征,使传统意义的水文频率计算方法遭受了极大的挑战[1]。

因此,研究变化环境下非一致性水文频率的计算方法,具有重要的理论意义及实际价值。

目前,我国对非一致性水文频率进行计算的传统方法源自20世纪80年代初我国的第一次全国水资源评价,是基于“还原”和“还现”的计算方法[2-3]。

该方法过程简单,计算容易,但存在还原精度较低、资料获取困难、有些项不能还原、成果使用不方便等问题。

因此,后来逐步转向直接基于概率论与数理统计的分析方法,如时变矩模型[4-5]、混合分布法[6]和条件概率分布法等[7]。

时变矩模型多用于分析计算具有趋势变异的水文频率,混合分布法及条件概率分布法多用于分析计算具有跳跃变异的水文频率。

如:谢平等[8]基于极值同频率法分析了东岭下站非一致性年径流序列,并确定了跳跃变异点;李伶杰等[9]以渭河流域具有跳跃变异的林家村、张家山和状头3个水文测站的年平均流量系列为研究对象,应用条件概率模型,计算效率系数,得出条件概率模型拟合效果较优,是对具有跳跃变异的非一致性水文序列频率计算较优的方法;王军等[10]应用混合分布法对淮河流域非一致性降雨系列进行分析,将设计标准值与经验标准值进行对比,得出混合分布模型的拟合度良好。

尽管针对不同方法的非一致性水文分析研究较多,但采用多种非一致性检验法综合识别序列非一致性成分并比较不同方法的实用性,从而进行非一致性水文频率计算的研究还较少。

因此,本文基于多种非一致性识别方法,对略阳站1940—2016年最
大日径流序列进行非一致性检验,综合不同检验方法所得结果,提高变异识别结果的可信度,进一步计算非一致性水文频率,以期为略阳站附近的水利工程规划与管理提供水文数据支撑。

1 研究区概况与数据来源
1.1 研究区概况
略阳县位于陕西省西南部,与勉县、宁强县于东南部接壤,是陕西、甘肃、四川3省的交界地带。

略阳水文站位于略阳县城关镇,于1939年建站,系嘉陵江干流上游控制性水文测站,集水面积19 206 km2,多年平均径流量36.3亿m3,年平均降水量860 mm。

嘉陵江上游白水江、西坝一带的年降水量700 mm左右,城关的860 mm左右,东南部的900 mm左右。

全县降雨量最高中心在何家岩,年降水量1 067.9 mm。

年内6—8月降水量较大,为洪水期；5、9、10月为丰水期；其余6个月为平水期。

历史记录大洪水:1981年8月22日,最大洪峰流量8 500
m3/s;1898年历史最大洪水,最大洪峰流量9 500 m3/s[11]。

嘉陵江略阳断面是陕西省引嘉入汉工程的重点取水口,通过隧道向汉江流域调水,以补充汉江径流量,优化丹江口水库入库过程,从而保证引嘉入汉工程作为南水北调中线的后续水源工程,对南水北调中线工程的调水起积极作用。

研究区位置及高程信息如图1所示。

图1 区域概况及站点位置
1.2 数据来源与处理
本文选取略阳站1940—2016年长序列日径流资料,径流过程中缺测3个月以上的年份予以剔除,其中包括1947—1950年及1967—1968年,择取每年年最大日径流数据,构成年最大日径流序列,实测研究序列长度共计71 a。

2 研究内容
2.1 年最大日径流趋势分析
根据略阳站1940—2016年长序列年最大日径流数据,绘制多年变化趋势如图2所示。

由图2可知,略阳站1940—2016年的年最大日径流呈下降趋势,递减率为-68.25 m3/(s·10a)。

取信度水平α=0.05,自由度n-2=70时,rα=0.223,由于实测线性拟合的R2值为0.021,r=0.146,则有r<ra,故未通过显著性检验,即递减趋势不显著。

1940—1970年,略阳站的年最大日径流值在多年均值上下均匀分布;而1971—2016年,该站年最大日径流值小于多年均值的年份比大于均值的年份多,推测在此段时间内可能存在径流量变异现象,导致年最大日径流量在多年均值附近分布不均,变化波动情况较大。

根据5 a滑动平均曲线,略阳站的年最大日径流的递减趋势基本遵循“增—减—增—减—增”的周期性规律。

其中:1940—1962年为第一个增势阶段,1963—1974年为第一个减势阶段,实测数据增减过程总历时29 a;1975—1981年为第二个增势阶段,1982—1997年为第二个减势阶段,实测数据增减过程总历时23 a。

总体来看,略阳站的年最大日径流的每一个增减变化周期持续约30 a。

图2 略阳站年最大日径流变化趋势
2.2 非一致性变异识别方法
本文主要应用Hurst系数法、Lee-Heghinian检验法、有序聚类法、滑动T检验法及Mann-Kendall检验法5种非一致性变异识别方法对略阳站1940—2016年的年最大日径流序列进行非一致性变异识别。

Hurst系数法是1951年由Hurst[12]提出,其原理是利用R/S法计算水文序列的H 指数,H偏离0.5的程度越大,表明持续的变化效应越明显,变异程度相应越大,由此可以根据H指数判定序列是否发生变异以及发生变异的程度。

谢平等[13]通过结合分数布朗运动增量相关函数与Hurst系数之间的关系,提出了基于Hurst系数的水文变异分级方法,并将变异程度划分为无变异、弱变异、中变异、强变异及巨变异
5个等级。

Lee-Heghinian检验法是基于贝叶斯理论的方法。

对序列xt(t=1,2,…,n),在假定总体为正态分布以及分割点τ的先验分布为均匀分布的情况下,推求可能分割点的后
验条件概率密度函数f(τx1,x2,…,xn)[14]。

若点τ满足则该点为最可能的变异点。

有序聚类法以有序分类来推求最可能干扰点τ0,实质上就是推求最优分割点,保证同类间的离差平方和较小,而不同类之间的离差平方和较大[15]。

假设分割点为τ,则
分割前、后的离差平方和分别为Vτ和Vn-τ,总离差平方和为Sn(τ)。

当τ满足时
即为最可能的变异点。

滑动T检验法设变异点τ前、后两序列总体的分布函数分别为f1(x)和f2(x),从总
体f1(x)和f2(x)两序列中分别抽取容量为n1和n2两个样本,并假设f1(x)=f2(x)，构造T统计量[16],T服从t(n1+n2-2)分布。

在给定信度水平下,当T>tT/2时,表明存在显著性差异;当T<tT/2时,则接受原假设,即f1(x)=f2(x)。

对于所有满足
T>tT/2的点τ,选择其中T统计量达到最大值的点,作为所求的最可能变异点τ0。

Mann-Kendall检验法是世界气象组织推荐的非参数检验方法,已广泛应用于气象、水文等领域。

该方法不需要样本遵从一定分布,少数异常值对其结果影响较小,且可
以用于检验序列的变异情况[17]。

计算统计量UF、UB,若UFk>0,表明序列呈上升趋势;反之,为下降趋势。

若UF、UB两条曲线超过临界直线,表明序列上升或下降趋势显著;若两曲线在临界区间内存在交点,则该点对应的时刻即为变异开始的时间,但若该点交于两年份之间,需要进一步采取相应方法选定两年份中具体一点为变异点。

表1 各方法识别结果及适用性情况方法名称变异结果判断标准适用性Hurst系数法根据H指数与分数布朗运动增量相关函数的大小判断变异等级仅能判定序列是
否发生变异及变异等级，不能识别出具体的变异点Lee-Heghinian检验法概率密
度函数最大的点为变异点较适用于均值变异有序聚类法离差平方和最小的点为变异点适用性较普遍，检验结果较好滑动T检验法滑动T统计量最大的点为变异点较
适用于均值变异，长序列检验结果较好Mann-Kendall检验法置信区间内UF、
UB交点为变异点适用性普遍，计算简便
2.3 非一致性变异识别结果分析
首先采用Hurst系数法计算长序列的H指数,结合分数布朗运动增量相关函数C(t)，可确定略阳站1940—2016年的年最大日径流序列有弱变异情况发生。

进一步采用Lee-Heghinian检验法进行变异检验,经验函数的最大值点为最可能变异点,检验结果如图3所示。

由图3可知,局部最大值对应的变异年份为1981年和1985年,因此,1981年和1985年为可能发生变异的年份。

图3 Lee-Heghinian检验法变异检验结果
采用有序聚类法进行变异检验,检验结果如图4所示,满足离差平方和最小的点为最可能发生变异的年份,即1994年。

图4 有序聚类法变异检验结果
采用滑动T检验法进行变异检验,结果如图5所示。

由图5可知,1963年为T统计
量达到最大值的点,且通过0.05显著水平检验,1987年及1994年为局部极值点,也可能是变异发生的年份。

采用Mann-Kendall检验法所得检验结果如图6所示,UF、UB曲线在0.05置信区间内有多个交点存在,初始交点为1994年,可判断变异从1994年开始。

图5 滑动T检验法变异检验结果
图6 Mann-Kendall检验法变异检验结果
不同变异识别方法所得结果有一定差异,多种检验结果汇总识别出的变异年份见表2。

其中,拟合程度最高的变异年份为1994年,有3种变异识别方法均得到该结果。

Lee-Heghinian检验法虽未识别出变异年份1994年,但该方法所识别出的变异年
份1981年和1985年与1994年相距较近,且检验结果在由趋势分析所得1971—2016年可能发生变异的时间段内,不排除相近变异点未能识别出的情况。

综合考虑
4种非一致性变异检验方法所得结果,1994年为变异年份被检验出的概率为75%。

因此将1994年确定为略阳站年最大日径流发生变异年份的可信度较高,可采用该结果进行非一致性水文频率的计算。

表2 各种检验方法检验结果汇总检验方法变异年份滑动T检验法1963Lee-Heghinian检验法1981Lee-Heghinian检验法1985滑动T检验法1987有序聚类法、滑动T检验法、Mann-Kendall检验法1994
3 非一致性条件下水文频率计算
3.1 非一致性水文频率计算原理
我国学者谢平[18]提出:实测水文变异序列由理想的平稳性即一致性成分与随机性成分线性组合而成,通过对变异分割点前后的一致性修正,可获得修正后总体的统计特征参数。

假设变异点为τ,序列被分割为x1,x2,…,xτ和xτ,xτ+1,…,xτ+n两部分,变异前、后序列的均值分别为Exa与Exb,变异的水文序列可表示为[19]:
Ex=α·Exa+(1-α)·Exb。

(1)
式中α为变异点前的权重,由Exa与Exb确定[19]:
(2)
Ex可通过式(1)、(2)联合求出,其中f(t)为确定性趋势成分,因此,随机性成分可表示为[19]:
y(t)=f(t)-Ex。

(3)
则经过一致性修正后的序列可以表示为[19]:
Y(t)=x(t)-y(t)。

(4)
采用Pearson-Ⅲ型概率密度函数[20]对水文序列进行频率分析,从而确定径流序列不同保证率下的水文设计值。

3.2 水文频率非一致性分析
采用本文2.3节确定的变异年份(1994年),将略阳站1940—2016年的年最大日径流序列分割成1940—1994年及1995—2016年两个时段的径流序列。

本文以1940年作为对应历史时期气候及下垫面条件的参照年,2016年作为对应近期相应
条件的参照年,对该序列进行频率分析。

首先,对变异的水文序列进行剔除随机性趋势的还原计算,所得的对应不同时期的气
候及下垫面条件下还原序列的趋势分布情况如图7所示。

由图7可知,修正后对应历史与近期的气候及下垫面条件下的略阳站年最大日径流
序列均已剔除趋势性成分。

图7 修正序列的趋势检验
然后,将修正后的序列重新采用Mann-Kendall检验法进行变异识别,所得结果如图8所示。

由图8可知,在0.05置信区间内,UF、UB曲线无交点,变异成分已经剔除,
略阳站年最大日径流序列满足频率计算的一致性要求,可进行传统的水文频率计算。

图8 修正序列的Mann-Kendall变异检验
最后,分别基于历史与近期的气候及下垫面条件,采用Pearson-Ⅲ法计算略阳站年最大日径流序列修正前、后在不同保障率下的水文设计值,并确定其设计参数,结果见
表3。

基于历史与近期的气候及下垫面条件下略阳站年最大日径流量序列修正前、后的频率曲线如图9所示。

根据表3中实测年最大日径流基于不同气候及下垫面
条件下的修正结果,百年一遇洪水的设计值变化显著,对应历史气候及下垫面条件下
洪水的设计值为5 181.06 m3/s,对应近期气候及下垫面条件下洪水的设计值为5 028.92 m3/s,变幅-2.94%。

近期气候及下垫面条件下不同保障率洪水的设计值较
历史值均有所降低,这一变化与年最大日径流序列呈递减趋势相符。

因此,对于略阳站附近的水利工程的防洪要求,需满足近期气候及下垫面条件下洪水的设计值,可采取降低防洪水位的措施,以减少对水库坝体的压力等。

表3 历史与近期的气候及下垫面条件下略阳站年最大日径流序列修正前、后的不同保障率设计值序列类别ExCvCs年最大日径流量设计值/(m3/s)P=1%P=2%P=5%未修正1493．450．662．245166．814433．123472．00历史条件下1732．310．561．905181．064528．303662．16近期条件下1257．000．772．625028．924233．483204．92
图9 修正前、后不同保证率的设计径流量
4 结语
本文选取嘉陵江上游具有代表性的略阳水文测站,对其年最大日径流序列进行趋势
演变分析及基于多种方法的非一致性识别,并根据其非一致性变异结果进行非一致
性水文频率分析,所得结论如下:
1)略阳站1940—2016年的年最大日径流呈递减趋势,递减速率为-68.25 m3/(s·10
a),但未通过0.05显著性水平检验,递减趋势不显著。

2)经5种非一致性变异识别方法的检验,综合不同方法的适用性及所得检验结果,确定1994年为略阳站的变异年份,变异程度属弱变异。

3)经一致性修正计算后,略阳站对应近期气候及下垫面条件下不同保障率洪水的设
计值均小于历史气候及下垫面条件下的修正值,百年一遇洪水的设计值对应历史与
近期的气候及下垫面条件下的变幅为-2.94%,所得结果符合年最大日径流呈逐渐递减的变化趋势。

4)本文对嘉陵江上游略阳站的年最大日径流序列趋势及非一致性进行了分析,得出
其年最大日径流量呈逐年递减趋势,百年一遇洪水的设计值相应降低。

充分考虑一
致性修正问题,可以降低防洪水位以缓解水库坝体承受的压力,同时也可以有效降低
新修水利工程设施的投资成本。

5)本文仅选取了嘉陵江上游略阳站作为代表站进行分析,样本及数据量较少,分析所得结果可能存在通用性较低的问题,需进一步收集更多的数据资料并开展更深入的研究。

参考文献
[1] 李伶杰,宋松柏,王银堂,等.变化环境下石羊河流域年径流频率计算[J].水
文,2016,36(4):19-25.
[2] 宋承新.径流还原计算的综合修正法[J].水文,1999(2):46-48.
[3] 陆中央.关于年径流量系列的还原计算问题[J].水文,2000(6):9-12.
[4] 杜涛,熊立华,江聪.渭河流域降雨时间序列非一致性频率分析[J].干旱区地
理,2014,37(3):468-479.
[5] 江聪,熊立华.基于GAMLSS模型的宜昌站年径流序列趋势分析[J].地理学
报,2012,67(11):1505-1514.
[6] 冯平,曾杭,李新.混合分布在非一致性洪水频率分析的应用[J].天津大学学报(自然科学与工程技术版),2013,46(4):298-303.
[7] 宋松柏,李扬,蔡明科.具有跳跃变异的非一致分布水文序列频率计算方法[J].水利学报,2012,43(6):734-739，748.
[8] 谢平,张波,陈海健,等.基于极值同频率法的非一致性年径流过程设计方法:以跳跃变异为例[J].水利学报,2015,46(7):828-835.
[9] 李伶杰,康艳,宋松柏,等.条件概率模型在非一致性水文序列频率计算中的应用[J].西北农林科技大学学报(自然科学版),2016,44(8):219-225.
[10] 王军,宁亚伟,胡义明,等.混合分布在非一致性水文频率分析中的应用[J].南水北调与水利科技,2017,15(3):1-4.
[11] 周扬.略阳县地质灾害气象预报预警[D].西安:长安大学,2013.
[12] Hurst H E.Long-term storage capacity of reservoirs[J].Trans
ASCE,1951,116:770-808.
[13] 谢平,陈广才,雷红富.基于Hurst系数的水文变异分析方法[J].应用基础与工程科学学报,2009,17(1):32-39.
[14] 王璨,周秀平,王文圣.窟野河洪水序列变异点综合诊断[J].水电能源科
学,2012,30(7):50-53.
[15] 周园园,师长兴,杜俊,等.无定河流域1956—2009年径流量变化及其影响因素[J].自然资源学报,2012,27(5):856-865.
[16] 雷红富,谢平,陈广才,等.水文序列变异点检验方法的性能比较分析[J].水电能源科学,2007,25(4):36-40.
[17] 魏凤英.现代气候统计诊断与预测技术[M].北京:气象出版社,2007.
[18] 谢平,陈广才,雷红富,等.水文变异诊断系统[J].水利发电学报,2010,29(1):85-91.
[19] 高玉丹.梅江松源河段非一致性水文频率分析[J].广东水利水电,2016(4):8-11.
[20] 张家鸣,陈晓宏,叶长青.Pearson-Ⅲ型频率曲线对负偏水文序列的计算[J].水利学报,2012,43(11):1296-1301.。