邵雪杰，顧圣平，曹愛武，何海祥，于婷婷

(河海大學(xué) 水利水電學(xué)院，南京 210098)

多沙河流水沙變化特征小波分析

邵雪杰，顧圣平，曹愛武，何海祥，于婷婷

(河海大學(xué) 水利水電學(xué)院，南京 210098)

為了給多沙水庫(kù)運(yùn)行調(diào)度提供依據(jù)，降低水沙序列長(zhǎng)度和隨機(jī)性的影響，將小波變換應(yīng)用于研究流域徑流與含沙量的變化特性。采用db3小波對(duì)標(biāo)準(zhǔn)化徑流與含沙量序列進(jìn)行多分辨率分析，研究徑流與含沙量變化的趨勢(shì)性；選用復(fù)Morlet小波繪制出小波方差圖，分析徑流過程與含沙量過程存在的周期性。以崖羊山水電站所在的李仙江流域?yàn)槔?，針?duì)電站壩址斷面的月平均流量與含沙量序列進(jìn)行小波分析，從低頻重構(gòu)序列的結(jié)果中可以看出該流域徑流與含沙量的變化趨勢(shì)；并結(jié)合當(dāng)?shù)氐慕涤炅颗c水土保持狀況分析，表明結(jié)果是合理的。根據(jù)小波方差圖可以看出，崖羊山水電站所在流域徑流過程與含沙量過程存在非常接近的顯著周期，均為2 a左右，且兩者變化周期具有同步性。研究結(jié)果表明小波分析是研究非平穩(wěn)隨機(jī)時(shí)間序列的有效方法。

小波變換；水沙序列；徑流；含沙量；時(shí)間尺度；多分辨率分析；小波方差

水沙序列一般是非平穩(wěn)的離散等間距序列，具有隨機(jī)性與周期性[1]。為了給多沙水庫(kù)運(yùn)行調(diào)度提供依據(jù)，需要掌握水庫(kù)所在流域水沙序列的變化規(guī)律。以往的水文學(xué)中常常用譜分析法和相關(guān)分析法來分析水沙序列的變化特性，前者受序列長(zhǎng)度影響較大，后者受數(shù)據(jù)的隨機(jī)性影響較大。由于具備較好的時(shí)頻局部化的特點(diǎn)，小波分析方法非常適合用于顯示水沙時(shí)間序列的細(xì)微結(jié)構(gòu)[2]，其核心為小波變換。小波變換通過計(jì)算數(shù)據(jù)中能反映變化特性的小波變換系數(shù)，降低序列長(zhǎng)度和隨機(jī)性的影響，能夠凸顯序列的變化規(guī)律。

1 小波分析原理

(1)

式中：ψa,b(t)為子小波；a為尺度伸縮因子，又稱頻域參數(shù)；b為時(shí)間平移因子，又稱時(shí)域參數(shù)。

(2)

當(dāng)研究離散函數(shù)f(iΔt)(i=1,2,…,N；N為樣本容量；Δt為取樣時(shí)間間隔)時(shí)，式(2)的離散形式表達(dá)為

(3)

當(dāng)基本小波滿足容許性條件時(shí)，即

(4)

(5)

2 小波理論在水沙序列分析中的應(yīng)用

水沙序列受季節(jié)變化的影響。為降低這類影響，在小波理論應(yīng)用到水沙序列分析之前，一般對(duì)原始序列進(jìn)行標(biāo)準(zhǔn)化處理[4]，即

(6)

2.1 水沙序列趨勢(shì)性分析

水沙序列過程一般認(rèn)為由3部分構(gòu)成，即趨勢(shì)項(xiàng)、隨機(jī)項(xiàng)以及周期項(xiàng)。其中趨勢(shì)項(xiàng)對(duì)應(yīng)序列的低頻成分，隨機(jī)項(xiàng)對(duì)應(yīng)序列的高頻成分。

應(yīng)用小波理論進(jìn)行水沙序列分析時(shí)，首先要選擇適當(dāng)?shù)男〔ê瘮?shù)，通常有2種方法：一種是通過多次試算來挑選誤差較小的小波函數(shù)；另一種方法是參照研究對(duì)象的分布形態(tài)特征[5]，依靠經(jīng)驗(yàn)選擇與其序列形態(tài)相近的小波函數(shù)，本文使用的是這種方法。dbN小波是一種具有N階消失矩的緊支撐正交小波[6]。當(dāng)N取3時(shí)，即db3小波，其尺度函數(shù)和小波函數(shù)存在明顯的波峰和波谷[7]，能較好地反映標(biāo)準(zhǔn)徑流與含沙量序列的豐枯變化，選用其作為基本小波進(jìn)行水沙序列趨勢(shì)性分析比較合適。db3小波的小波基(尺度函數(shù)和小波函數(shù))沒有明確的解析表達(dá)式，可以通過給定的濾波器生成。

在實(shí)際研究過程中，水沙時(shí)間序列大多是離散的，運(yùn)用相應(yīng)的離散小波變換對(duì)水沙序列進(jìn)行趨勢(shì)性分析，通常從多分辨率分析(MRA)的方法入手，其基本原理為：對(duì)?f(t)∈L2(R)，都可以分解為

(7)

設(shè)f1,f2∈L2(R)，它們的內(nèi)積定義為

(8)

低頻成分在更低的分辨率下還可以進(jìn)一步分解成高頻成分與低頻成分，分解公式為

(9)

其中：

圖1 多分辨率分析過程Fig.1 Process of Multiple Resolution Analysis (MRA)

相應(yīng)地，分辨率為M的低頻成分重構(gòu)公式為

(10)

多分辨率分析逐步地去除高頻成分，使得低頻成分得以保留，水沙序列的總體趨勢(shì)性越來越明顯。

2.2 水沙序列周期性分析

在周期分析方面，本文選用復(fù)Morlet小波[8]作為水沙序列周期分析的基本小波，其小波函數(shù)為

(11)

式中:ω0為常數(shù)；i為虛數(shù)單位。與實(shí)數(shù)小波相比，復(fù)數(shù)小波更能準(zhǔn)確地反映水沙序列各時(shí)間尺度的周期大小及其在時(shí)域中的分布。

徑流與含沙量序列經(jīng)過復(fù)Morlet小波變換后可得到小波變換系數(shù)的模及實(shí)部。分析小波變換系數(shù)的實(shí)部，可判斷同一時(shí)間尺度的信號(hào)在不同時(shí)域的強(qiáng)弱，信號(hào)較強(qiáng)的時(shí)間段即為大致的周期范圍。

在上述所得到的周期范圍的基礎(chǔ)上，借助小波方差，可以進(jìn)一步確定序列隱含的顯著周期。方差的大小反映不同時(shí)間尺度自身信號(hào)的強(qiáng)弱，方差越大，信號(hào)越強(qiáng)[9]。當(dāng)小波方差達(dá)到最大值時(shí)，小波函數(shù)所選的尺度與隱含周期匹配度最高。離散時(shí)間序列的小波方差可以根據(jù)方差公式計(jì)算得出，其公式為

(12)

式中W(a,xj)代表頻域參數(shù)為a、時(shí)域參數(shù)為xj時(shí)的小波變換系數(shù)。由式(12)可知，在尺度a處，方差V(a)等于每個(gè)樣本對(duì)應(yīng)的小波變換系數(shù)平方W2(a,xj)的均值。對(duì)于復(fù)小波函數(shù)的W2(a,xj)而言，即為小波變換系數(shù)模的平方。

3 實(shí)例分析

本文以崖羊山水電站所在的云南省思茅地區(qū)的李仙江流域?yàn)槔?，針?duì)電站壩址斷面的月平均流量與徑流含沙量序列，進(jìn)行小波分析，研究水沙序列的趨勢(shì)性和周期性特征。

3.1 數(shù)據(jù)標(biāo)準(zhǔn)化

根據(jù)電站壩址處的月平均流量與對(duì)應(yīng)的平均含沙量序列(1976—1995年，以月為單位，各240個(gè)樣點(diǎn))，運(yùn)用式(6)進(jìn)行標(biāo)準(zhǔn)化處理，得到標(biāo)準(zhǔn)化的月平均流量序列和平均含沙量序列(見圖2)。

圖2 標(biāo)準(zhǔn)化處理序列與低頻重構(gòu)序列對(duì)比Fig.2 Comparison between standardized monthly series and reconstruction series of the lowest frequency

3.2 水沙序列趨勢(shì)性小波分析

采用多分辨率分析方法，以db3小波為基本小波，對(duì)水沙序列進(jìn)行分辨率為7(分辨率取log2N，N為樣本容量)的小波分解，提取序列的低頻成分。

徑流與含沙量序列的小波分解重構(gòu)結(jié)果如圖2中的點(diǎn)劃線所示。為直觀起見，作經(jīng)過低頻重構(gòu)序列突變點(diǎn)(2個(gè)序列的突變點(diǎn)均為樣點(diǎn)133處，對(duì)應(yīng)于1987年1月)的水平直線，可以看出在突變點(diǎn)兩側(cè)，2個(gè)序列的趨勢(shì)明顯不同，徑流呈先下降后上升趨勢(shì)，含沙量呈先上升后下降趨勢(shì)。

圖3 年降雨量序列Fig.3 Yearly rainfall series

經(jīng)分析，徑流的趨勢(shì)性改變可能與該流域降雨量變化有關(guān)[10]。根據(jù)李仙江流域?qū)崪y(cè)的降雨量序列(見圖3)可以看出，年降雨量在1976—1987年大致呈下降趨勢(shì)，在1987—1995年大致呈上升趨勢(shì)，徑流的趨勢(shì)性變化基本與降雨量變化同步，二者是正相關(guān)關(guān)系。而含沙量的趨勢(shì)性改變可能與該流域水土保持狀況有關(guān)。在20世紀(jì)80年代中期之前的一段時(shí)期內(nèi)，該流域上游地區(qū)的森林覆蓋率較低，土壤侵蝕較強(qiáng)，且一度有逐年加重的趨勢(shì)，導(dǎo)致徑流的含沙量不斷升高。自80年代中期前后起，該地區(qū)逐漸認(rèn)識(shí)到水土流失的危害性，開始重視水土保持建設(shè)，采取退耕還林、植樹造林的措施，逐漸恢復(fù)并增加森林植被的覆蓋面積，使得徑流含沙量呈下降趨勢(shì)。

3.3 水沙序列周期性小波分析

對(duì)徑流與含沙量序列，經(jīng)過復(fù)Morlet小波變換后可得到小波變換系數(shù)的模及實(shí)部，其實(shí)部的時(shí)頻分布如圖4所示。圖4中圖例是色標(biāo)，亮度最大為1，最小為0，實(shí)部越大，對(duì)應(yīng)的亮度也越大，相應(yīng)的時(shí)間尺度信號(hào)越強(qiáng)。從圖中可看出徑流序列和含沙量序列分別存在時(shí)間尺度為42～54個(gè)月與40～58個(gè)月的周期。

圖4 標(biāo)準(zhǔn)月徑流序列與含沙量序列復(fù)Morlet 小波變換系數(shù)實(shí)部的時(shí)頻分布Fig.4 Real part of CWT for standardized monthly flow series and sediment series using complex valued Morlet wavelet

根據(jù)式(12)計(jì)算得到2個(gè)序列的小波方差,見圖5。

圖5 徑流與含沙量小波方差Fig.5 Flowandsedimentwaveletvariance

由圖5可知，2個(gè)序列小波方差的最大峰值均出現(xiàn)在時(shí)間尺度為48處，表示該時(shí)間尺度下的周期振蕩最強(qiáng)，即可判斷徑流與含沙量序列的大周期應(yīng)該為48個(gè)月，約4a；第2大峰值均出現(xiàn)在時(shí)間尺度為24處，可判斷小周期應(yīng)該為24個(gè)月，約2a。

從實(shí)部的時(shí)頻分布中抽取2個(gè)序列大周期對(duì)應(yīng)的縱坐標(biāo)(時(shí)間尺度)的序列，繪制對(duì)應(yīng)時(shí)間尺度下的小波變換系數(shù)(實(shí)部)圖，結(jié)果如圖6所示。

圖6 徑流與含沙量序列大周期(4a)對(duì)應(yīng)的復(fù)Morlet小波變換系數(shù)的實(shí)部時(shí)頻分布Fig.6 RealpartofcomplexvaluedMorletwaveletcoefficientcorrespondingtothedominantperiod(4years)offlowseriesandsedimentseries

圖6中復(fù)小波變換系數(shù)實(shí)部為正數(shù)時(shí)，代表流量與含沙量在此時(shí)段偏多，對(duì)應(yīng)豐水期；為負(fù)數(shù)時(shí)，代表流量與含沙量在此時(shí)段偏少，對(duì)應(yīng)枯水期?？梢悦黠@看出2個(gè)序列的大致變化規(guī)律是在2年一個(gè)小周期內(nèi)均呈1年大1年小的交替變化，且二者的周期大小變化在時(shí)域上基本同步。

4 結(jié) 語(yǔ)

采用db3小波對(duì)水沙序列進(jìn)行多分辨率分析，研究徑流與含沙量變化趨勢(shì)；選用復(fù)Morlet小波分析水沙序列的隱含周期，結(jié)果合理可信。表明小波分析作為時(shí)頻分析的重要手段能夠較為準(zhǔn)確地反映水沙序列的變化特性，降低序列長(zhǎng)度和隨機(jī)性的影響，進(jìn)而為多沙水庫(kù)運(yùn)行調(diào)度決策提供科學(xué)依據(jù)。

[1] 丁 晶. 隨機(jī)水文學(xué)[M]. 北京: 中國(guó)水利水電出版社, 1997.

[2] 于 浩，張曉萍，李 銳. 延河流域徑流和輸沙周期變化特征的小波分析[J]. 中國(guó)水土保持科學(xué), 2008, 6(4): 18-22.

[3] 王文圣，丁 晶，向紅蓮. 小波分析在水文學(xué)中的應(yīng)用研究及展望[J]. 水科學(xué)進(jìn)展, 2002, 13(4): 515-520.

[4] 劉建梅，王安志，裴鐵璠，等. 雜谷腦河徑流趨勢(shì)及周期變化特征的小波分析[J]. 北京林業(yè)大學(xué)學(xué)報(bào), 2005, 27(4): 49-55.

[5] 劉素一，權(quán)先璋，張勇傳. 不同小波函數(shù)對(duì)徑流分析結(jié)果的影響[J]. 水電能源科學(xué), 2003, 21(1): 29-31.

[6] 孫延奎. 小波分析及其工程應(yīng)用[M]. 北京：機(jī)械工業(yè)出版社，2009.

[7]BRADSHAWGA,MCINTOSHBA.DetectingClimate-inducedPatternsUsingWaveletAnalysis[J].EnvironmentalPollution, 1994, 83(1): 135-142.

[8]MORLETJ,ARENSG,FOURGEAUE,etal. Wave Propagation and Sampling Theory—Part I: Complex Signal and Scattering in Multilayered Media[J]. Geophysics, 1982, 47(2): 203-221.[9] 藺秋生，范北林，黃 莉. 宜昌水文站年徑流量演變多時(shí)間尺度分析[J].長(zhǎng)江科學(xué)院院報(bào),2009,26(4):1-3.[10]薛小杰，王 煜. 小波分析在水文序列趨勢(shì)分析中的應(yīng)用[J]. 應(yīng)用科學(xué)學(xué)報(bào),2002, 20(4): 426-428.

(編輯：陳 敏)

Wavelet Analysis on Flow and Sediment Variation in Sandy Rivers

SHAO Xue-jie, GU Sheng-ping, CAO Ai-wu, HE Hai-xiang, YU Ting-ting

(College of Water Conservancy and Hydropower Engineering, Hohai University,Nanjing 210098, China)

In the aim of providing basis for the operation of sandy reservoir and reducing the influences of length and randomness of flow and sediment series, wavelet transform was adopted to analyze the variation characteristics of flow and sediment series. Standardized monthly flow and sediment series were decomposed by Multi-resolution Analysis using the db3 wavelet function, and continuous wavelet transform was used to evaluate the periodic variations of the two standardized series using the complex valued Morlet function. The flow and sediment series at the dam site of Yayangshan Hydropower Station located in Lixianjiang watershed were taken as an example. The reconstruction of the lowest frequency part revealed the trend of the flow and sediment series. According to the local rainfall and soil and water conservation, the result is considered reasonable. Wavelet variances were obtained to identify the dominant period as 2-year approximately for both flow series and sediment series and reveal the synchronization between them. The results indicate that wavelet transform is effective for nonstationary stochastic series analysis.

wavelet transform; flow and sediment series; runoff; sediment concentration; time scale; multi-resolution analysis; wavelet variances

2016-05-19；

2016-06-26

國(guó)家“十二五”科技支撐計(jì)劃項(xiàng)目(2013BAB06B01)

邵雪杰(1992-)，男，江蘇南通人，碩士研究生，研究方向?yàn)樗娤到y(tǒng)規(guī)劃與工程經(jīng)濟(jì)，(電話)15850603872(電子信箱)1289770074@qq.com。

顧圣平(1957-)，男，江蘇泰州人，教授，碩士，研究方向?yàn)樗娤到y(tǒng)規(guī)劃與工程經(jīng)濟(jì)，(電話)13951893043(電子信箱)spgu@hhu.edu.cn。

10.11988/ckyyb.20160482

2017,34(5):5-8

P333.9

A

1001-5485(2017)05-0005-04