变化环境下汀江流域非一致性水文频率分析

吕福水

(广东省水文局梅州水文分局，广东 梅州 514000)

1 概述

水文频率计算是水文分析计算的主要内容之一，是各类涉水工程规划、设计确定工程规模和管理决策的主要依据[1]。现行的水文频率分析方法的一个基本前提是水文系列满足一致性假设，即水文极值的概率分布或统计规律在过去、现在和未来保持不变[2-3]。当前，国内外进行非一致性水文频率计算主要基于两个途径[4]：一是通过还原(还现)方法将水文极值系列还原(还现)到过去(现在)某一状态，使其满足一致性要求；二是直接基于非平稳极值系列进行水文频率分析，如时变矩模型[5-6]和混合分布模型[7-8]等。混合分布模型认为整个非一致性的样本系列可以分成若干个一致性的子系列，先估计各子系列的分布函数(子分布)，再对各子分布进行加权综合，以获得一个可综合考虑各子系列特征的分布函数，用以描述极值系列的总体分布特征[9-10]。本文采用基于趋势分析和基于两P-Ⅲ混合分布的非一致性水文频率，分析汀江流域溪口站最大流量及最枯月平均流量设计值，分析传统水文频率计算与非一致性水文频率计算之间的差异，为汀江流域防洪调度提供参考依据。

2 流域概况

汀江为韩江的一级支流，发源于福建宁化木马山，在永定棉花滩流入梅州境内的大埔青溪，于茶阳汇小靖河和漳溪河，在三河坝与梅江、梅潭河汇合流入韩江，流域集雨面积为11 802 km2，河长为323 km，平均坡降为0.24%。溪口水文站位于广东省梅州市大埔县青溪镇溪口村，设立于 1958年，集水面积为9 228 km2，占汀江流域的78%，历史最大洪峰流量为8 140 m3/s，历史最枯月平均流量为21.2 m3/s。汀江曾经是闽粤主要水路交通线，大埔茶阳至福建永定石下坝可通小机船，其天然落差大，水力资源丰富，在上游永定石下坝已建有大型棉花滩水库，总库容为20亿m3，装机60万kW；在下游大埔青溪，已建中型青溪水库，总库容为7 468万m3，装机14.4万kW；在青溪水库下游约10 km、茶阳镇上游约1 km处建有茶阳水电站，总库容为3 723万m3，装机2.24万kW。汀江流域水系及站网分布如图1所示。

图1 汀江流域水系及站网分布示意

3 研究数据与方法

3.1 数据来源

本文采用溪口水文站1959—2021年最大洪峰流量序列及1963—2021年最枯月平均流量序列分析变化环境下水文频率变化。受水利工程建设影响，该站年最小流量序列值极小且较为接近，无法反应水利工程建设影响，因此采用最枯月平均流量，通过时间尺度上的拉长可以更好地反应溪口站的枯水变化规律。数据来源于广东省水文局梅州水文分局。

3.2 基于趋势分析的非一致性水文频率分析

将序列进行突变检验得出分割点τ，将样本容量为n的序列分割成x1，x2…，xτ和xτ+1，xτ+2…，xn，计算前后序列的均值EXa和EXb，则平稳状态下的震动中心Ex可表示为：

Ex=A·EXa+(1-A)·EXb

(1)

A为两序列样本权重，可通过前后序列样本均值EXa和EXb来确定，并考虑近期洪水对未来洪水影响较大的因素进行权重分配。权重A公式表示为：

(2)

由于Ex可以求出，则确定性成分可表示为：

y(t)=f(t)-Ex

(3)

经过一致性修正后的系列可表示为：

y′(t)=x(t)-y(t)

(4)

其中，f(t)序列包含趋势成分对应值，x(t)为原始序列。

3.3 基于两P-Ⅲ混合分布的非一致性水文频率分析

以τ为分割点将样本容量为N的序列分割成两段，变异点之前的序列为x1，样本长度为N1=τ，其概率密度函数为f(x1)。变异点之后的长度为N2=N-τ，其概率密度函数为f(x2)。序列x服从混合分布f(x)。其概率密度函数表示为[9-11]：

f(x)=αf(x1)+(1-α)f(x2)

(5)

其中α为权重系数；在两P-Ⅲ混合分布中，将f(x1)、f(x2)代入式(5)。

(6)

(7)

4 年最大洪峰流量及最枯月平均流量变异特性分析

4.1 Mann-Kendall检测法

采用Mann-Kendall检测法[12]对溪口水文站年最大洪峰流量及最枯月平均流量进行突变点检验(见图2、图3)。从图2可以看出，溪口水文站年径流量UF统计量和UB统计量相交于2001年，并通过显著性0.05 的检验，呈显著下降趋势。从图3可以看出，溪口水文站最枯月平均流量近53 a未发生明显变异。

图2 溪口水文站年最大洪峰流量Mann-Kendall检验示意

图3 溪口水文站最枯月平均流量Mann-Kendall检验示意

4.2 滑动t检验法突变检验

滑动t检验法[13]是通过t检验法的原理，逐一对序列进行t检验。t服从t(n1+n2-2)分布，显著性水平α采用0.01，通过查阅t分布表获得临界值tα/2。当t>tα/2时，拒绝原假设，说明其存在显著性差异；当t

图4 溪口水文站年最大洪峰流量滑动t检验示意

图5 溪口水文站最枯月平均流量滑动t检验示意

5 非一致性水文频率分析

根据Mann-Kendall检测法及滑动t检验法突变检验结果，溪口水文站年最大洪峰流量序列在2001年发生变异，最枯月平均流量未发生明显变异。因此，本文仅对溪口水文站年最大洪峰流量序列采用基于趋势分析的非一致性水文频率分析及基于两P-Ⅲ混合分布非一致性水文频率计算，并分析两种非一致性水文频率计算方法成果的合理性。

5.1 基于趋势分析的非一致性水文频率计算

年最大洪峰流量序列最优分割点为2001年，将序列分为两个子样本；样本均值为EX1=3 900 m3/s和EX2=1 840 m3/s，采用公式(3)和(4)进行序列的一致性修正，溪口水文站年最大洪峰流量修正前后设计值见表1，序列修正后Mann-Kendall检验见图6，修正前后频率曲线见图7～图8，序列修正后Mann-Kendall检验近45 a未发生明显变异。

表1 基于趋势分析的非一致性水文频率修正前后设计值

图6 溪口水文站年最大洪峰流量修正后Mann-Kendall检验示意

图7 溪口水文站年最大洪峰流量修正前频率曲线示意

图8 溪口水文站年最大洪峰流量修正后频率曲线示意

从表1可以看出，基于趋势分析的非一致性水文频率计算修正前后均值较为接近，修正后100 a一遇设计值减小6.51%，50 a一遇设计值减小6.59%，20 a一遇设计值减小6.51%，10 a一遇设计值减小6.13%，5 a一遇设计值减小5.10%。受2001年棉花滩水电站调蓄影响，溪口站年最大洪峰流量设计值修正后较修正前均减小。

5.2 基于两P-Ⅲ混合分布非一致性水文频率计算

采用L—矩对汀江溪口水文站年最大洪峰流量进行参数计算，并用信息熵进行权重系数α的估算，得出各参数见表2。将各参数代入公式(6)和公式(7)生成两P-Ⅲ混合分布新密度函数，绘制频率曲线(见图9)，分析在混合分布下的重现期。从表3可以看出，基于两P-Ⅲ混合分布非一致性水文频率计算，修正后100 a一遇设计值减小2.88%，50 a一遇设计值减小2.24%，20 a一遇设计值减小2.51%，10 a一遇设计值减小2.51%，5 a一遇设计值减小2.75%。

表2 两P-Ⅲ混合分布各参数值

表3 基于两P-Ⅲ混合分布非一致性水文频率修正前后设计值

图9 年最大洪峰流量基于两P-Ⅲ混合分布非一致性水文频率曲线示意

对比两种不同的非一致性水文频率计算修正方法，图9水文频率拟合结果比图8更好，更接近高洪期设计值。因此，溪口水文站年最大洪峰量序列设计值宜采用基于两P-Ⅲ混合分布非一致性水文频率计算修正。

6 结语

1)通过Mann-Kendall突变检验及滑动t检验法突变检验，溪口水文站年最大洪峰流量序列在2001年发生突变，呈显著下降趋势，而最枯月平均流量序列未发生明显变异。2001年，溪口站上游棉花滩水电站建成并开始蓄水，对洪峰流量的削峰作用明显。平水期及枯水期，电站蓄水发电较为规律，对最枯月平均流量序列影响较小。

2)采用基于趋势分析的非一致性水文频率计算并对比修正前后溪口水文站年最大洪峰流量设计值差异，结果表明，修正后100a一遇设计值减小6.51%，50a一遇设计值减小6.59%，20a一遇设计值减小6.51%，10a一遇设计值减小6.13%，5a一遇设计值减小5.10%。受2001年棉花滩水电站调蓄影响，溪口水文站年最大洪峰流量设计值修正后较修正前均明显减小。

3)采用基于两P-Ⅲ混合分布非一致性水文频率计算并对比修正前后溪口水文站年最大洪峰流量设计值差异，结果表明，修正后100a一遇设计值减小2.88%，50a一遇设计值减小2.24%，20a一遇设计值减小2.51%，10a一遇设计值减小2.51%，5a一遇设计值减小2.75%。

4)基于两P-Ⅲ混合分布非一致性水文频率曲线比基于趋势分析的非一致性水文频率曲线拟合成果较好，前者在高洪期的设计值更具有代表性。

5)溪口水文站年最大洪峰量序列受上游棉花滩水电站调蓄影响，已失去了水文资料的一致性，本文采用两种不同方式对序列进行修正，可为汀江流域水利工程建设及防洪调度提供参考依据。