一种考虑水文周期性的水库汛期分期方法

王大洋，杜 懿，王大刚，龙铠豪，安 程

(中山大学地理科学与规划学院，广东 广州 510275)

0 引 言

汛期是水库发挥防洪功能的主要时期；因此，汛期的科学划定是保障水库更好发挥防洪调度功能的前提[1-2]。以往学者的研究多为直接对收集的水文资料进行整体分析，筛选出反映汛期特性的要素，如降水、蒸发、径流等，然后直接对总体数据进行汛期研究划分[3- 8]。这种方法虽然较为直接，但未能考虑水文要素的周期性变化规律特性。笼统地对全部要素进行分析可能会掩盖部分因周期性变化产生的现象和趋势，不能全面反映汛期各分期的变化过程。因此，本文在分析水文要素周期性变化规律的基础上，以周期长度为时间尺度，借鉴滑窗理论[9]，研究水库汛期分期问题，旨在为现有的汛期分期研究方法提供新的思路和参考。

1 研究思路及方法

1.1 研究思路

研究首先采用小波分析方法对水文要素进行周期性分析，从而确定水文要素变化的主周期；其次，以谱聚类为基础，对长系列水文要素进行分析，确定非汛期与汛期的界限；之后，以主周期为时间尺度，参考滑窗理论，研究主周期长度时间尺度下，汛期各分期的变化规律；最后，通过各个分期的起止时间，确定最终的汛期各分期。

1.2 小波分析方法

小波分析(Wavelet Analysis)是以Fourier变换为基础，通过利用一簇小波函数系列来近似逼近某一信号序列。其常用于分析水文气象要素的多时间演变特征[10-12]。在周期性研究中，Morlet小波函数应用较为广泛，其小波函数和连续变换表达式为

(1)

(2)

式中，Wf(a，b)为小波变换函数；a、b分别为尺度因子和平移因子。小波方差值可以反映水文要素序列能量波动随时间尺度的变化状态。因此，小波方差图中极大值点对应的时间尺度为水文序列变化的主周期。小波方差表达式为

(3)

1.3 谱聚类算法

谱聚类(Spectral Clustering)算法来源于谱图划分理论。该算法的主要思想是把所有聚类的对象看作是某个维度空间中的点，这些点用边线进行连接。边线上具有权重，距离较远的两个点之间的边线权重较小，距离较近的两点之间的边线权重较大，通过对所有数据点进行切割，使得切割后不同的子图间边线权重尽可能小，而子图内的边线权重尽可能大，从而达到聚类的目的。目前，谱聚类算法在计算机视觉、机器学习、图像识别等领域具有较为广泛的用途[13]。

如果需要将图G划分成两个子图A和B，划分的原则为：子图内的边线权重最大化，而各子图间的边线权重最小化，则可以建立目标函数Ncut。有

(4)

(5)

(6)

式中，cut(A，B)表示子图A和B的边线，也称作“边切集”；此外，若能实现目标函数的最小化，则分割后的结果自然能满足不同类样本间的相似度最小，同时也能满足同一类样本相识度最大。如果需要同时划分成多个子图，则上式可以调整为

(7)

1.4 集对分析

集对分析(Set Pair Analysis)是一种将中国哲学中的对立统一和普遍联系进行有机统一，从全新的角度分析变量之间不确定关系的方法。该方法的核心思想是对不确定系统的两个有联系的集合构建集对，通过对集对的特性进行统一性、差异性、对立性分析，建立集对的同、异、反联系度，从而定量分析变量之间的不确定关系[14-16]。

若有集合X和Y，其均有n项属性表征其特性，则有X(a1，a2，…，an)和Y(b1，b2，…，bn)，则可以构建集对H(X，Y)，那么表征集对H(X，Y)联系度的表达式为

(8)

式中，S为同一性个数；F为差异性个数；P为对立性个数，三者总和为n；I为差异度系数，其取值在(-1，1)区间，根据集合反映的具体现象确定，但有时仅作为表示差异的符号；J为对立系数，取值为-1；若令a=S/n，b=F/n，c=P/n，则有μX～Y=a+bI+cJ；式中，a，b，c分别为集对H(X，Y)的同一度，差异度和对立度。

对于差异度，其可以进行进一步挖掘扩展，如令bI=b1I1+b2I2+…，则可得到更多元的联系度

μX～Y=a+b1I1+b2I2+…+bk-2Ik-2+cJ

(9)

式中，a+b1+b2+…+bk-2+c=1；b1，b2，…，bk-2为差异度分量，表示分量的不同程度，如轻度、中度、重度等；I为不同等级的差异性系数。

2 实例研究

2.1 研究对象及数据

研究以珠江水系红水河上游龙滩水库为对象，龙滩水库位于广西壮族自治区天峨县境内，处中国西南部桂、黔、滇三省交界地带，兼具发电、防洪、灌溉等多项功能。坝址以上流域面积为9.85万km2，占据整个红水河流域面积的71%，有效库容为205.3亿m3，具有多年调节能力。流域内多年平均降水量为760～1 860 mm，多集中在汛期。其中，汛期降水占据全年降水总量的85%以上。

研究的水文要素采用龙滩水库下游天峨水文站1960年～2012年逐日径流量数据。考虑到汛期洪水在时段上的延续性，以“候”为时间单位[17](古代气候学中，一候是五天，三候为一节气，一年二十四节气，对应七十二候)。通过对53 a的逐日径流量数据进行整理分析，得到影响汛期的5个因子，分别为多年平均候流量、多年平均候内最大1日洪量、多年平均候内最大3日洪量、多年平均候内最大3日洪量变差系数Cv、候内年最大洪峰出现次数，共同组成分期指标体系。

2.2 径流周期性分析

径流周期分析是对龙滩水库53 a的年径流序列进行分析，寻找径流周期性变化的规律，确定年径流变化的主周期，以此作为滑动窗口的时间尺度。通过采用Morlet小波分析方法，对龙滩水库年径流序列变换分析，得到径流Morlet小波实部图和小波方差图(见图1、2)。

图1 年径流Morlet小波实部

图2 年径流Morlet小波方差

从图1、2可知：①年径流以22～26 a为尺度的周期震荡信号最为明显，20～25 a为尺度的小波实部图中曲线高度闭合，整个值域内的特征明显；②在小波方差图中，共出现了4个极大值点，分别对应5、10、15 a和24 a。其中，24 a对应的小波方差峰值最大，表明震荡能量最为剧烈，其余3个峰值相差不大。因此，年径流的第一主周期为24 a，以此时间尺度作为汛期分期滑动窗口的“窗口宽度”。

2.3 汛期和非汛期确定

水文要素的特性在汛期和非汛期时域内表现特征是不同的，汛期要素相较非汛期会发生明显的变化。对于水库流域而言，为了能够从容地应对汛期洪水灾害，实现适应性调度，一般采用年度长序列进行汛期和非汛期的划分，汛期起止时间相对固定。汛期和非汛期的水文要素应当归属为两个不同的类别。为此，利用谱聚类的独特优势，对指标体系进行聚类分析，进而确定汛期和非汛期的时间界限。

首先，对分期指标体系进行归一化处理，消除指标之间因数据量级差异对分类结果产生的影响；其次，以一候为一个样本，每个样本由5个指标组成，对应5维空间上的一个向量。之后，分别计算5维空间中相邻两个向量之间的边线权重，并构建目标函数。最后，通过对目标函数进行优化分析，使其达到最小值来确定样本的归属类别。指标体系中共有72候，将此72候进行分割后，对于同一类的样本用相同编号标识，汛期类别序号为2，非汛期为1，见图3。

图3 谱聚类划分结果

从谱聚类划分结果可知，72个候可以归为两类，标识非汛期的候序号为{1～21，63～72}，其对应的时间段分别为1月1日～4月15日，11月11日～12月31日；相应地，标识汛期的候序号为{22～62}，其对应的时间段为4月16日～11月10日。从汛期时间结果来看，相比传统的简单以月为单位的划分情况(常以4月～10月为汛期)，以候为单位则更加精确和科学。从汛期的起止时间和延续长度分析，汛期的起始时间延后了15 d，结束时间延后了10 d，汛期总体持续时间约为7个月，和传统汛期持续时间长度基本保持一致。

2.4 汛期分期划定

通过谱聚类分析，确定了汛期的起止时间，即第22候～第26候。通过小波变换周期分析，确定年径流序列的主周期为24 a。以两者计算结果为基础，引入滑窗理论，从1960年开始，以24 a年为尺度，分别将53 a的候序列进行滑动分组，如第1组为1960年～1983年，第2组为1961年～1984年，依次类推，共计30组，逐组进行汛期时间段(第22～62候)分期指标体系构建，然后以集对分析方法对汛期各个分期进行划定。

针对每一组，计算过程大致可分为3个步骤：

(1)建立划分标准集合Y={Y1，Y2}，其中Yi的形状与分期指标体系一致(汛期有41个候，共5个指标，为41×5的矩阵)，且令Y1中所有元素均为1，则有Y1={1，1，…，1}，令Y2中所有元素均为0，则有Y2={0，0，…，0}；本研究中Y1标识主汛期，Y2标识非主汛期，依据时间确定为前汛期或者后汛期。

(2)建立分期指标体系符号集合M，M的形状与分期指标体系一致，其中元素的确定方法为，以各项指标的均值为标准，对分期指标体系中的各个元素逐项进行判别归类，当元素大于均值时，符号集合中的与该元素相同位置标记为1，反之标记为0。由此可确定所有分期指标体系中的符号元素类型，也就是指标体系符号集合。

表1 1960年～1983年汛期各个候联系度及所属类

(3)构建集对H=(Mi，Yj)，并计算集对的联系度。其中，Mi为分期指标体系符号集合，Yj为分期标准集合。通过将集合Mi中的元素逐一与Y1，Y2进行对照比较，统计符号相同的个数，记作S；统计符号相差为1的个数，记作F。根据下式可以得到2元联系度。式中，因F为差异度系数，根据经验取值为0，通过计算符号集合M的联系度，可以判断每候的归类。即

(10)

式中，S1，S2分别为集合Mi中的元素逐一分别与Y1，Y2进行对照比较时符号相同的个数；F1，F2为相差为1的个数。

通过计算可得每组汛期的41个候的联系度，以第一组1960年～1983年为例，其联系度见表1。

从表1划分结果可知，候序号{22～31}和{49～62}的所属类均为Y2，候序号{32～48}所属类为Y1，因此该组可以将汛期划分为前汛期、主汛期和后汛期共3期。前汛期时域为4月16日～6月5日；主汛期时域为6月6日～8月31日；后汛期时域为9月1日～11月10日。

通过时间尺度进行滑动，可得到每组汛期的分期情况，为了更直观地表现划分结果，以候序列为横坐标，以起始年份为纵坐标(如第1组年份为1960年，代表范围1960年～1983年)，对每组汛期分期情况进行统计，结果如图4所示。图中方框中○表示所属类为Y2的候，方框中●表示所属类为Y1的候。

图4 各组汛期划分结果汇总

从30个组的整体情况看，每组的前汛期结束时间点变化不大，基本在第30、31和32候结束。对于主汛期，情况则较为复杂，存在类似两个Y1中间夹杂一个Y2的情况，如{○，●，○}。为避免单个不同类别掺杂影响划分的情况，规定只有当Y2前后面均连续出现2次及以上的Y1时，如{●，●，○，●，●}，此Y2方可转化为Y1，其余情况均以Y2首次出现的时间为主汛期截止时间。因此，按照此法重新进行修正后得到图5。

图5 各组汛期划分结果修正后汇总

从图5可知，随着年份的向后滑动，主汛期的时间范围有变大的趋势，从第1组的17个候(约102 d)，变化到第21组的23个候(约132 d)，整体而言，主汛期结束的时间点多为第48～58候。对前汛期和主汛期结束时间进行统计分析，可得到表2。

表2 前汛期和主汛期结束候序号统计比例

由表2可知，前汛期在第30候结束时的比例最大，为73.33%，平均候序号为30.37(约为第30候)，其对应时时间为4月30日；主汛期在第53候结束的比例最大，为43.33%，平均候序号为51.06(约为第51候)，其对应的时间为9月15日。因此可以确定龙滩水库的汛期各分期起止时间，前汛期为4月16日～5月31日、主汛期为6月1日～9月15日、后汛期为9月16日～11月10日。

3 结果分析

结合之前有关学者对龙滩水库汛期分期的相关研究，将结果进行汇总比较(见表3)。

表3 不同划分方法研究结果对比

由表3可知，史文海[18]采用分形方法对龙滩水库进行划分时，划定汛期起讫时间的依据是样本序列中第一场暴雨和最后一场暴雨出现的时间，如此进行划定略显粗糙。且该法在进行分期时段划分时是以旬(10 d)为单位，不够精确。模糊集合分析法[19]在进行汛期和非汛期的划分时则是直接采用珠江全流域的划分标准，以4月开始10月结束，未能考虑和捕捉到各子流域在汛期时域上存在的微小差异性。考虑水文周期的划分方法兼顾考虑了划分时间的精准性和子流域间汛期时域存在的差异性情况，先对汛期时域进行划定，在此基础上考虑水文要素的周期性变化规律，进行汛期各分期的划分，从而使得划分结果更科学可信。

4 结 语

(1)本研究首次将水文要素的周期性变化规律引入到水库汛期分期的研究中，以水文要素的周期为时间尺度，借鉴滑窗理论，对时间尺度内的汛期进行划分。对于汛期各子分期划分时，采用能同时反映要素的宏观性和微观联系的集对分析方法。最后，对划分结果进行统计分析，并与之前有关研究进行比较，从而确定龙滩水库的汛期分期情况。

(2)科学的汛期分期可以为更好地计算各分期的设计洪水，提出适应性调度规则，制定科学的水库管理方案提供有力支撑。