水文时间系列的周期分析方法探讨
1 前言
水文现象由于受多种气候因素的影响，具有随机性，其时间系列所隐含的震荡周期的识别、判定是一个重要但又较为困难的问题。解决方法从最初的离散周期图、方差周期图、方差分析过渡到谱分析、小波分析，经历了从时间域到时频域的发展过程，特别是1965年傅立叶变换的出现以及上世纪80年代计算机技术的飞速发展，使时频域分析走向实用并迅速拓展。
2 方法原理
2.1 最大熵谱分析方法
谱分析方法把时间系列看成是多种不同频率的规则波（正弦波或余弦波）迭加而成，比较不同频率波的方差大小，从而找出波动的主要周期，对一时间系列的谱分析，主要有功率谱、最大熵谱分析方法等。其中最大熵谱具有突出的分辨能力，峰值偏离小，提取的主次周期更符合实际。
对一组系列值，熵定义为：
ωωd S H ?+∞
∞
-=)(ln （1）
式中：ω为频率，)(ωs 为谱密度。
由于功率谱和自相关函数)(j r 互为傅立叶变换，即：
ωωωd e s j r j i ?+∞
∞
-=)()( （2）
因为熵值极大输出功率也达到极大，现在要解决的关键问题就是如何利用)(j r 去估计谱密度)(ωs ，使其满足式（2），且式（1）中的熵谱为最大。利用泛函分析中拉格朗日乘子法可以证明，满足以上条件的熵谱密度为：
2
1
)(2
00
1)(∑=--=
k k k
i k k k H e
a S ωσω （3）
式中0k 为自回归的阶数；)
(0k k
a 为自回归系数；2
0k σ为预报误差方差估计。
将计算得到的不同波长（可换算成周期）的最大熵谱密度绘成图，如果谱密度有尖锐的峰点，其对应的周期就是序列存在的显著周期。 2.2 小波分析方法
小波分析是一种时、频多分辨率分析方法，具有时频局部化功能，可以对函数和信号系列进行多尺度细化分析，以分析不同尺度（周期）随时间的演变情况。
小波分析不仅能将水文时间系列的频率特征在时间域上展现出来，分析出其主要周期，而且能清晰地给出各种时间尺度（周期）的强弱和分布情况以及旱涝变化趋势和突变点。方法原理简述如下。
小波函数是指具有震荡特性、能够迅速衰减到零的一类函数，即：
0)(=?
-∞
∞
+dt t φ （4）
目前有很多函数可以选用，本文采用“墨西哥帽”小波，)(t φ通过伸缩和平移构成一族函
数系：
)(
)(2
1
,a
b
t a
t b a -=-φφ （5） 式中：R b ∈，R a ∈，0≠a 。)(,t b a φ为连续小波，称为母小波，对于时间系列函数)(t f ，其小波变换的连续形式为：
dt a
b
t t f a
b a W f )(
)(),(2
1
-=?
+∞
∞
--φ （6） 其中)(t φ为)(t φ的共轭函数，),(b a W f 称小波变换系数。在实际工作中，时间系列常常是离散的，如t ?为取样时间间隔，n 为样本量，则式（6）的离散形式为：
)(
)(),(1
2
1
a
b
t i t i f t a
b a W n
i f -???=∑=-φ （7） ),(b a W f 随参数a 和b 变化，作出以b 为横坐标，a 为纵坐标的关于),(b a W f 的二维
等值线图，正的小波系数对应于偏多期，负的小波系数对应于偏少期，小波系数为零对应着突变点；小波系数绝对值越大，表明该时间尺度变化越显著。通过小波系数的分析，可识别降水量时间系列多时件尺度演变特性和突变特征。
3 实例分析
3.1 最大熵谱分析方法
根据上述原理，运用最大熵谱研究1956-2000年共45年降水量时间序列，判断1x ，
2x ，…，n x 显著周期。
(1) 建立适当阶数的自回归模型，并根据最终预报误差（FPE ）最小原则，确定模型的最佳阶数0k ，最终预报误差计算见式（8）
2
)(k k
n k n k FPE σ-+=
（8）
(2) 用递推模型计算最终的自回归系数。
∑∑=-=++=n
t t t n
t t t x x
x
x a 2
2122
1
11
)(2
（9）
∑∑∑∑∑∑+===+-----+===+-----++??
????-+-?
?
?
???--=n k t k j k j j k t kj k t j t kj t n k t k j k j j k t kj k t j t kj t k k x a x x a x x a x x a x a 2112
1122111111
)()())((2 （10）
(3) 利用（3）式计算最大熵谱，由于降水系列为离散值，只能采用离散形式计算熵谱值，因此计算出的最大熵谱的离散形式：
2
1
)(2
1)(20
00
0)sin()cos(1)(??????+??????-=
∑∑==k k k k
k k k k k H m lk a m lk a l S ππσ （11）
上式中，频率取m
l
l 22πω= ( l =0，1，2，…，m)，m 为选取的最大波数，在序列样本量不大时，m 通常取
2
n
。 根据上述原理，最终计算出最佳阶数为6以及不同周期的最大熵谱值，如表1，并绘制最大熵谱图，见图1。
表1 年降水最大熵谱表
图1 年降水最大熵谱图
根据降水系列的最大熵值计算图、表，可以得出谱密度曲线有三个明显尖锐峰点，为
7.18)4(=H S ，其次为76.9)19(=H S ，16.7)11(=H S ，对应的周期11)4(=主T 年，
32.2)19(=次T 年，4)11(=次次T 年。因此，降水主要存在11年左右和2~4年左右的周期
震荡。
3.2 小波分析方法
由于采用小波变换来分析时间系列的多时间尺度演变规律，需较长资料系列，而年降水系列较短，某雨量站却有较完整的长系列（1910-2000年）月降水资料，同时与区域降水有较好的一致性，所以本文对某雨量站月降水进行小波变换分区域降水的多时间尺度演变特征。
为消除边界效应，将1910-2000年共91年1092个月降水资料进行距平处理，同时采用对称延伸法将两端数据外延；再进行数据处理后选用“墨西哥帽”母小波，对上述资料小波变换分析，得出各不同时间尺度下的小波系数，见图2所示。
图2 某雨量站降水系列的小波变换系数时频图
由图2分析得到：区域降水周期存在持续时间最长、震荡最强烈的是6~12年，主要的多雨期在上世纪20年代、50年代和80年代的中期，少雨期主要在上世纪10年代中期、40年代初期、60年代后期。另外还存在持续时间较短、较为明显的20~30年、2-6年周期震荡规律。
4 结论
本文应用最大熵谱和小波变化两种方法分析区域降水时间系列的周期变化规律，两者分析结果基本一致，主要结论如下：
（1）从二十世纪初到二十世纪末，我省年降水主要存在阶段性的20-30年、6-12年、2-5年左右的主要周期振荡。
（2）20-30年左右的周期变化从上世纪初持续到60年代，振荡最强；5-12年周期明显存在于上世纪20年代至80年代末期；2-5年周期振荡主要存在于在上世纪10年代初至20年代末、50年代初、90年代末期。
（3）上世纪50年代以来，我省降水主要表现为10年左右的周期振荡，60年代末、80年代初、90年代中后期表现出较强烈5-7年的周期振荡，值得注意的是从80年代开始，周期振荡更趋频繁。
由于降水的变化规律复杂多变，以上分析结果虽然比较接近广东省降水变化的客观规律，运用功率谱和小波分析的方法进行降水变化规律的研究还需进一步完善，有待进一步深入。
参考文献
〔1〕张吉峰. 谱分析在测定时间系列周期中的应用[J]，预测. 1994，13（4）：40-45.
〔2〕崔锦泰.小波分析导论[M]，西安：西安交通大学出版社，1995.
〔3〕胡昌华，张军波，夏军等. 基于Matlab的系统分析与设计——小波分析[M]. 西安：西安电子科技大学出版社，2000.
〔4〕纪忠平，谷德军.近百年来气候变化的多时间尺度分析[J]，热带气象学报，1999（2）：48-55.
〔5〕Farge M. Wavelet transforms and their applications to turbulence. Annu Rev Fluid Mech, 1992, 24: 395-457