完全互補小波噪聲輔助集總經驗模態(tài)分解

何 劉， 丁建明， 林建輝， 劉新廠

(西南交通大學 牽引動力國家重點實驗室， 成都 610031)

完全互補小波噪聲輔助集總經驗模態(tài)分解

何 劉， 丁建明， 林建輝， 劉新廠

(西南交通大學 牽引動力國家重點實驗室， 成都 610031)

經驗模態(tài)分解(EMD)是一種自適應非線性非平穩(wěn)數據處理方法。噪聲輔助的EMD方法能克服EMD方法在處理間歇信號時出現的“模態(tài)混疊”現象。在這些噪聲輔助方法中，互補集總經驗模態(tài)分解(CEEMD)和完全噪聲輔助噪聲集總經驗模態(tài)分解(CEEMDAN)恢復了EMD分解的完整性。在現有分析方法上提出了完全互補小波噪聲輔助集總經驗模態(tài)分解(CCWEEMDAN)算法。該算法能用更小的集總數、更少的迭代次數和極小的計算消耗獲得更好的光譜分離效果和數目較少的篩選模態(tài)。

經驗模態(tài)分解；集合經驗模態(tài)分解；噪聲輔助；模態(tài)混疊；互補集總經驗模態(tài)分解

經驗模式分解(Empirical Mode Decomposition，EMD)[1-2]是一種非線性非平穩(wěn)信號自適應分解的信號處理方法。該方法能根據數據自身的時間尺度自適應的將信號分解成具有實際物理意義的各個本征模態(tài)函數(Intrinsic Mode Function，IMF)分量和分解殘差，將所有模態(tài)函數和最終趨勢相加能完美重構原始數據。

EMD方法本身也存在一些問題，例如同一模態(tài)中出現不同尺度或頻率的信號，或者同一尺度或頻率的信號被分解到多個不同的模態(tài)中，該問題被稱為“模態(tài)混疊”[3]。為了克服EMD方法的“模態(tài)混疊”問題，集總經驗模態(tài)分解(Ensemble Empirical Mode Decomposition，EEMD)方法被提出。通過在分解前添加高斯白噪聲至原始信號中，利用高斯白噪聲的頻譜均勻分布特性和EMD方法的二進濾波特性填充整個時頻空間[4]，有效的解決了模態(tài)混疊問題。EEMD方法在解決模態(tài)混疊問題的同時卻產生了新的問題。①由于EEMD方法添加了不同量級的高斯白噪聲，最終分解可能得到不同數量的模態(tài)；②由于添加高斯白噪聲的原因，使得分解得到的各個模態(tài)和最終趨勢中存在噪聲，這使得該方法不能精確重建原始數據。③由于高斯白噪聲的引入需要相應集總次數的EMD分解求平均以抵消噪聲影響，使得迭代次數增加、計算效率下降。

針對EEMD方法不能精確重建的問題，互補集總經驗模態(tài)分解(Complementary Ensemble Empirical Mode Decomposition，CEEMD)[5]通過向待分析信號添加成對相反的高斯白噪聲減小由白噪聲引起的重構誤差，并且提高了EEMD的計算效率。然而該方法依然不能解決添加不同噪聲的信號產生不同數量模態(tài)分解的問題。完全噪聲輔助集總經驗模態(tài)分解(Complete Ensemble Empirical Mode Decomposition with Adaptive Noise，CEEMDAN)[6]的提出是對EEMD方法的重要發(fā)展，該方法實現了分解信號的近似完美重構以及解決了不同加噪信號產生不同數量模態(tài)的缺陷。盡管如此，CEEMDAN依然存在迭代次數多，計算效率低的問題。本文主要解決問題就是減少CEEMDAN方法的迭代次數，提高計算效率。

1 EMD, EEMD, CEEMD and CEEMDAN理論

EMD將信號分解為若干個基本模式分量IMF和一個余項?；灸Ｊ椒至勘仨殱M足兩個基本：①在數據列中，數據的極值點數量和過零點數量必須相等或最多相差一個；②任何一點，極大值和極小值的包絡的平均值為0。該算法的詳細過程為

步驟1 令k=0，計算信號r0=x的極大值和極小值。

步驟2 采用三次樣條對信號rk所有極大值點、極小值點進行插值得到上下包絡emax(emin)。

步驟3 計算并得到上下包絡的平均值m=(emax+emin)/2。

步驟4 根據式(1)計算固有模態(tài)函數

rk-m=hk+1

(1)

步驟5 判斷hk+1是否滿足IMF的條件,如果不滿足，對hk+1重復步驟2～步驟4得到滿足條件的包絡平均值。

步驟7 重復上述步驟知道得到滿足停止條件的最終趨勢項rk。

EEMD分解的計算框架和EMD分解基本相同，將原始信號x篩分為IMF。在信號x篩分過程中添加不同白噪聲ω(i)(i=1,…,I)，I為集總次數，為每次EMD分解后剩余分量的時域分布提供一致的參考結構。其具體分解步驟為

步驟1 對原始數據添加白噪聲，形成總體分解信號，其(i)中ε>0。

x(i)=x+εω(i)

(2)

(3)

(4)

為解決EEMD存在的重構誤差，YEH等提出了CEEMD方法。噪聲被成對(一正一負)的添加到原始數據中進行兩次集總分解。

(5)

(6)

步驟1 在每次集總數i=1,…,I對信號x(i)=x+β0ω(i)進行EMD分解，得到第一階模態(tài)。

(7)

步驟2 利用式(6)計算第一階(k=1)殘差r1。

步驟3 給第一階殘差添加噪聲β1E1(ω(i))，得到分解輸入信號r1+β1E1(ω(i))，(i=1,…,I)，再進行EMD分解得到該信號的第一階模態(tài)也就是原信號的第二階模態(tài)。

(8)

步驟4 對于k=2,…,K的k階殘差為

(9)

步驟5 計算信號rk+βkEk(ω(i))，(i=1,…,I)，的第一階模態(tài)得到原始信號的第k+1階模態(tài)。

(10)

步驟6 回到步驟4計算下一第k階模態(tài)。

該方法停止規(guī)則為殘差不能再進行EMD分解。

最終趨勢項為

(11)

其重構公式為

(12)

該方法由于每一階分解的初始信號均相同 (rk)，并且只需要得到一階模態(tài)，所以確保了不同集總數下均有相同數量的模態(tài)(模態(tài)數為1)，并且該方法具有完美的重構誤差。

但是該方法在處理一些信號時會出現同一模態(tài)的偽分量出現在固有模式前后，并且分解模態(tài)中依然包含殘留噪聲。基于這些問題文獻[7]又對CEEMDAN方法做了進一步的改進。從原始CEEMDAN中可以發(fā)現，首先是第一階模態(tài)求解時直接向信號中添加了高斯白噪聲，而第二階求解中卻添加的是經過EMD分解后的噪聲，這使得每一階求解中出現噪聲的交叉干擾，改進方法將第一階模態(tài)求解時信號添加噪聲修改為一階EMD分解得到的噪聲，其它階模態(tài)求解添加的噪聲是相應階次EMD分解后的噪聲。其次原方法每次求解模態(tài)分量時均對模態(tài)函數求平均以消除噪聲的影響，但是每次添加的噪聲頻帶均和分解出的固有模態(tài)頻帶重合，實際添加的噪聲被大量分解在了模態(tài)函數中。然而，信號分解殘差的頻帶小于分解出的模態(tài)函數頻帶，相應的也小于添加噪聲的頻帶。因此分解到殘差中的噪聲最小，最容易被平均，所以改進方法是將原來求解模態(tài)函數均值修改為求解殘差均值，而模態(tài)函數有分解信號減去殘差均值得到。定義算子M(·)為計算EMD第一階模態(tài)后的殘差信號。改進后的方法具體算法步驟為

步驟1 計算信號x(i)=x+β0E1(ω(i))(i=1,…,I)第一階殘差信號

(13)

步驟3 第一階殘差中添加噪聲β1E2(ω(i))，按照步驟1和步驟2求解第二階模態(tài)

(14)

步驟4 對于k=3,…,K的k階殘差為

(15)

步驟6 回到步驟4計算下一階模態(tài)。

其分解的每一階噪聲信噪比均能選擇βk=εkstd(rk)。

2 完全小波噪聲輔助集總經驗模態(tài)分解和完全互補小波噪聲輔助集總經驗模態(tài)分解

詳細分析CEEMDAN算法發(fā)現，在每階殘差信號中均添加特定的噪聲信號βkEk(ω(i))，該噪聲是通過EMD分解得到，那么CEEMDAN方法實際的迭代計算次數Nreal為

Nreal=Nx+Nnoise

(16)

式中:Nx為計算各階模態(tài)的迭代次數;Nnoise為計算各階噪聲的迭代次數。利用CEEMDAN分析如下仿真信號x。

(17)

CEEMDAN方法分析該仿真信號，計算實際模態(tài)的迭代次數箱型圖為圖1，計算各階噪聲的迭代次數箱型圖為圖2。圖中可以看出計算各階噪聲花費的計算量也非常巨大。所以該方法的實際迭代次數為前兩者迭代次數的總和，其箱型圖為圖3所示。

圖1 計算各階模態(tài)的迭代次數Fig.1 The iteration number of each mode

圖2 計算各階噪聲的迭代次數Fig.2 The iteration number of each noise

圖3 完整CEEMDAN的實際迭代次數Fig.3 The actual iteration number of complete CEEMDAN

在本組分解中計算各階模態(tài)函數花費的總迭代次數為Nx=19 735次，計算各階噪聲的總迭代次數為Nnoise=7 376次，該方法實際總迭代次數為Nreal=27 111次，噪聲迭代次數約占比實際總迭代次數的27%，相應的計算各階模態(tài)函數的迭代次數約占比實際總迭代次數的73%，并且分解層數越少則花費在噪聲分解中的迭代次數占比會急劇上升。

單獨分析EMD方法對高斯白噪聲的分解性質，設定20 000個獨立的高斯白噪聲數據(其中均值為0，方差為1)，對每個高斯白噪聲數據進行EMD分解得到不同的IMF分量，計算每個分量的功率譜。計算20 000個信號對應功率譜的平均值并進行幅值歸一化，計算結果為EMD的等效濾波特性。具體如圖 4所示，其橫坐標化為了角頻率并做了歸一化。

圖4 EMD的等效濾波特性Fig.4The equivalent filter characteristics of EMD

分析圖4可得出，白噪聲被EMD規(guī)律性的分解為能量均勻分布的各頻率分量，EMD表現出有效的二元濾波器組特性。而小波分解也具有相似的特征，由Mallat算法可知，二進小波分解對信號頻帶劃分也表現為二元濾波器組特性。本文用coif3小波分解高斯白噪聲得到小波分解的等效濾波特性如圖5所示。對比圖4和圖5發(fā)現coif3小波和EMD方法對白噪聲的頻帶劃分同樣表現為二元濾波器組特性。

圖5 小波的等效濾波特性Fig.5 The equivalent filter characteristics wavelet

正是基于小波分解和EMD分解均表現為相同濾波特性的特點，可將CEEMDAN分解中求解各階添加噪聲的EMD算法改為小波分解代替，以達到在繼承CEEMDAN優(yōu)良特性的基礎上減少算法迭代次數提高計算效率的目的。新方法為完全小波噪聲輔助集總經驗模態(tài)分解，定義算子Wj(·)是通過小波計算給定信號j層小波分解并單子重構第j個細節(jié)信號的算子，算子M(·)為用EMD求解一階模態(tài)后的殘差信號的算子，同樣高斯白噪聲為ω(i)(i=1,…,I)，I為集總次數。CWEEMDAN方法詳細計算步驟流程為

步驟1 在信號x添加I個噪聲，即x(i)=x+β0W1(ω(i))(i=1,…,I)，對每組信號進行一層EMD分解得到I個殘差，平均后得到第一階殘差分量。

(18)

步驟3 第一階殘差r1疊加I個噪聲信號β1W2(ω(i))，并計算第二階殘差。

(19)

步驟5 對于k=3,…,K的k階殘差為

(20)

步驟7 回到步驟5計算下一階模態(tài)。

此外，詳細分析上訴算法發(fā)現CWEEMDAN還可以根據CEEMD方法中減小添加噪聲影響的方法，將CWEEMDAN中每次添加的特定的噪聲信號βkEk(ω(i))正負成對的添加到分解信號中以達到進一步的減小添加噪聲的影響，最終達到在求解相同分解結果時，減少集總次數。所以最終改進算法為完全互補小波噪聲輔助集總經驗模態(tài)分解,算法的具體步驟如下(見圖6)。

圖6 CCWEEMDAN算法流程圖Fig.6 Flowchart describing the CCWEEMDAN

同樣高斯白噪聲為ω(i)(i=1,…,I/2)，I(偶數)為集總次數。

步驟1 在信號x添加I/2對噪聲，對每組信號進行一層EMD分解得到I個殘差，平均后得到第一階殘差分量

(21)

(22)

(23)

(24)

步驟5 對于k=3,…,K的k階殘差為

(25)

(26)

步驟7 回到步驟5計算下一階模態(tài)。

3 EMD, EEMD, CEEMD ,CEEMDAN，CWEEMDAN and CCWEEMDAN分解性能比較

為對比EMD，EEMD，CEEMD，CEEMDAN CWEEMDAN和CCWEEMDAN方法分解性能，構造仿真信號:s=s1+s2。

(27)

(28)

該仿真信號及各組成成分的時域波形圖如圖7所示。用以上涉及的六種方法對該組信號進行分解，噪聲輔助EMD分解方法的集總數為50，噪聲水平為0.2，分解結果如圖8所示。

圖7 仿真信號及其各組成成分波形圖Fig.7 The waveforms of simulation signal and its components

圖8 仿真信號EMD,EEMD,CEEMD,CEEMDAN,CWEEMDAN和CCWEEMDAN分解結果Fig.8 Decomposition of artificial signals by EMD，EEMD，CEEMD，CEEMDAN，CWEEMDAN and CCWEEMDAN

對比圖8的結果發(fā)現，EMD分解存在嚴重的模態(tài)混疊，未能將間歇信號s1和諧波信號s2完整分解開。而剩下的四種噪聲輔助方法均將兩個分量信號完整分解開來，避免了模態(tài)混疊。但是EMD，EEMD和CEEMD方法分解出了9個或以上的模態(tài)，雖然從第三個模態(tài)開始具有很小的能量，但是并不代表原始信號的信息。而CEEMDAN CWEEMDAN和CCWEEMDAN方法可以在分解的每一階驗證全局的停止準則，一旦滿足IMF條件，整體分解很快停止。噪聲輔助的EMD分解具有一定隨機性，對于同一組信號同樣兩次的分解結果稍有不同。此外，隨著集總數的增加，分解的殘余噪聲能量應當進一步減小，為量化該性能指標，在集總數I=50，100，200，400，800時用上述方法中的五種噪聲輔助方法(EEMD，CEEMD，CEEMDAN，CWEEMDAN，CCWEEMDAN)對仿真信號進行100次分解，其中噪聲幅度ε=0.2。

圖9 五種噪聲輔助EMD方法對仿真信號s的分解性能Fig.9 Performances of the five noise-assisted methods on artificial signal s

為比較分析方法對原始信號各分量的恢復能力，定義恢復信號a和基準信號b的相對平方根誤差為

(29)

對已知分量信號的恢復性能分析(見圖10和圖11)，CWEEMDAN比EEMD，CEEMD，CWEEMDAN的結果均要好，而CWEEMDAN在集總次數<100時的結果也優(yōu)與CEEMDAN。當集總次數>200時CEEMD，CEEMDAN，CWEEMDAN和CCWEEMDAN的恢復性能不隨集總次數的提高而提高，而CEEMDAN方法在集總次數大于200時的恢復性能略優(yōu)于CCWEEMDAN，這是因為在集總次數足夠大時大部分噪聲幾乎被平均，而由小波濾波特性與EMD濾波特性的差異引起的噪聲未被完全平均。

對比CEEMDAN和CCWEEMDAN的恢復性能可以發(fā)現，雖然在集總次數非常大時，CEEMDAN具有比CCWEEMDAN較好的恢復性能，但是其相差不大，幾乎不會對結果有影響，并且CCWEEMDAN具有在低集總次數得到高恢復性能的優(yōu)點。最后，對五種方法的信號重構性能做對比分析(見圖12)，CCWEEMDAN重構性能均優(yōu)于EEMD，CEEMD和CWEEMDAN，并且重構誤差和CEEMDAN達到了同一數量級10-17。CCWEEMDAN和CEEMD比較發(fā)現隨著集總次數的增加重構誤差變大，其原因是加入不同噪聲會使信號產生不同數量的模態(tài)，隨著集總次數的增加產生不同數量模態(tài)的概率增加，所以重建誤差增大。

五種噪聲輔助方法對仿真信號每次EMD分解的平均迭代次數如圖13所示。比較四副圖發(fā)現，CCWEEMDN方法比其它四種方法具有更少的迭代次數，更低的計算消耗。為進一步比較五種噪聲輔助EMD方法的計算開銷，使用一組真實數據對五種方法進行驗證，該數據來源于MIT-BIH正常竇性心律的心電圖(ECG)信號。對該組心電信號進行五種噪聲輔助EMD方法的分解(集總數為100，噪聲幅度為0.2)，其中EEMD、CEEMD、CEEMDAN和CWEEMDAN方法篩選各階模態(tài)需要的迭代次數如圖14所示，而本文最終提出的CCWEEMDAN方法篩選各階模態(tài)的迭代次數如圖15所示。比較圖14和圖15可知，EEMD和CEEMD方法

圖10 五種噪聲輔助EMD方法對仿真信號s1的恢復性能Fig.10 Performances of the five noise-assisted methods on artificial signal s1

計算各階模態(tài)的迭代次數最多，并且這兩種方法得到的模態(tài)也最多；CEEMDAN方法和CWEEMDAN、CCWEEMDAN方法得到的模態(tài)數量均較少，但是CEEMDAN方法篩選各階模態(tài)的迭代次數均大于后兩種新的算法。五種噪聲輔助EMD方法單次EMD分解的迭代次數如圖16所示。由圖16可知，迭代次數最多的EEMD和CEEMD，它們單次EMD分解最大概率的迭代次數分別為543次和529次，CEEMDAN單次EMD分解最大概率的迭代次數為356次，迭代次數最少的是CWEEMDAN和CCWEEMDAN，其迭代次數分別為238次和240次。通過計算可知，CEEMDAN方法計算效率大約是EEMD和CEEMD的1.5倍；CWEEMDAN和CCWEEMDAN計算效率大約是CEEMDAN的1.5倍；總體而言新噪聲輔助EMD方法的計算效率大約是原始EEMD方法的2.25倍。如果考慮到單次EMD分解中的異常迭代次數，集總次數為100次的EEMD迭代次數為65 729次，CEEMD總迭代次數為62 918，CEEMDAN總迭代次數為39 937次，CWEEMDAN和CCWEEMDAN總迭代次數分別為25 594次和26 978次(見圖 17)。分析圖17可知，CCWEEMDAN總迭代次數僅僅是CEEMDAN的67.55%，是CEEMD的42.84%，是EEMD的41.04%，幾乎和CWEEMDAN持平。

圖12 五種噪聲輔助EMD方法重構性能Fig.12 Performances of the five noise-assisted methods on artificial signal s: reconstruction error

圖13 五種噪聲輔助EMD方法平均迭代次數Fig.13 The average number of iterations of five noise-assisted methods

圖14 EEMD，CEEMD，CEEMDAN和CWEEMDAN刷選各階模態(tài)的迭代次數Fig.14 The iteration number of each mode by EEMD, CEEMD, CEEMDAN and CWEEMDAN

圖15 CCWEEMDAN篩選各階模態(tài)的迭代次數Fig.15 The iteration number of each mode by CCWEEMDAN

圖16 五種噪聲輔助EMD方法單次EMD分解的迭代次數Fig.16 The sifting iterations for each mode of the five noise-assisted methods

圖17 五種噪聲輔助EMD方法總迭代次數Fig.17 The total iterations of five noise-assisted methods

五種噪聲輔助EMD方法的重構誤差如圖18所示，除EEMD外其余四種噪聲輔助EMD方法的重構誤差均小于5×10-13。比較后四種方法的重構誤差平方和(如圖19)，CCWEEMDAN的重構誤差平方和最小。

圖18 五種噪聲輔助EMD方法重構誤差Fig.18 The reconstruction error of five noise-assisted methods

圖19 四種噪聲輔助EMD方法重構誤差平方和Fig.19 The quadratic sum of reconstruction error of four noise-assisted methods

4 結 論

本文提出了非線性、非平穩(wěn)信號處理新算法，新算法對仿真信號和實測信號具有很好的測試效果。新方法(CCWEEMDAN)概括起來具有以下特點：

(1)CCWEEMDAN同EEMD一樣具有抗頻率混疊的特性。

(2)CCWEEMDAN比EEMD和CEEMD具有更少的偽模態(tài)數量。

(3)CCWEEMDAN比EEMD和CEEMD有更小的殘余噪聲和重構誤差，以較小集總次數獲得比CEEMDAN更好的分解效果和重構精度。

(4)CCWEEMDAN較其它噪聲輔助方法具有更小的迭代次數和計算開銷。

[1] HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for non-linear and non-stationary time series analysis [J]. Proceedings of the Royal Society a Mathematical Physical & Engineering Sciences, 1998,454:903-995.

[2] HUANG N E，SHEN Z，LONG R S． A new view of nonlinear water waves the Hilbert spectrum [J].Annual Review of Fluid Mechanics，1999，31(1): 417-457．

[3] WU Z, HUANG N E. Ensemble empirical mode decomposition: a noise-assisteddata analysis method [J]. Advances in Adaptive Data Analysis, 2009,1 (1):1-41.

[4] FLANDRIN P, RILLING G, GONCALVS. Empirical mode decomposition as a filter bank [J]. IEEE Signal Processing,Letters,2004,11(2):112-114.

[5] YEH J R, SHIEH J S, HUANG N E, Complementary ensemble empirical modedecomposition: a novel noise enhanced data analysis method [J]. Advances in Adaptive Data Analysis,2010,2 (2):135-156.

[6] TORRES M E, COLOMINAS M A, SCHLOTTHAUER G, et al. A complete ensemble empirical mode decomposition with adaptive noise [J]. Brain Research Bulletin, 2011,125(3):4144-4147.

[7] COLOMINAS M A, SCHLOTTHAUER G, TORRES M E. Improved complete ensemble EMD: a suitable tool for biomedical signal processing [J]. Biomedical Signal Processing and Control,2014,14(1):19-29.

A complete complementary wavelet ensemble empirical mode decomposition with adaptive noise

HELiu,DINGJianming,LINJianhui,LIUXinchang

(State Key Laboratory of Traction Power, Southwest Jiaotong University, Chengdu 610031, China)

Empirical mode decomposition (EMD) is a self-adaptive method and suitable to analysing the non-stationary and nonlinear signals. Noise-assisted versions have been proposed to alleviate the so-called “mode mixing” phenomenon, which may appear when an EMD algorithm is used to deal with a signal with intermittency. Among them, the complementary ensemble EMD (CEEMD) and complete ensemble EMD with adaptive noise (CEEMDAN) recover the completeness property of EMD. In this work a new algorithm named complete complementary wavelet ensemble empirical mode decomposition with adaptive noise(CCWEEMDAN) was presented based on those existing techniques, obtaining better spectral separation of the modes with fewer sifting iterations and less noise of components with small ensemble number and extremely low computational cost.

empirical mode dercompostion; ensemble empirical mode decomposition; noise-assisted data analysis; mode mixing; complementary ensemble empirical mode decomposition

國家自然科學基金項目(61134002; 51305358)

2015-09-09 修改稿收到日期：2016-01-04

何劉 男，碩士生，1990年生

丁建明 男，博士，助理研究員，1981年生

U211；U270

A

10.13465/j.cnki.jvs.2017.04.037