基于小波消噪的无定河径流趋势变化分析

基于小波消噪的无定河径流趋势变化分析——以赵石窑站为例

王化齐^[1]1，2，张茂省1，党学亚1

（1中国地质调查局西安地质调查中心，西安710054；2长安大学环境科学与工程学院，西安710054）

【摘　要】根据无定河赵石窑站1954～2001年径流实测资料，采用Morlet函数作为小波函数，对序列结构进行详细剖析，探讨径流序列多时间尺度变化以及丰枯突变特征。结果表明，径流时间序列存在近似334个月的主周期，141个月、84个月左右的次周期；进行小波消噪后的径流序列的信息数减少，降低了径流序列的复杂性，使其尽可能的展现径流序列的内在特征；在不同的尺度下，对应小波系数的波形都相似，随着尺度的细化，突变点有所增加，说明在更小的尺度上更具体的反应出径流序列的局部变化特征。

【关键词】小波分析法，径流量，径流时间序列

The tendency analysis on runoff of Wuding River based on waveletdenoising

——A case study in Zhaoshiyao station

Wang Huaqi^1，2　Zhang Maosheng¹　　Dang Xueya¹

（1.Xi’an Center of Geological Survey，CGS，Xi’an，710054；

2.School of Environment Science＆Engineering，Chang’an University，Xi’an，710054）

【Abstract】According to the runoff data（1954～2001）at Zhaoshiyao station of Wuding River，the series structure was analyzed by means of the Morlet function as wavelet function，and the multi scale and oscillation characteristics of the time series were discussed.The results indicate that，runoff series exist multi－scale period characteristic，namely，the main period is 334 months and the secondary period is141 and 84 months；the information deceased and became simple by denoising，thus it can represent the inner characteristic of the runoff series；for different time scale，the wave form of wavelets coefficient were similar，the saltation points increased with scale refinement，which means the local variation characteristics of runoff series couldbe identify in smaller scale.

【Key words】wavelet analysismethodrunoffrunoff time series

1.引　言

由于实测水文序列受到多因素影响，含有系统噪声和测量噪声，因而采用含有噪声的数据进行分析计算就不能真实反映水文系统的本质，故需要做消噪处理。^［1］传统的维纳滤波、卡尔曼滤波等消噪技术仅能适用于线性系统且严格依赖于状态空间函数的建立，而水文系统是非线性的且难于建立合适的状态空间模型，极大地限制了其应用范围。对径流序列多尺度分析的传统方法有滑动平均、滤波、Fourier分析等。但上述方法或多或少的具有自身的缺陷:在时域和频域上不具有局部化性质；对突变点的诊断缺乏数学上的严谨性。20世纪80年代中期发展起来的小波分析法（wavelet Analysis Method），亦称为多分辨率分析（Multi－resolution analysis），是近些年国际上十分热门的前沿领域，被认为是傅里叶分析方法的突破性进展。^［2］小波分析是一种多尺度分析方法，能有效分离高频成分与低频成分。根据不同信号在小波变换后表现的不同特性对小波分解序列进行处理，将处理后的序列加以重构即可实现信噪分离，^［3］因此采用小波分析消噪效果较为理想。

无定河是黄河中游较大的一级支流，发源于陕西省定边、靖边和吴起县交界处的白于山，穿越陕西、内蒙两（区）13个旗县，于清涧县河口附近入黄河，干流全长491km，流域面积约3026×10⁴ km²，^［4］年径流量15.3×10⁸ m³，夏秋两季径流量占全年的65%以上。北岸是毛乌素沙漠，南岸是黄土沟壑区，水土流失严重，河水含沙量大，多年平均每立方米浑水中含泥沙144kg，平均每年输入黄河的泥沙达2.23×10⁸ t。

本文以榆林地区无定河赵石窑站1954～2001年实测径流资料为依据，通过小波变换将径流序列分解成高频、低频序列，采用stein无偏阈值方法对分解后的序列进行消噪；通过对消噪前后径流序列的小波变换系数所含信息数大小的分析，确定消噪效果如何；对消噪后的径流序列进行连续小波变换，通过小波变换系数的实部、模平方、方差的分析确定径流序列所包含的主要尺度，以及这些尺度下径流序列的丰枯变化。

2.小波分析方法

小波分析是一种窗口大小固定但形状可变（时宽和频宽可变）的时频局部变化分析方法，它具有自适应的时频窗口:高频段时，频域窗口增大，时间窗口减小；低频段时，时间窗口增大，而频域窗口减小。小波分析的关键在于引入了满足一定条件的基本小波函数ψ（t）以代替Fourier变换中的基函数e^－iωt。ψ（t）经伸缩和平移得到一族函数:

式中:ψ_a，b称分析小波或连续小波；a为尺度（伸缩）因子，在一定意义上1/a对应于频率ω；b为时间（平移）因子，反应时间上的平移。

因此，小波分析实际上是一族小波函数表示的时间函数f（t），通过小波变换可得到小波变换系数，对这些系数的分析可以显示出f（t）的时频特性和其重要的局部变化特性。当a减小时，ψ_a，b（t）的时域波形在时间轴方向上收缩，分析信号细节，得到信号的高频信息；当a增大时，ψ_a，b（t）的时域波形在时间轴方向上伸展，分析信号的概貌，获得信号的低频信息。也就是说，通过调整a的大小，改变时频窗口的时宽和频宽，从而实现了信号时频局部不同分辨率的分析。

式中，W_f（a，b）既包含了f（t）的信息，又包含了ψ_j，k（t）的信息。因此，小波函数的选择十分重要。目前广泛使用的有Haar小波、墨西哥帽（Marr）小波、Morlet小波、Daubechies小波等。本文采用Morlet小波函数对榆林地区的降水、气温时间序列进行连续小波变换，Morlet小波是高斯包络下的单频率复正弦函数，其形式如下:

式中，ω₀为常数，i表示虚数，Morlet小波的时、频域局部性都比较好。

3.数据分析及处理

3.1　日径流序列的小波消噪分析

对1954～2001年的无定河赵石窑站径流量资料进行分析。首先对径流序列进行分辨率分解，选择小波函数为db4，经分析确定最大分解层数为4层。经过四层分解后的径流序列图如图1所示。

由图中可以看出，小波分析可以清晰的将径流序列分为低频序列和高频序列。选择Stein（史坦）的无偏风险阈值方法和软阈值处理方法对原始序列进行消噪，消噪效果图如图2所示。

图1　四层小波分解后的径流序列

Fig.1　Runoff series of four layer wavelet decomposition

图2　小波消噪结果图

Fig.2　Map of Wavelet denoising results

将低频系数和消噪后的高频系数进行小波重构得到消噪后的径流序列。消噪后的径流序列如图3所示。由图1可见，消噪后的年径流序列较平滑，趋势更明显，深入分析不同时间段的径流变化趋势，发现从大时间尺度上，整体上呈现减小趋势。

图3　消噪后的径流序列

Fig.3　Denoised runoff series

3.2　径流序列的小波信息数

径流序列消噪后，无法直接确定其消噪效果的好坏。计算径流的小波信息数可以确定径流序列的复杂程度，因此可以通过比较消噪前后径流序列的小波信息数的大小来确定消噪的效果。

对径流序列f（t）进行离散小波变化，其中尺度因子a和平移因子b的离散化公式为:a＝2^j，b_jk＝2^j，k，j，kZ。a_j（j＝1，2，…，N；N为最大分解层数）尺度下的离散小波变换系数d_j，k可由Mallat算法求出:

利用小波变换系数d_j，k可求得径流序列f（t）在各种尺度j下的能量E_j:

则径流序列f（t）的总能量为:

则径流序列f（t）各个尺度上的概率分布P _j为:

尺度a＝2^j下的小波变换系数d_j，k对应的频带ΔF为:

将P_j的Shanon熵定义为信息量系数（information cost function，ICF），即:

信息量系数ICF计算出了系统能量分布的相应信息，表征了被研究系统的复杂程度。径流序列f（t）越无规则，越复杂，ICF越大；反之，能量分布集中于某一频带，时间序列性越强，ICF越小。

根据基于小波分析的信息数的计算办法计算消噪前后序列的信息数ICF，结果如下表所示:

表1　消噪前后径流序列的信息数

Table 1　The information number of runoff series before and after denoising

由上表可以看出消噪后的序列的信息数减少，小波消噪方法可以降低径流序列的复杂性，使其尽可能的展现径流序列的内在特征。

3.3　日流量序列的多尺度分析

小波分析在时域和频域上都有很好的局部化功能，因而可以用来分析径流序列的局部特征。所谓多时间尺度是指系统变化在时域中存在多层次时间尺度结构和局部化特征。多时间尺度分析可以为系统提供不同时间尺度下的演变规律和发展趋势，从而能够为系统中长期预测提供背景分析依据。

选择Morlet小波函数对消噪后的径流序列进行连续小波变换，计算小波变换系数的模平方、实部和方差，并绘制它们随尺度变化的分布图，如图4、图5所示。

图4　小波变换实部图

Fig.4　Real Wavelet transforms graph

图5　小波变换模平方图

Fig.5　Wavelet transform modulus square graph

图6　小波变化方差图

Fig.6　The variance diagram of wavelet transforms

根据图6可以看出，赵石窑站径流序列存在334个月、141个月以及84个月的主周期。分析主要尺度下的小波变换实部变化过程，可以找到各种尺度下的径流序列的丰枯变化过程。

分析主要尺度下的小波变换实部变化过程，可以找到各种尺度下的径流序列的丰枯变化过程。图7所示，在334个月尺度，丰枯的突变点在1970年（丰到枯）、1986年（枯到丰），从未来一段时间内径流是偏枯的。

图7　在334个月尺度下小波变换实部变化过程图

Fig.7　The real part of the change process of wavelet transform in 334 months

从图8可以看出，在141个月尺度，丰枯的突变点在1960年（丰到枯）、1967年（枯到丰）、1973年（丰到枯）、1980年（枯到丰）、1987年（丰到枯）、1994年（枯到丰）、2001年（丰到枯）。图9为在84个月尺度下小波变换系数的实部变化过程，也反映了径流序列的丰枯变化过程。随着尺度的细化，突变点越多。

图8　在141个月尺度下小波变换实部变化过程图

Fig.8　The real part of the change process of wavelet transform in 141 months

图9　在84个月尺度下小波变换实部变化过程图

Fig.9　The real part of the change process of wavelet transform in 84 months

4.结　论

本文利用Morlet小波分析方法对榆林地区降水和蒸发进行了时间变化特征分析，得到如下结论:

（1）通过无定河赵石窑水文站1954～2001年间径流序列的Morlet小波分析，发现径流序列存在多时间尺度变化特征，径流时间序列存在主周期为334个月，次周期为141个月、84个月左右等。

（2）对径流量的趋势分析，消噪系列与实测序列的变化趋势吻合，但在变化趋势的显著性分析上存在差异，消噪后的序列的信息数减少，小波消噪方法可以降低径流序列的复杂性，使其尽可能的展现径流序列的内在特征。采用Stein阈值确定法进行消噪，且应对消噪后的径流序列进行合理性分析，适当进行人工修正，避免采用纯数学方法而使分析结果失真。

（3）通过对不同尺度的径流时间序列的小波变换过程进行分析，在不同的尺度下，对应小波系数的波形都相似，随着尺度的细化，突变点有所增加，说明在更小的尺度上更具体的反应出径流序列的局部变化特征，即在短时期内径流量的变化存在波动性，但在较长时期内还是具有整体变化趋势。

（4）小波分析的时频局部化特性可展现径流时间序列的精细结构，不仅可以将隐含在序列中的随时间变化的周期震荡显现出来，还能确定出径流丰枯突变点的位置，也可以对时间序列的演变趋势进行定性的估计，为分析多时间尺度变化特征及短期水文预报提供了一条新途径。

（中译英:张艳，王金凤检）

参考文献（Reference）

［1］刘冀，董晓华，李帅.基于小波消噪的径流趋势变化对比分析.［J］.水电能源科学，2009.27（5）:4－7.

［2］王化齐，张茂省，党学亚.榆林地区降水蒸发时间序列的多尺度特征和突变分析.［J］.水电能源科学，201028（9）:1－4.

［3］郭慧芳.鄱阳湖水资源优化配置的智能算法及其应用研究.2010，河海大学:南京.

［4］杨新，延军平，刘宝元.无定河年径流量变化特征及人为驱动力分析.［J］.地球科学进展，2005.20（6）:637－642.

【注释】

[1]作者简介:王化齐，中国地质调查局西安地质调查中心，博士后。张茂省，中国地质调查局西安地质调查中心研究员、博导。