基于CEEMDAN和小波變換的地震信號(hào)隨機(jī)噪聲壓制新方法

郭曉菲 歐同庚 馬武剛 吳林斌 趙義飛 劉 軍 徐春陽(yáng)

1 中國(guó)地震局地震大地測(cè)量重點(diǎn)實(shí)驗(yàn)室，武漢市洪山側(cè)路40號(hào)，430071

2 武漢地震科學(xué)儀器研究院有限公司，湖北省咸寧市青龍路11號(hào)，437099

隨機(jī)噪聲壓制或去噪是地震數(shù)據(jù)分析與應(yīng)用的重要前提?；诮?jīng)驗(yàn)?zāi)B(tài)分解(EMD)[1-3]的信號(hào)分解技術(shù)在地震信號(hào)去噪領(lǐng)域應(yīng)用廣泛，但EMD存在模態(tài)混疊現(xiàn)象。為了改進(jìn)該問(wèn)題，基于白噪聲平衡的集合經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD)[4-5]逐漸受到重視。但單純應(yīng)用EMD或EEMD進(jìn)行信號(hào)去噪，往往直接剔除首條高頻IMF分量、最后一條低頻IMF分量或剩余分量，對(duì)其他分量未作更加精細(xì)的去噪分析，這樣可能遺漏地震數(shù)據(jù)關(guān)鍵信息，且重構(gòu)后的信號(hào)還原效果不理想。 CEEMDAN[6-8]是一種基于EMD分解并借鑒EEMD中添加白噪聲來(lái)平衡信號(hào)隨機(jī)干擾的新型信號(hào)分解方法，能有效抑制EMD分解時(shí)模態(tài)混疊的出現(xiàn)，克服殘余噪聲在分解過(guò)程中定向傳遞的問(wèn)題，并保證信號(hào)不失真，具有良好的自適應(yīng)性、還原性和完備性。

本文基于皮爾斯相關(guān)系數(shù)(Pearson correlation coefficient，PCCs)，提出一種結(jié)合CEEMDAN和小波變換(wavelet transform，WT)的地震信號(hào)噪聲壓制新方法。并通過(guò)模擬實(shí)驗(yàn)和2021年青?，敹?.4級(jí)地震數(shù)據(jù)驗(yàn)證其有效性。

1 信號(hào)去噪、重構(gòu)與效果評(píng)估

1.1 CEEMDAN-WT方法去噪流程

CEEMDAN-WT方法的去噪流程如圖1所示，具體步驟如下：

1)設(shè)原始信號(hào)為x(n)，對(duì)其進(jìn)行CEEMDAN分解，得到m個(gè)IMFi(i=1,2,…,m)分量和1個(gè)剩余信號(hào)R，上述分量按照頻率高低依次排列。

2)計(jì)算各IMFi分量及余量R與初始信號(hào)x的皮爾森相關(guān)系數(shù)PCCsj(j=1,2,…,m+1)。

3)根據(jù)相關(guān)性強(qiáng)弱劃分原則，分別設(shè)定α、β為弱相關(guān)閾值和顯著相關(guān)閾值，將|PCCsj|∈(α,β)的K個(gè)分量分別作小波變換[9]處理，剔除高頻成分，得到wk(k=1,2,…,K)。

4)|PCCsj|∈[β,1]的分量與初始信號(hào)相關(guān)性顯著，且基本為高頻分量，為了更多地保留有效成分，不作任何處理，將P個(gè)該分量記為Gp(p=1,2,…,P)；|PCCsj|∈[0,α]的分量與初始信號(hào)相關(guān)性較低，且基本為低頻分量，予以刪除。

5)信號(hào)線性重構(gòu)：

(1)

1.2 去噪效果定量評(píng)估

分別從去噪有效性(信噪比、樣本熵變化量)、信號(hào)還原度(皮爾森相關(guān)系數(shù)、互信息、均方誤差)、算法效率(去噪速度)等3個(gè)方面設(shè)計(jì)指標(biāo)定量評(píng)估信號(hào)去噪的效果。各指標(biāo)物理意義如下：

1) 信噪比(SNR)。SNR是衡量地震信號(hào)去噪質(zhì)量的主要指標(biāo)，是信號(hào)中有效成分與噪聲的相對(duì)程度，去噪目標(biāo)之一即是提高SNR值。

2) 互信息(mutual information，MI)。MI是信息論中表示2個(gè)信息之間相關(guān)性的指標(biāo)。

3) 樣本熵[10]變化量(variety of sample entropy，VSE)。樣本熵可用來(lái)描述信號(hào)的復(fù)雜狀態(tài)，值越大，信號(hào)復(fù)雜程度越高，而VSE表示原始含噪信號(hào)與去噪后信號(hào)的樣本熵值之差，其值為正，則表示隨機(jī)噪聲得到一定抑制，信號(hào)復(fù)雜性有所降低及波形更加平滑；反之，去噪效果不理想。

4) 皮爾森相關(guān)系數(shù)(PCCs)。PCCs用來(lái)反映兩個(gè)對(duì)象之間的相關(guān)程度，其值介于-1與+1之間，表征去噪后信號(hào)的還原程度。

5) 均方誤差(MSE)。MSE表示原始含噪信號(hào)與去噪后信號(hào)間的差異程度。

6) 去噪時(shí)長(zhǎng)。去噪時(shí)長(zhǎng)為信號(hào)分解、去噪與重構(gòu)的計(jì)算耗時(shí)。

2 仿真實(shí)驗(yàn)

為了驗(yàn)證CEEMDAN去噪效果的有效性，首先利用MATLAB 2019A平臺(tái)生成一個(gè)簡(jiǎn)單的仿真地震波形，加入6 dB的高斯白噪聲(圖2)，然后對(duì)含噪波形信號(hào)分別進(jìn)行EMD、EEMD、CEEMDAN，以及EMD-WT、EEMD-WT、CEEMDAN-WT去噪處理。前3種方法直接對(duì)首階高頻IMF分量進(jìn)行剔除后線性重構(gòu)，而EMD-WT、EEMD-WT的去噪原理類(lèi)似于CEEMDAN-WT去噪。

實(shí)驗(yàn)參數(shù)設(shè)置為：EMD分解維數(shù)為8；EEMD的白噪聲標(biāo)準(zhǔn)差為0.15，白噪聲添加數(shù)目為20；CEEMDAN的噪聲標(biāo)準(zhǔn)差為0.2，噪聲數(shù)量為20，最大迭代次數(shù)為3 600；α、β分別設(shè)置為0.20、0.85；小波基函數(shù)采用db4；樣本熵參數(shù)值為0.15(STD)。

圖3為EMD、EEMD及CEEMDAN分解得到的各IMF分量與初始純凈信號(hào)間相關(guān)系數(shù)值的分布?？梢钥闯?，CEEMDAN得到的與初始信號(hào)高度相關(guān)的IMF分量數(shù)目較多，分解效率較高，便于后期的高低頻成分分析與有效成分提取。

圖4給出了6種方法的去噪結(jié)果，表1為6種去噪方法的評(píng)價(jià)指標(biāo)結(jié)果。可以看出，EEMD和EEMD-WT去噪效果一般，其余4種方法去噪效果較明顯，特別是CEEMDAN-WT方法，將初始含噪信號(hào)的信噪比由6 dB提升到了15.696 0 dB，且其MI、PCCs和MSE指標(biāo)結(jié)果也最佳。

需要特別注意的是，EEMD方法雖然在波形復(fù)雜度降低方面效果較為突出(VSE=1.215 6)，但去噪后信號(hào)失真較嚴(yán)重，信噪比結(jié)果未達(dá)預(yù)期，其對(duì)信號(hào)與噪聲的分離效果不佳。同時(shí)，從圖表信息可觀察到，EMD-WT方法幾乎在每一個(gè)去噪性能指標(biāo)數(shù)值上都接近于CEEMDAN-WT方法，且其去噪時(shí)間更短。因此，在對(duì)數(shù)據(jù)處理效率要求較高的場(chǎng)合，EMD-WT方法可以作為CEEMDAN-WT方法的可靠替代方案。

圖5為仿真地震波信號(hào)的去噪結(jié)果頻譜圖?？擅黠@看出，CEEMDAN-WT與EMD-WT去噪頻譜結(jié)果的高頻噪聲成分(10～50Hz)得到有效抑制，其他4種去噪方法相對(duì)表現(xiàn)較差，仍存在大量的高頻噪聲，尤其是EEMD-WT方法，這與表1的分析結(jié)論一致。

表1 仿真地震信號(hào)去噪效果對(duì)比

3 地震信號(hào)去噪實(shí)驗(yàn)

為了進(jìn)一步驗(yàn)證CEEMDAN-WT算法的去噪效果，應(yīng)用2021-05-22發(fā)生于瑪多的M7.4地震強(qiáng)震動(dòng)記錄進(jìn)行去噪實(shí)驗(yàn)，該震為我國(guó)近幾年震級(jí)最大、余震較多的一次典型地質(zhì)破壞事件，地震烈度X度區(qū)域半徑達(dá)69 km，具有顯著的研究?jī)r(jià)值。地震數(shù)據(jù)由中國(guó)地震局工程力學(xué)研究所負(fù)責(zé)觀測(cè)記錄，包括西寧、湟源、民和等16個(gè)地震臺(tái)站(圖6)記錄的三分量加速度數(shù)據(jù)，共48條。

表2為48個(gè)地震樣本數(shù)據(jù)實(shí)施CEEMDAN-WT去噪的總體結(jié)果。從VSE來(lái)看，有97.916 7%(47/48)的樣本的波形隨機(jī)性得到抑制，PCCs與MI指標(biāo)皆在0.930 0與0.830 0以上，其均值分別達(dá)到0.982 7和0.872 4，去噪重構(gòu)后的有效成分保留較好。

圖7為48個(gè)地震數(shù)據(jù)樣本在不同區(qū)間下CEEMDAN各分量與初始樣本相關(guān)系數(shù)的分布結(jié)果，統(tǒng)計(jì)區(qū)間閾值分別為0.2、0.5、0.8。由圖可知，分解出顯著相關(guān)(PCCs>0.8)分量的概率是27.083%(13/48)，分解出較大相關(guān)(PCCs>0.5)分量的總體概率是100%(48/48)；對(duì)于單個(gè)地震樣本得到的全部9個(gè)分量而言，CEEMDAN能準(zhǔn)確分離出該地震樣本的有效成分，其數(shù)量在1～4之間。地震信號(hào)有效分量的提取除了利于進(jìn)行后續(xù)的CEEMDAN-WT高頻成分去噪外，還有助于進(jìn)行地震波形多尺度特征的提取與地震事件屬性辨識(shí)，對(duì)其他基于信號(hào)處理的地震學(xué)研究也有一定的參考價(jià)值。

表2 基于CEEMDAN-WT的瑪多地震信號(hào)去噪效果

圖8為地震樣本1(表2)的CEEMDAN-WT去噪結(jié)果?？梢钥闯觯卣鹦盘?hào)經(jīng)過(guò)CEEMDAN分解得到9個(gè)IMF分量，IMF1～I(xiàn)MF9的中心頻率依次下降，波形光滑性與周期性逐漸增強(qiáng)，即高頻隨機(jī)噪聲的影響持續(xù)減弱。根據(jù)IMF分量波形特征以及皮爾森相關(guān)系數(shù)計(jì)算結(jié)果，IMF6～I(xiàn)MF8的PCCs值分別為0.731 1、0.585 6、0.563 6，可以推測(cè)其為地震波真實(shí)成分，其余PCCs值較低的IMF分量可能為高頻噪聲成分或低頻噪聲成分。經(jīng)過(guò)CEEMDAN-WT去噪、線性重構(gòu)后，地震波形至少在8個(gè)局部區(qū)域(圖中虛線圓圈標(biāo)注)受到明顯的光滑濾波處理，并與原信號(hào)保持較為一致的波形形狀，該樣本的去噪指標(biāo)MI、VSE、PCCs、MSE分別為0.878 1、0.045 4、0.985 7、0.023 4(表2)。

4 結(jié) 語(yǔ)

1)相對(duì)于EMD-WT、EEMD-WT等5種去噪方法，CEEMDAN-WT方法在各個(gè)去噪效果評(píng)價(jià)指標(biāo)上皆取得最好的噪聲壓制效果，SNR從6 dB提升到了15.696 dB，并最大程度地保留了原始地震仿真信號(hào)的有效成分，MI、PCCs和 MSE 分別為 0.878 8、0.988 3和0.0461；在保證去噪效果良好的前提下，如果對(duì)去噪時(shí)效性要求較高，EMD-WT方法可作為CEEMDAN-WT方法的代替方案。

2)2021青?，敹郙7.4地震的天然地震信號(hào)去噪實(shí)驗(yàn)結(jié)果整體符合預(yù)期目標(biāo)：97.916 7%的樣本波形復(fù)雜程度得到一定抑制，信號(hào)還原性較好，PCCs與MI指標(biāo)均值分別達(dá)到0.982 7和0.872 4；CEEMDAN-WT方法分離出PCCs>0.5的IMF分量的概率達(dá)到100%，這有助于對(duì)指定IMF分量進(jìn)行高頻去噪；在保證地震波形整體變化趨勢(shì)無(wú)明顯變化的同時(shí)，該方法實(shí)現(xiàn)了一定的局部光滑去噪效果。

3)考慮到臺(tái)站場(chǎng)地條件、震源深度和震級(jí)等因素對(duì)隨機(jī)噪聲的影響，下一步將探索分析CEEMDAN-WT方法在各類(lèi)地震案例樣本中的表現(xiàn)。

致謝：感謝中國(guó)地震局工程力學(xué)研究所提供數(shù)據(jù)。