CEEMD與KSVD字典訓練相結合的去噪方法

樂友喜 楊 濤 曾賢德

(①中國石油大學(華東)地球科學與技術學院，山東青島 266580；②黃河勘測規劃設計研究院有限公司,河南鄭州 450003；③中國能源建設集團新疆電力設計院有限公司，新疆烏魯木齊 830001)

0 引言

地震信號中混入隨機噪聲時，會大大降低地震資料的信噪比，嚴重影響對有效信號的檢測和準確估計。研究表明，由于在各個時刻具有不確定性，決定了隨機噪聲呈現在零點處自相關函數值最大，而在其他點處自相關函數值迅速衰減的特點。對于理想高斯白噪聲，其自相關函數值在零點處應表現為一個尖脈沖；而對于其他一般信號，由于存在固有的關聯度，其自相關函數值的變化規律明顯不同于隨機噪聲[1]。

現有的地震信號去噪方法已相當成熟，主要包括Wiener濾波[2]、譜減法[3]、小波變換[4-5]等傳統的(去噪)算法及其改進方法。

KSVD(K Singular Value Decomposition)字典稀疏去噪法屬于基本壓縮感知理論[6-7]的改進算法，該方法通過大量樣本信號訓練獲得過完備冗余字典，能更自適應地以稀疏基表征[8-9]信號，從而更精確地重構信號，實現去噪目的。盡管KSVD字典[10-11]稀疏去噪可使信號中大部分隨機噪聲得到有效抑制或消除，但在信號的各個波峰和波谷處總會殘留一些毛刺[12]。

經驗模態分解(Empirical Mode Decomposi- tion，EMD)可對信號局部時變特征做自適應時頻分解[13-14]，非常適用于非平穩、非線性信號分析，但它存在模態混疊和端點效應問題[15]。

完備總體經驗模態分解(Complete Ensemble Empirical Mode Decomposition，CEEMD)[16]彌補了EMD存在的缺陷，在分解信號的每一階段添加一種特定白噪聲，同時計算唯一的殘差以獲取每個符合定義的固有模態分量(IMF)，這是一個完備的分解過程。CEEMD提供的模態分離效果明顯好于EMD，同時對原始信號的重構更為精確，且計算效率相對更高[17]。

本文針對傳統地震信號的隨機噪聲去除方法的不足及KSVD字典稀疏去噪后信號中總會殘留毛刺的問題，提出一種CEEMD與KSVD字典訓練相結合的去噪方法，以期能夠在實際地震資料中更好地消除隨機噪聲的影響。

1 理論基礎

1.1 CEEMD基本原理

總體經驗模態分解(Ensemble Empirical Mode Decomposition，EEMD)在一定程度上解決了EMD的模態混疊效應，提高了EMD的分解質量；同時也造成了新問題，如加入不同高斯白噪聲會導致分解出的IMF不是唯一解，與IMF的原始定義不能很好吻合，這樣就失去它所代表的物理含義。EEMD算法的完備性不強，對處理后結果做累加難以完美重構原來信號；同時受到加總平均[12]的影響，在很大程度上，噪聲添加越多分解結果越好，但其計算成本也會相應地增加。

Torres等[18]于2011年提出了一種噪聲自動適應的完備總體經驗模態分解(CEEMD)，該算法與EEMD相同點在于也利用了噪聲進行輔助分析，但具體實現過程與EEMD區別很大。此算法不僅解決了模態混疊問題，同時提高了原信號重構的精度、提升了計算效率也降低了成本，總體上明顯優于EMD、EEMD方法。其具體步驟描述如下。

(1)將不同的白噪聲分別與原來信號合成T個混合信號，應用EMD對其處理，計算集合平均并將其作為原來信號的第一個固有模態函數IMF1

(1)

式中：Fj(·)為EMD處理后獲取的第j階模態；ωi為i個高斯白噪聲；εk為每一階段白噪聲加入的比例大小；x(t)為初始信號。

(2)令r0(t)=x(t)，對k=1、…、K計算第k階殘差rk(t)

rk(t)=rk-1(t)-IMFk(t)

(2)

(3)對rk(t)+εkFk[ωi(t)]做EMD處理，獲取對應的IMF1；計算總體平均并將其作為IMFk+1

(3)

(4)重復步驟(2)和步驟(3)，直到殘差信號不能分解為止，得到最終殘差

(4)

1.2 KSVD字典構造

KSVD算法[19]主要包括字典更新和稀疏編碼兩個步驟，獲取字典的同時也獲得稀疏系數矩陣。

在稀疏編碼過程中，常用的稀疏分解方法為正交匹配追蹤算法(Orthogonal Matching Pursuit，OMP)。字典更新的常用方法有MOD(Method of Optimal Directions)和KSVD，兩者差異之處在于字典更新原子時所用更新機制不同： MOD是對字典中所有原子進行更新，而KSVD只針對字典中某個原子進行。

在字典[20]更新過程中，KSVD學習字典使用了奇異值分解更新字典原子，與固定基字典相比，訓練的字典具有更能反映信號特征的優點。KSVD算法具體流程(圖1)如下。

圖1 KSVD稀疏去噪方法框圖

參數初始化： 給定訓練樣本X=[x1,x2,…,xN]， 初始超完備字典D， 正則化參數λ和迭代次數J。

(5)

(3)達到迭代次數J即結束，否則返回步驟(1)。

經過上述步驟(1)～步驟(3)，即可得到學習后的KSVD字典，在整個KSVD運算中，每次字典更新時都能保證總誤差單調下降，因此該算法具有較好收斂特性。此外，KSVD算法還具有一定的去噪功能。

2 CEEMD與KSVD聯合去噪方法

CEEMD去噪方法盲目舍棄高頻IMF模態分量有可能導致高頻有效信號損失或高頻噪聲壓制不完全。KSVD字典稀疏去噪盡管去除了大部分噪聲，但弱隨機噪聲尚存在于局部信號中，尤其在信號各個波峰和波谷處尚殘留一些較明顯毛刺。鑒于此，本文將CEEMD與KSVD相結合，通過互補方式彌補上述方法存在的兩點不足。

CEEMD-KSVD聯合去噪具體實現步驟(圖2，其中“過渡模態”是指介于噪聲主導模態與信號主導模態之間的模態)如下：

(1)采用CEEMD對觀測含噪信號y模態分解獲得本征模態函數集。

(2)求取每個IMF分量的峰度，識別噪聲主導模態，過渡模態和信號主導模態。

(3)噪聲主導模態分量：直接舍棄。

(4)過渡模態分量： ①對過渡模態分量累加合成新信號x，并對其二次CEEMD分解；②利用各個模態分量求取的峰度，舍棄噪聲主導模態；③對步驟②中的剩余IMF分量累加合成新信號x1，并與觀測信號y作為訓練樣本，獲得KSVD學習字典；④利用上步字典及OMP算法對信號x1稀疏去噪獲得去噪結果y1。

(5)信號主導模態分量： ①各個信號主導模態累加重構得到新信號z；②將信號z和y作為訓練樣本，獲取KSVD學習字典；③利用學習字典及OMP算法對信號z稀疏去噪得到y2。

圖2 本文去噪方法流程框圖

(6)將步驟(4)和步驟(5)中的去噪結果y1和y2重構合并，得到最終去噪結果。

峰度又稱峰態系數，是用于表征數據的概率密度分布曲線在平均值處峰值高低的特征數[21]。對于服從高斯分布的隨機噪聲信號，其峰度F接近0；對于服從超高斯分布的隨機噪聲信號，其峰度大于0；反之， 服從亞高斯分布的隨機噪聲信號，其峰度小于0。

下文中峰度評判標準： 噪聲主導模態是|F|≤0.1， 過渡模態是0.1<|F|≤1.0， 信號主導模態是|F|>1.0。

3 模型試算

設計如圖3所示的褶積模型，該模型共有三層，上下兩層為水平層，中間層為傾斜層，并與上覆層斜交。橫向共有25道，縱向上有300個采樣點，采樣頻率為500Hz，雷克子波主頻為25Hz，模型中加入兩種不同比例高斯白噪聲，分別采用五種方法對含噪模擬地震數據進行去噪，通過分析驗證五種方法去噪效果。

圖3 無噪模擬地震記錄

特別說明： 本文所使用的剖面數據的字塊為8×8，字典原子的維數是64×256；由于過渡模態部分合成信號的隨機噪聲較強，使用的迭代次數為20，信號主導模態部分合成信號的隨機噪聲相對弱小，為提高計算效率，使用的迭代次數為10；同時分割子塊尺寸和字典維數不變時，迭代次數越多，字典就會變得越稀疏，計算效率也就越高；反之迭代次數不變時，子塊尺寸和字典維數越大，計算效率越低。本文加入的高斯白噪聲頻譜范圍分布廣泛，包含地震頻帶(圖4)。

3.1 信噪比RS/N=0.88

加入較高比例的高斯白噪聲構成信噪比為0.88的加噪模擬地震記錄(圖5a)， 圖5b～圖5f為五種方法的去噪結果。從圖中可看到：F-X域去噪和CEEMD分頻去噪效果最差；小波軟閾值去噪效果也不太理想；KSVD稀疏去噪效果較好，但尚存一些毛刺現象；本文方法基本保存了原始信號的有效信息，去噪效果最好。表1為五種方法去噪結果評價指標對比。從各方法去噪結果及其評價指標的對比得知，五種方法中后兩種去噪效果較好，尤其本文方法相比其他方法去噪后信噪比最高，有效信息損失最少，是一種有效保幅的去噪方法。信噪比計算公式如下

圖4 高斯白噪聲(a)及其頻譜特征(b)

表1 五種方法去噪結果評價指標對比

圖5 含較高比例噪聲模擬地震記錄五種方法去噪效果對比

(6)

式中：S1為去噪后的數據；S為無噪數據；RS/N為去噪后信噪比。

3.2 信噪比RS/N=3.13

加入較低比例的高斯白噪聲構成信噪比為3.13的加噪模擬地震記錄(圖6a)，圖6b～圖6f為五種方法的去噪結果。從圖中可知： 五種方法都能將大部分噪聲壓制，但F-X域去噪和CEEMD分頻去噪有明顯的噪聲殘余；小波軟閾值去噪沒有毛刺現象，但在局部改變了信號的幅值；KSVD稀疏去噪效果較好，與本文方法去噪效果相差無幾，但是仔細觀察仍存在輕微的毛刺現象。表2為五種方法去噪結果評價指標對比。因此，綜合對比五種方法的去噪結果及其評價指標可看出：對高信噪比資料，后三種方法去噪效果均較好，尤其KSVD稀疏去噪與本文方法去噪效果只有些許差別。

圖6 含較低比例噪聲模擬地震記錄五種方法去噪效果對比

表2 五種方法去噪結果評價指標對比

通過上述分析可得到以下認識： F-X域去噪和CEEMD分頻去噪不管是在高信噪比還是低信噪比資料中盡管能壓制一定程度的噪聲，但都未獲取理想去噪結果；小波軟閾值去噪和KSVD稀疏去噪效果質量與信噪比呈線性關系，信噪比越高，效果越好；本文方法實用性強，在不同信噪比資料中都能獲取較好的去噪效果。

4 實際地震資料應用

為了進一步測試本文去噪方法對實際地震資料的處理，將其應用于實際地震資料。

圖7a為某一工區的疊后地震記錄，該地震記錄的CDP數目為250，采樣點個數為500，采樣時間為2ms，從圖中可見地震剖面上由于存在大量的隨機噪聲致使剖面同相軸不連續或錯斷，剖面顯示不清晰，對地質解釋帶來諸多不便。為了消除隨機噪聲的影響，分別利用上述五種方法進行去噪處理，對比這五種方法去噪后的成果剖面(圖7b～圖7f)可知： F-X域去噪與CEEMD分頻去噪效果最差，去噪后還存在許多噪聲殘留；小波軟閾值去噪方法除掉了大量隨機噪聲，仍存有少許隨機噪聲，去噪效果一般；KSVD稀疏去噪與本文方法去噪效果相比前三種較好，但KSVD稀疏去噪在一些細節(弱幅度有效信號)上不如本文方法，由圖7d和圖7e可知，本文方法充分壓制了隨機噪聲，同相軸相對更加連續，有效信息得以保留，地震剖面清晰度提高。圖8為各種去噪方法獲得的殘差剖面圖。由殘差剖面可知，F-X域、KSVD稀疏、小波等去噪法都或多或少地去除了部分有效信號，尤其F-X域去噪損失有效信號最嚴重；CEEMD分頻去噪，當舍棄模態IMF1時，去除的基本上是隨機噪聲，但是去除IMF1+IMF2時，導致一部分有效信號的丟失，而本文方法去噪不僅隨機噪聲得以完全壓制，同時又不傷害有用信號，取得了很好的去噪效果。

圖7 針對實際資料五種方法去噪效果對比

圖8 五類去噪方法對應的殘差剖面

5 認識和結論

(1)本文研究表明： F-X域去噪方法能去除部分隨機噪聲，但會導致假同相軸、有效波波形畸變；由于閾值很難精準選取，小波去噪方法在去除大量噪聲的同時也丟失了部分有效信號；CEEMD去噪方法因盲目舍棄IMF分量可能造成高頻有效信息缺失；KSVD稀疏去噪總有殘留毛刺，且當資料信噪比較低時，去噪效果也隨之變差；本文提出的CEEMD-KSVD聯合去噪方法的去噪效果優于其他四種方法，能顯著提高地震資料品質。

(2)隨機噪聲的頻率大多較高，CEEMD分頻去噪方法通過舍棄前若干個模態，主要用于去除高頻部分的隨機噪聲，KSVD稀疏去噪方法是以隨機噪聲與信號在整個定義域上稀疏分解系數分布差異較大為基礎，兩種方法都立足于隨機噪聲的特性。因此，本文去噪方法更適用于隨機噪聲的壓制，對于其他類型噪聲，其適用性有待測試。