沙 磊

1 安徽省地質(zhì)實驗研究所(國土資源部合肥礦產(chǎn)資源監(jiān)督檢測中心),合肥市阜陽北路318號,230041

變形監(jiān)測數(shù)據(jù)一般包含真值和噪聲2個部分,對變形監(jiān)測數(shù)據(jù)進行去噪處理、提取變形特征、分析變形規(guī)律,是變形分析研究的重要內(nèi)容[1]。小波分析具有良好的時頻特性,被廣泛應用于變形分析中,特別是小波閾值去噪在非平穩(wěn)變形監(jiān)測數(shù)據(jù)處理中取得了較好的效果[2-3]。眾多學者從影響小波去噪的關鍵因素出發(fā),如小波分解層數(shù)[4-5]、去噪評價[6-7]、閾值函數(shù)[8-10]、閾值估計[8-11]等對傳統(tǒng)小波閾值去噪進行研究。在閾值估計方面,常用的4種閾值分別為通用閾值、啟發(fā)式閾值、自適應閾值和極大極小閾值。其中,自適應閾值和極大極小閾值選取規(guī)則比較保守,僅將部分系數(shù)置0;而通用閾值和啟發(fā)式閾值選取規(guī)則的去除噪聲效果更佳,但易將有用信號的高頻部分當作噪聲信號除去[9]?？紤]到噪聲在小波變換下隨尺度增大而減小的特性,學者們給出了改進閾值估計算法,取得了一定效果,但也存在一些不足。如李明峰等[8]給出的閾值算法首層閾值估計過小,去噪不徹底;甘若等[10]、趙瑞珍等[11]給出的閾值算法首層閾值估計過大,易產(chǎn)生“過扼殺”小波系數(shù)的現(xiàn)象,且上述文獻都沒有考慮帶噪小波系數(shù)含噪幅度隨尺度增加而減小的函數(shù)關系。基于以上問題,本文提出一種新的小波閾值估計算法,并通過理論分析、工程實例驗證其優(yōu)越性和有效性。

1 小波閾值去噪法

小波閾值去噪是基于小波變換的“集中”或“聚焦”能力。通常認為,信號在經(jīng)過小波變換后,由信號產(chǎn)生的小波系數(shù)包含有信號的重要信息,其幅值大、數(shù)目少,而噪聲對應的小波系數(shù)幅值小、數(shù)目多。通過在不同分解尺度上選取合適的閾值,并將小于該閾值的小波系數(shù)置0,保留大于閾值的小波系數(shù),使信號噪聲得到有效抑制,最后經(jīng)過小波逆變換,得到去噪后的重構(gòu)信號。

1.1 ZR閾值法

考慮到噪聲在小波變換下隨尺度增大而減小的特性[10],給出小波閾值法的改進法(稱為ZR閾值法)[11]:

(1)

1.2 1/2閾值法

對于均值為0、方差為σ2的高斯白噪聲n(t),其二進小波變換為WTj,kn(t),若小波函數(shù)ψj,k(t)為實函數(shù)時,則有[4]:

|WTj,kn(t)|2=|?Rn(t)ψj,k(t)dt|2

(2)

由于E[n(t)n(t)]=σ2,小波函數(shù)ψj,k(t)為小波基ψ(t)伸縮平移所得,小波函數(shù)中每一個都有規(guī)范化的范數(shù)‖ψj,k(t)‖=‖ψ‖=1,故:

(3)

即隨著小波分解尺度的增加,含有噪聲的觀測序列在低頻信號中噪聲影響幅度按1/2的規(guī)律逐漸減小[12]。隨著分解尺度的增加,高一級尺度上的低頻信號中含有的噪聲影響幅度會越來越小,曲線也會也來越光滑。

對變形觀測數(shù)據(jù)序列進行如式(4)所示的小波分解與重構(gòu):

an+1=an+dn=an-1+dn-1+dn=

…=a0+d1+d2+…+dn

(4)

式中,an為低頻信號,表現(xiàn)為趨勢項,dn為高頻信號,含有噪聲、粗差等局部變形信息。

設an噪聲影響幅度為σ2,由以上分析規(guī)律可知,an-1噪聲影響幅度為σ2/2,由式(4)可知,dn-1=an-an-1,即an噪聲影響幅度σ2在高一級尺度上一分為二,也就是an-1分解到的噪聲影響幅度為σ2/2,dn-1分解到的噪聲影響幅度為σ2/2。同理,dn-2=an-1-an-2,dn-2分解到的噪聲影響幅度為σ2/4。因此,隨著分解尺度的增加,高頻信號分解到的噪聲影響幅度也是按1/2的規(guī)律逐漸減小。

綜上所述,噪聲在小波變換下隨著分解尺度的增大而減小,且高頻信號分解到的噪聲影響幅度也是按1/2的規(guī)律逐漸減小。據(jù)此,本文給出新閾值計算公式(稱1/2閾值法):

(5)

2 仿真實驗

多期的變形監(jiān)測信號一般包括趨勢項、周期項及多種噪聲,為保證與實際觀測數(shù)據(jù)相符合,構(gòu)造復合信號進行仿真模擬實驗。如圖1所示,仿真信號時長為120 s,其中包括2個趨勢項、2個周期項及隨機噪聲項,實驗采用MATLAB中自帶的awgn函數(shù),添加信噪比為db10的高斯白噪聲,其表達式為:

圖1 仿真數(shù)據(jù)Fig.1 Simulation data

Y=0.000 6t2-0.185 2t-

0.649 9sin(2πt/12+0.153 7)+

0.147 7sin(2πt/20+0.424 0)-

2.038 0+noise

為驗證本文1/2閾值法的降噪效果,采用分解能力較好的sym8小波基,使用文獻[4-5]的方法確定該信號最佳分解尺度為3層。使用硬閾值方法、軟閾值法、ZR閾值法進行對比,使用信噪比(SNR)、均方根誤差(RMSE)作為小波去噪的評價指標。

消噪的信號如圖2所示?？梢钥闯?軟閾值法去噪后曲線過于光滑,未能有效提取周期性信號,去噪效果最差;硬閾值法由于閾值函數(shù)的不連續(xù)性,在采樣序號75處有較大折線現(xiàn)象,即偽吉布斯現(xiàn)象,對周期性信號的提取有所提高,但在采樣序號5～10處未能提取到信息;ZR閾值法去噪信號曲線光滑,無多余折線現(xiàn)象,在采樣序號5～10處也未能提取到信息;1/2閾值法去噪信號曲線光滑程度最好,無多余振蕩,對周期性信號提取最好。從表1可知,1/2閾值法去噪的SNR最大,RMSE最小,去噪效果最佳,說明1/2閾值法去噪是可行的。

表1 4種閾值法的去噪指標Tab.1 Denoising indicators of four threshold methods

圖2 4種方法去噪結(jié)果Fig.2 Denoising results of four methods

3 工程實例

為進一步驗證本文小波閾值法去噪的有效性和優(yōu)越性,對文獻[10]中某江岸危巖上H監(jiān)測點10 a的沉降監(jiān)測數(shù)據(jù)進行實例分析,如圖3所示。可以看出,該危巖H點的高程隨時間的增加呈現(xiàn)出下降趨勢,同時含有周期性變形信息。由于各種不確定的因素,該信號含有很多毛刺,嚴重影響變形監(jiān)測的準確分析和后續(xù)預測,因此需要對其進行去噪處理。

圖3 原始變形監(jiān)測數(shù)據(jù)Fig.3 Original measuring deformation data

為了驗證1/2閾值法去噪的優(yōu)越性,采用分解能力較好的sym8小波基,使用文獻[4-5]的方法確定該信號最佳分解尺度為3層,同時使用硬閾值法、軟閾值法、ZR閾值法進行對比,去噪結(jié)果如圖4所示?？梢钥闯?軟閾值法由于估計小波系數(shù)和帶噪小波系數(shù)間有恒定的偏差,導致重構(gòu)的信號失真比較嚴重,去噪后曲線過于光滑,未能有效保留周期性細節(jié)變形信息,去噪效果最差。硬閾值法由于閾值函數(shù)的不連續(xù)性,在采樣序號43處有較大折線現(xiàn)象,即偽吉布斯現(xiàn)象,去噪后曲線光滑性最差,對周期性細節(jié)變形信息保留有提高,如在采樣序號20～40處保留了周期性細節(jié)變形信息。ZR閾值法克服了軟、硬閾值法的缺點,也考慮到噪聲在小波變換下隨尺度增大而減小的特性,去噪后曲線的光滑性好,對周期性細節(jié)變形信息的保留程度有所提高,如在采樣序號50～60處保留了周期性細節(jié)變形信息。1/2閾值法克服了軟、硬閾值法的缺點,不僅與噪聲在小波變換各尺度上的傳播特性相一致,也更符合帶噪小波系數(shù)隨分解尺度增加而減小的函數(shù)關系,去噪后曲線的光滑性好,沒有多余振蕩,特別是在采樣序號110～120處,對周期性細節(jié)變形信息保留最好。從表2可知,1/2閾值法去噪的SNR最大,RMSE最小,與其他閾值法去噪對比,效果更好,對周期性細節(jié)變形信息提取更準確。

表2 4種閾值法的去噪指標Tab.2 Denoising indicators of four threshold methods

圖4 4種方法去噪結(jié)果Fig.4 Denoising results of four methods

4 結(jié) 語

通過理論分析與工程實例可知,本文1/2閾值法去噪效果優(yōu)于軟、硬閾值法及ZR閾值法,能更有效地恢復被噪聲湮沒的周期性細節(jié)變形信息,進一步提高去噪的精度與可靠性,更好地反映危巖變形規(guī)律。1/2小波閾值法不僅與噪聲在小波變換各尺度上的傳播特性相一致,也更符合帶噪小波系數(shù)隨分解尺度的增加而減小的函數(shù)關系。需要指出的是,變形監(jiān)測數(shù)據(jù)具有受噪聲干擾相對較小的特點,對噪聲干擾較大的復雜信號,本文方法的去噪效果有待進一步研究。