小波分析整理

《小波分析整理》由会员分享，可在线阅读，更多相关《小波分析整理（11页珍藏版）》请在金锄头文库上搜索。

1、青海湖地区近青海湖地区近 50 年降水量周期变化分析年降水量周期变化分析 杨沈斌，张弥，吕开龙 1问题问题 青海湖位于青海省东北部的青海湖盆地内，既是中国最大的内陆湖泊，也是中国最大 的咸水湖。近年来有多篇关于青海湖面积受气候变化影响的研究报告，认为温度升高，降水 减少和蒸发量大是造成面积下降的几个重要原因。 为此， 本实验拟采用连续小波分析方法对 青海湖地区年降水量周期变化进行分析，探讨降水量变化与青海湖面积变化的关系。 2. 资料资料 考虑到刚察气象观测站正好处于青海湖边上， 所得数据对反映青海湖周边降水变化具有 较好的代表性。因此，以青海湖西北位置的刚察气象台站 1961-2007 年年

2、降水量资料为例， 对该站年降水资料进行小波分析，获取其周期变化特征。同时，获取了青海湖 1961-2007 年 的面积变化资料。从资料发现，青海湖面积近 50 年呈现下降趋势，平均下降 3 km2/a。图 1 和图 2 分别显示了 1961-2007 年刚察气象站年降水量和青海湖面积的变化， 以及对应的趋势 线。 图 1 1961-2007 年刚察气象台站年降水量及趋势线 图 2 1961-2007 年青海湖面积及趋势线 3. 小波分析方法小波分析方法 小波变换具有多分辨率分析的特点，并且在时频两域都具有表征信号局部特征的能力。 小波变换通过将时间系列分解到时间频率域内， 从而得出时间系列的显

3、著的波动模式， 即周 期变化动态，以及周期变化动态的时间格局（Torrence and Compo, 1998） 。小波（Wavelet） ， 即小区域的波，是一种特殊的、长度有限，平均值为零的波形。它有两个特点：一是“小” ， 二是具有正负交替的“波动性” ，即直流分量为零。小波分析是时间（空间）频率的局部化 分析，它通过伸缩平移运算对信号（函数）逐步进行多尺度细化，能自动适应时频信号分析 的要求，可聚焦到信号的任意细节。小波分析将信号分解成一系列小波函数的叠加，而这些 小波函数都是由一个母小波（mother wavelet）函数经过平移与尺度伸缩得来的。用这种不 规则的小波函数可以逼近那些

4、非稳态信号中尖锐变化的部分， 也可以去逼近离散不连续具有 局部特性的信号， 从而更为真实的反映原信号在某一时间尺度上的变化。 小波分析这种局部 分析的特性使其成为对非稳态、 不连续时间序列进行量化的一个有效工具 （Stoy et al., 2005） 。 小波是一个具有零均值且可以在频率域与时间域内进行局部化的数学函数（Grinsted et al., 2004） 。一个小波被称为母小波（mother wavelet） ，母小波可沿着时间指数经过平移与尺 度伸缩得到一系列子小波。子小波可以通过尺度（s，频率的反函数）函数和时间（n）位置 或平移来描述。 利用一系列子小波， 一个信号可以在不同的

5、时间尺度上进行计算并显示出详 细的特征尺度。 拉伸更大的小波窗口， 使其宽度更大便可以分析时间系列中波动较大的部分 并捕捉大尺度（低频）事件的特征。相反，压缩较小的窗口将包含小尺度（高频）的事件信 息。当信号被子小波相乘，被 s 与 n 唯一的表达，我们可以计算出信号在时间频率域一个具 体位置的系数。如果信号在时间 n 上的谱成分可以与小波 s 比较，那么计算的小波系数具有 相对较大的值。在其它 n 与 s 的组合（如其它的子小波）上都进行这样的计算，那么将会产 生一系列系数（小波变换）来表达信号在时间频率域内的分解。通过这样的变化便可得到时 间系列的波动模式 （周期变化模式） 以及这些模式随

6、时间的变化 （Furon et al., 2008; Jevrejeva et al., 2003） 。 小波变换可以分为连续小波变换（the Continuous Wavelet Transform, CWT）与离散小波 变换（Discrete Wavelet Transform, DWT） 。离散小波变换 DWT 是数据的紧凑表示，常用于 降噪与数据压缩。 连续小波变换 CWT 更适合于信号特征的提取 （Grinsted et al., 2004） 。 CWT 作为时间系列间歇式波动特征提取的工具被广泛的应用的地球物理学研究中 （Grinsted et al., 2004; Furon e

7、t al., 2008） 。 （1）连续小波变换）连续小波变换 CWT 可以将具有等时间步长 t 的离散时间系列 xn (n=1, N)的连续小波变换定义为小波函 数 0尺度化以及转换下的 xn的卷积： 1 0 * )( )( N n n X n s tnn x s t sW （1） 式中*表示共轭复数，N 是时间系列的总数据个数，(t/s)1/2是一个用于小波函数标准化的因 子从而使得小波函数在每个小波尺度 s 上具有单位能量。 通过转换小波尺度 s 并沿着时间指 数 n 进行局部化， 最终可得到一幅展示时间系列在某一尺度上波动特征及其随时间变化的图 谱，即小波功率谱（Torrence an

8、d Compo, 1998; Torrence and Webster, 1999; Grinsted et al., 2004） 。 对一个时间系列进行小波转换时，母小波的选择显得尤为重要，Farge（1992）曾经讨 论过母小波选择时需要考虑的因素，例如正交与非正交、负值与实值、母小波的宽度与图形 等等。 正交小波函数一般用于离散小波变换， 非正交小波函数即可用于离散小波变换也可用 于连续小波变换（Torrence and Compo, 1998） 。通常在对时间系列进行分析时，希望能够得 到平滑连续的小波振幅，因此非正交小波函数较为合适。此外，要得到时间系列振幅和相位 两方面的信息，就要

9、选择复值小波，因为复值小波具有虚部，可以对相位进行很好的表达 （Torrence and Compo, 1998） 。Morlet 小波不但具有非正交性而且还是由 Gaussian 调节的指 数复值小波。 2/4/1 0 2 0 )( tti eet （2） 式中 t 为时间，0是无量纲频率。当 0=6，小波尺度 s 与傅里叶周期（period）基本相等（, = 1.03s） （Torrence and Webster, 1999） ，所以尺度项与周期项可以相互替代。由此可见， Morlet 小波在时间与频率的局部化之间有着很好的平衡 （Grinsted et al., 2004） 。 此外，

10、 Morlet 小波中还包含着更多的振动信息， 小波功率可以将正、 负峰值包含在一个宽峰之中 （Torrence and Compo, 1998） 。 （2）小波功率谱）小波功率谱 为使计算更为快捷，公式（1）的卷积在傅里叶域内执行（Torrence and Compo, 1998; Grinsted et al., 2004） 。 2 )(sW X n 定义为小波功率谱(wavelet power spectrum)，该功率谱表达 了时间系列在给定小波尺度和时间域内的波动量级（Lafrenire and Sharp, 2003） 。由于我们 采用的 Morlet 母小波为复值小波，因此)(s

11、W X x 也为复数，其复值部分可以解释为局部相位 （Torrence and Compo, 1998） 。将小波功率谱在某一周期上进行时间平均，我们可以得到小 波全谱（global wavelet spectrum） ， 1 2 2 )( 1 )( N on n sW N sW （3） 小波全谱能够表明时间系列真实功率谱的无偏、一致估计（Torrence and Compo, 1998） 。 由于小波全谱可以显示出背景谱量度，所以局部小波谱的峰值可以得到验证。因为该特性， 通过小波全谱图中可以清晰的辨别时间系列的周期波动特征及其强度。 （3）小波功率谱边缘效应及影响锥）小波功率谱边缘效应及影

12、响锥 由于小波变换假设数据是循环的， 所以当我们处理有限长度的时间系列时， 在小波功率 谱中会出现边缘效应，即在功率谱的起始及末端部分出现误差。由于该原因，需要我们在时 间系列的末尾补零从而使得分析的时间系列的总长度 N 大于 2m而小于 2m+1。但是，当我们 采取这样的措施时会在小波功率图谱边缘引起端点不连续以及谱振幅下降的现象。 在这种情况下，需要明确一个概念，即影响锥（Cone of Influence, COI） ，影响锥 COI 表示小波谱区域以及相应的边缘效应。 在COI的边缘小波谱值会下降e-2（Torrence and Compo, 1998; Grinsted et al.

13、, 2004; Furon et al., 2008） 。 （4）小波功率谱的显著性检验）小波功率谱的显著性检验 小波功率谱的统计显著性可以对照一个原假设进行评价， 该原假设为假设信号由一个给 定背景功率谱 （Pk） 的稳定过程产生， 通常背景功率谱为白噪声或红噪声 （Torrence and Compo, 1998; Lafrenire and Sharp, 2003） 。 由于许多地球物理时间系列具有红噪声特征 （即方差随着 尺度的增加或频率的下降而增加） ，所以常采用红噪声作为背景谱对小波谱进行检验。红噪 声过程可以很好的由一阶自回归过程（AR1）来模拟（Torrence and Com

14、po, 1998; Grinsted et al., 2004） 。一个由 lag-1 自相关 处理的 AR1 的傅里叶功率谱可以定义为： 2 2 2 1 1 ki k e P （4） 式中 k 为傅里叶频率指数。 通常在研究中，每个尺度上用 COI 以外的值以 5%的显著水平进行估计。 （5）尺度选择）尺度选择 在进行小波变换时，还需要选择一系列尺度 s。本研究使用非正交小波变换，我们可以 使用任意一组的尺度来构建较完整的图像。一系列尺度可以用 2 的分数幂来表达： jj j ss 2 0 ，j = 0, 1, J （5） )/(log 02 1 stNjJ （6） 式中，s0为可分辨的最小

15、尺度，J 为确定的最大尺度。s0应该被选择恰当以便使相等的 傅里叶周期近似于 2t。 一个足够小的 j 的选择依赖于小波方程谱空间的宽度。 4. 小波分析步骤小波分析步骤 选用开源小波分析软件（http:/paos.colorado.edu/research/wavelets/），以 47 年年降水数 据为例，操作步骤如下（本实验中使用的数据和部分小波分析代码见附录 A）： （1）数据预处理数据预处理。进行小波分析的年降水量时序数据必须是连续等时间步长的。另外， 需要进行标准化处理。标准化处理可选择 SPSS 软件，对降水量时序数据进行 Z 值标准化， 使其成为无量纲时序数据； （2）母小波选

16、取母小波选取。由于小波分析能够很好的反应序列的尺度变化、周期变化和趋势变 化，将其应用在降水序列的分析中，不仅能反应降水序列的尺度变化，还可以显示出降水序 列的周期变化以及变化的时间位置。可选择 Mexican hat 小波，或 Morlet 小波。通常在对时 间系列进行分析时，希望能够得到平滑连续的小波振幅，因此非正交小波函数较为合适。本 实验中选择 Morlet 小波，它比 Fourier 分析更能反映出序列的局部特征。 （3）尺度选择尺度选择。由于年降水时间序列为 47 年，时间系列长度 N=47，为了减小功率谱 的边缘效应，在进行交互小波变换时选择 26个数据。时间步长 dt=1，即一年一个数据。j 可选择 0.125。 （4）显著性检验显著性检验。由于许多地球物