李梦楚,胡春宏
(中国水利水电科学研究院 流域水循环模拟与调控国家重点实验室,北京 100038)
小波分析是20世纪80年代初发展起来的一种信号时、频局部化分析新方法。目前在信号处理、图像压缩、语音编码、模式识别、地震勘探、大气科学以及许多非线性科学领域内取得了大量的研究成果。随着小波理论的形成和发展,它的优势逐渐引起许多水科学工作者的重视并且引入水文学科中[1-2]。水文时间序列是一个观测到的样本,是随时间变化的时变序列,几乎都是非平稳的。对于这一类非平稳序列,需要提取某一时段(或极小瞬间)的频域信息,或某一频段所对应的时间信息。而小波变换是一种具有时频多分辨率的分析方法,能把时间序列在时间和频率域上同时展开,进行多尺度细化分析,研究它们在不同周期随时间变化的特性[1]。
从1993年Kumar和Foufoular-Gegious[3]将小波变换介绍到水文学中,到目前已取得了一定的研究成果,主要表现在水文系统多时间尺度分析、水文时间序列变化特性分析、水文系统预测预报和水文系统随机模拟等方面。关于黄河水沙变化的小波分析多是关于径流量的分析[4-5],也有一些关于上游输沙量的分析[6],但对于黄河下游研究较少,也没有人把小波分析引进到河床演变研究中来。
本文选取黄河下游小浪底、花园口、高村和利津四站,利用小波理论研究此4站的水量、沙量序列多时间尺度变化,并以花园口、高村、利津三站平滩流量序列作为河道断面形态特征参数与水沙序列变化进行对比,对于黄河下游水沙变化及水沙变化对断面形态的影响研究都有重要意义。
1.1 小波函数 小波函数是具有震荡特性、能够迅速衰减到零的一类函数,即也称基小波,其伸缩和平移构成一簇函数系:
式中:Ψa,b(t)称子小波;t为时间,a为尺度因子或频率因子,反映了小波的周期长度;b为时间因子或平移因子,反映了在时间上的平移。
1.2 小波变换 对于能量有限信号f(t)∈L2(R),其连续小波变换定义为
当研究离散信号f(iΔt)(i=1,2,…,N;N为样本容量;Δt为取样时间间隔)时,式(2)的离散形式表达为
Wf(a,b)能同时反映时域参数b和频域参数a的特性,它是时间序列f(t)或f(iΔt)通过单位脉冲响应的滤波器的输出。当a较小时,对频域的分辨率低,对时域的分辨率高;当a增大时,对频域的分辨率高,对时域的分辨率低。因此,小波变换实现了窗口的大小固定、形状可变的时域局部化。
Wf(a,b)随参数a和b变化,可作出以b为横坐标,a为纵坐标的关于Wf(a,b)的二维等值线图,称为小波变换系数图。通过小波变换系数图可得到关于时间序列变化的小波变化特征。在尺度a相同情况下,小波变换系数随时间的变化过程反映了系统在该尺度下的变化特征;正的小波变换系数对应于偏多期,负的小波变换系数对应于偏少期,小波变换系数为零对应于突变点;小波变换系数绝对值越大,表明该时间尺度变化越显著。
1.3 小波方差 将时间域上的关于a的所有小波变换系数的平方进行积分,即为小波方差
小波方差随尺度a的变化过程称小波方差图。它反映了波动的能量随尺度的分布。通过小波方差图,可以确定一个水文序列中存在的主要时间尺度,即主周期。
水沙序列资料采用小浪底、花园口、高村、利津四站1952—2009年年水量、年沙量、汛期、非汛期水量、沙量资料。小浪底站年水沙量过程如图1。
采用Morlet小波为母函数。Morlet小波为复数小波,定义为φ(t)=eiω0te-t22,式中ω0为常数且当ω0≥5时Morlet小波就能近似满足允许性条件,其小波变换系数的模和实部是两个重要的变量。模的大小表示特征时间尺度信号的强弱,实部表示不同特征时间尺度信号在不同时间上的分布和位相两方面的信息。将各序列和Morlet小波函数代入式(3),取不同的a和b,计算小波变换系数的模平方和实部,作它们的二维等值线图,并求出小波方差,作图。Morlet小波变换系数的模值表示能量密度,
由于运用小波分析时间序列需较长资料序列,因此采用以端点为中心对称延拓的方法将两端数据外延,消除边界效应。同时为应用方便将数据进行距平处理。小波变换系数的模值越大,表明其所对应的时段和尺度的周期性越明显,模平方图即能量图;小波系数的实部包含着给定时间和尺度信号相对于其他时间和尺度信号的强度和位相两方面的信息;小波方差的峰值则对应主要周期,小波方差越大即表明在该尺度下周期变化明显。
图2(a)、(b)分别给出了小浪底站年水量距平序列Morlet小波变换系数的模平方和实部的时频分布,(a)中虚线表示负位相,对应偏少期,实线表示正位相,对应偏多期,交替处即为突变点。由图可以看出,18~30a时间尺度十分突出,周期中心在23a,在1952—2009年都表现明显,正负位相交替出现。12~18a时间尺度也比较强,周期中心在14a,在1952—2009年都表现明显,正负位相交替出现。2~6a时间尺度在1952—1980年表现比较明显,周期中心主要在3a左右,正负位相交替出现。6~11a时间尺度在1968—1995年表现比较明显,其周期中心在9a左右,正负位相交替出现。
由图3(a)、(b)可以看出年沙量周期变化、强弱、突变点分布及其位相结构。8~16a时间尺度表现明显在1952—1990年左右,周期中心在11a,1985—2009年在此时间尺度附近体现为11~19a及5~10a时间尺度,周期中心分别为14a和8a。17~30a时间尺度在1952—2009年全时段都体现比较明显,周期中心在23a。3~7a周期在1952—1995年表现比较明显,周期中心主要在5a左右。30a以上周期也有所体现,但由于采用的时间序列较短,不很明显。
用式(4)计算小波方差,并作出小浪底站年水沙量距平序列的小波方差图如图4、图5。由图可清楚看出,小浪底年水量具有3、9、14、23a的变化周期,其中23a是主要周期;小浪底年沙量具有5、11、23、36a的变化周期,其中11a是主要周期。1a自然周期能量也较大,但忽略之。
小浪底站年水量与年沙量序列总体来说变化规律是一致的,都体现出了4a左右、10a左右、23a左右的周期变化。同样方法分析小浪底站汛期、非汛期水量、沙量以及花园口、高村、利津三站的年水沙量、汛期、非汛期水沙量,可得类似结果,将四站水沙变化周期列于表1(表中黑体表示为此序列主周期)。可知各站相同序列变化基本一致,年水量与汛期、非汛期水量变化周期也基本一致,非汛期沙量与年沙量、汛期沙量的周期稍有差异。水沙总体变化周期都是一致的。
中国大陆华北地区年降水存在2、4、7、15、26a左右的主周期,15a左右为第一主周期,西北地区存在2、4、9、20a左右的主周期,20a左右为第一主周期[11]。而太阳黑子存在9a、41a的周期变化特征[10]。故可知黄河下游水沙变化周期与年降水变化相关,可能也与太阳黑子活动有关。
表1 四站水沙变化周期 (单位:a)
采用花园口、高村、利津三站1950—2008年平滩流量资料(来源于黄委会)。花园口站平滩流量序列如图6。
同样的,对平滩流量序列也进行延拓、距平,进行如上分析。以花园口为例,分析如下:
图7(a)、(b)分别给出了花园口站平滩流量距平序列mor⁃let小波变换系数的模平方和实部的时频分布,由图可以看出花园口站平滩流量周期变化、强弱、突变点分布及其位相结构。18~30a时间尺度十分突出,周期中心在24a,在1950—2008年都表现明显,正负位相交替出现。12~20a时间尺度在1950—1975年表现明显,周期中心在16a,正负位相交替出现。7~12a时间尺度在1975—2009年表现比较明显,周期中心主要在10a左右,正负位相交替出现。5~10a时间尺度在1950—1975年表现比较明显,其周期中心在7a左右,正负位相交替出现。作小波方差图如图8。可看出24a时间尺度表现非常明显,为花园口站平滩流量主要周期,其他周期表现较弱。
高村、利津两站同理分析,将平滩流量变化周期列于表2(表中黑体表示为此序列主周期)。
由此可见,平滩流量的变化与水沙变化周期基本一致,平滩流量变化受到水沙变化周期影响。
为进一步分析,选取小尺度3a、中尺度10a、大尺度23a的各序列小波变换系数实部变化作图进行分析(以小浪底年水沙量、花园口平滩流量为例)。
图9给出了小浪底站年水沙量及花园口站平滩流量3a的小波变换实部变化过程,正相位对应偏多时期,负相位对应偏少时期,小尺度下变异点不具有代表性,有许多伪突变点,应在中尺度和大尺度下判断,故变化点不一一列出。
表2 花园口、高村、利津平滩流量变化周期(单位:a)
由图9可以看出,在20世纪70年代中期以前年水沙量序列与平滩流量序列变化基本一致;70年代中期到80年代中期这一时段内,3个序列变化不同步,影响因素比较复杂,由于多时间尺度是从整体序列分析出来的,而同一尺度下的丰枯变化又是相对于本身序列的数值变化得出的[6],初步分析认为可能与1974年、1986年三门峡水库调水调沙及1986年龙羊峡水库建成运用有关;80年代后期以后3个序列变化逐渐趋于一致,90年代后期到2003年左右平滩流量序列与水沙序列变化稍有紊乱后归于一致,可能与小浪底水库建成运用有关。
图10给出了小浪底站年水沙量及花园口站平滩流量10a的小波变换实部变化过程,可以看出:小浪底年水量序列在1953年以前、1958—1963年、1968—1973年、1977—1982年、1986—1991年、1996—2002年、2007年以后各时段为负相位,表示年水量偏少;而1953—1958年、1963—1968年、1973—1977年、1982—1986年、1991—1996年、2002—2007年各时段为正相位,表示年水量偏丰。由此可推断在此尺度下2009年以后两三年内仍然处于流量偏少时段,2012年以后会处于流量偏丰时段。小浪底年沙量序列在1954年以前、1960—1965年、1970—1976年、1981—1986年、1991—1994年、1998—2003年、2007年以后各时段为负相位,表示年沙量偏少;而1954-1960年、1965—1970年、1976—1981年、1986—1991年、1994—1998年、2003—2007年各时段为正相位,表示年沙量偏多。由此可推断在此尺度下2009年以后两三年内仍然处于年沙量偏少时段,2012年会处于年沙量偏多时段。花园口平滩流量序列在1953年以前、1957—1963年、1969—1974年、1979—1984年、1989—1994年、2000—2006年各时段为负相位,表示平滩流量偏小;而1953—1957年、1963—1969年、1974—1979年、1984—1989年、1994—2000年、2006年以后各时段为正相位,表示平滩流量偏丰。由此可推断在此尺度下2009年以后两年内仍然处于平滩流量偏丰时段,2011年以后会处于平滩流量偏丰时段。
在中尺度下,水沙变化及平滩流量变化是一致的。沙量变化及平滩流量变化始终跟随水量变化,符合实际力学规律,且能够体现出河床演变滞后于水沙变化的现象,没有明显的紊乱的时段。
图11给出了小浪底站年水沙量及花园口站平滩流量23a的小波变换实部变化过程,可以看出:小浪底年水量序列在1957年以前、1968—1980年、1991—2003年各时段为负相位,表示年水量偏少;而1957—1968年、1980—1991年、2003年以后各时段为正相位,表示年水量偏丰。由此可推断在此尺度上2009年以后5、6a内仍然处于流量偏丰时段,2015年以后会处于流量偏少时段。小浪底年沙量序列在1957—1969年、1980—1992年、2003年以后各时段为负相位,表示年沙量偏少;而1957年以前、1969—1980年、1992—2003年各时段为正相位,表示年沙量偏多。由此可推断此尺度上2009年以后5、6a内仍然处于年沙量偏少时段,2015年以后会处于年沙量偏多时段。花园口平滩流量序列在1956年以前、1967—1979年、1991—2003年各时段为负相位,表示平滩流量偏小;而1956—1967年、1979—1991年、2003年以后各时段为正相位,表示平滩流量偏丰。由此可推断2009年以后5、6a内仍然处于平滩流量偏丰时段,2015年以后会处于平滩流量偏少时段。
在大尺度下可看到,水沙序列呈相反趋势,水量偏多期对应于沙量偏少期,反之亦然。而平滩流量的变化基本与年水量变化一致,略有所超前,在图9上补充了小浪底站汛期水量的变化,发现平滩流量的变化与汛期水量的变化更为吻合,是跟随汛期水量变化而变化的,所以在长时间尺度上汛期水量对于造床的贡献可能更大。
综上可知,在大、中、小不同时间尺度上,水沙序列、平滩流量序列的小波变换实部变化基本一致,在小尺度和大尺度上体现出的不同步可能有两个主要原因:一是黄河流域水沙异源的特性,二是人类活动特别是大型水库建设的影响。
将小浪底年水沙序列及花园口平滩流量序列在中、大尺度上得变异点列于表3。同尺度下各序列变异点基本吻合,比较稳定。序列的奇异点可能在不同尺度上是不一致的,与时间尺度有关。一般在较大尺度上寻找奇异点,从而判断奇异性的存在,但也要综合各个尺度进行系统分析。此结果基本吻合1960年三门峡水库建成投运,1968年刘家峡水库蓄水运用,1969年开始的三门峡水库第二期改建,2002年小浪底水库建成完工等大的历史事件,故人类活动在很大程度上影响了水沙序列的变化,下游河道也因此产生很大影响。
表3 小浪底年水沙序列及花园口平滩流量序列在10a、23a时变异点
从1952年到2009年,下游来水来沙经历了3个阶段:1960年以前的天然阶段、1960—1999年三门峡水库单独运用阶段以及1999年以后小浪底水库投入运用。本文选取全时段进行研究,先在多时间尺度层面对整个黄河下游水沙及平滩流量变化周期进行了整体规律的研究,发现变化规律比较明显,又在不同时间尺度上进行了分析,从宏观到微观,总结了在不同自然条件及人类活动影响下黄河下游的变化情况。通过对黄河下游小浪底、花园口、高村、利津4站水沙序列及花高利三站平滩流量序列进行小波分析,发现了在多时间尺度下水沙与平滩流量的变换特征,结果表明:(1)黄河下游水量主要存在3a、9a、14a、23a变化周期,其中23a为最显著周期;沙量主要存在5a、11a、23a、36a变化周期,其中11a为最显著周期;平滩流量主要存在24a变化周期。水沙序列存在相似的多时间尺度变化特征,平滩流量序列受到水沙序列显著周期的影响。水沙变化周期与华北、西北地区降水变化周期基本吻合,也可能与太阳黑子变化有关。(2)在小尺度(3a)下,年水量、年沙量、平滩流量曲线变化基本一致,在某些时段发生紊乱;在中尺度(10a)下,年水量、年沙量、平滩流量曲线变化有一定相关性及同步性,相位基本一致;而大尺度(23a)下,年水量、平滩流量曲线变化基本一致,平滩流量与汛期流量变化相关性更强,与年沙量序列相位相反。各尺度下水量变化、沙量变化与平滩流量变化有所不一致可能是由于受到人类活动及自然情况的影响,也可能是自然规律的反映,需进一步研究。(3)各尺度下各序列的突变点基本一致,随着时间尺度的增大,小时间尺度上的一些突变点逐渐退化为平常点,相反,随着时间尺度的降低,除原有的突变点外,还不断增加突变点,大时间尺度上的一些平常态,在低层次看来则是突变态。因此,突变点与时间尺度有密切关系。黄河下游水沙及平滩流量序列的突变点主要为1957、1968、1980、1991、2003年左右。
[1]王文圣,丁晶,李跃清.水文小波分析[M].北京:化学工业出版社,2005.
[2]王文圣,丁晶,向红莲.水文时间序列多时间尺度分析的小波变换法[J].四川大学学报:工程科学版,2002,34(6):14-17 .
[3]P-Kumar,E-Foufoula-Georgiou.A multi-component decomposition of spatial rainfall fields 1.Segregation of large and small scale features using wavelet transforms[J].Water Resources Research.1993,29(8):2515-2532.
[4]蒋晓辉,刘昌明,黄强.黄河上中游天然径流多时间尺度变化及动因分析[J].自然资源学报,2003,18(2):142-147.
[5]刘俊萍,田峰巍,黄强.黄河上游河川径流变化多时间尺度分析[J].应用科学学报,2003,21(2):117-121.
[6]王秀杰,练继建.近43年黄河上游来水来沙变化特点[J].干旱区研究,2008,25(3):342-347.
[7]姜世中,梁川.龙川江年输沙率时间序列的小波特征[J].哈尔滨工业大学学报,2009,41(11):197-200.
[8]高鹏,穆兴民,王炜.基于Morlet小波的太阳活动与黄河输沙量关系分析[J].水土保持研究,2010,17(3):181-183 .