基于二階包絡逆EMD算法改進與實現

何經偉， 胡維平， 莫家玲

(廣西師范大學 電子工程學院 廣西多源信息挖掘與安全重點實驗室，廣西 桂林 541004)

基于二階包絡逆EMD算法改進與實現

何經偉， 胡維平， 莫家玲

(廣西師范大學 電子工程學院 廣西多源信息挖掘與安全重點實驗室，廣西 桂林 541004)

在分析經典EMD(Empirical Mode Decomposition)和高階極值點逆向篩分EMD方法的基礎上，提出了一種基于二階包絡逆EMD的改進方法。該方法結合標準EMD與逆EMD的良好屬性，采用一階、二階包絡，通過給定條件在篩分過程選擇最優包絡均值進行篩分。實驗結果表明：該方法適用于分解頻率差值較低的多分量信號且減小高頻分量誤差；有效地減小估計誤差及抑制多余IMF分量；以最優包絡作為篩分包絡，EMD分解效果更佳。

標準EMD；逆EMD；二階包絡逆EMD

Hilbert-Huang變換是 HUANG等[1]提出的一種用于非平穩、非線性信號的自適應處理方法。它突破了傳統的基于短時傅里葉變換的信號處理方法，是信號處理方法上的一次飛躍。HHT主要有兩個過程：經驗模態分解(EMD)和Hilbert變換。EMD目的是為了得到IMF，以便進一步利用Hilbert變換進行分析。在EMD算法中，信號的極值點信息在分解過程中至關重要。極值點信息的有效利用及EMD分解算法的改進研究，一直吸引著學者們。徐冠雷等[2]通過分析研究二級極值點結構，提出了從多分量到單分量用EMD分解的極值點判斷準則。 KOPSINIS等[3]提出了一種基于差值點選擇的EMD分解方法。 KOPSINIS 等[4]提出了一種基于特定最優極值點的啟發式搜索的改進EMD分解法。XU等[5]提出了限鄰域均值的EMD分解法。LI等[6]通過利用高階極值點信息，提出了一種由最優包絡均值逆向篩分的EMD分解法。

本文首先介紹兩種包絡方式的EMD分解方法，分析其存在的優勢及局限性，然后提出一種新的EMD方法并通過大量實驗驗證該方法的有效性。

1 經典EMD算法的基本原理

對給定信號x(t)的EMD分解算法[7]具體流程如下：

(1) 首先確定給定信號x(t)所有的局部極大值點以及極小值點，通過利用三次樣條插值法通過極大值點(極小值點)分別擬合出上u(t)(下1(t))包絡曲線。計算上下包絡均值m(t)=(u(t)+1(t))/2，將原始信號x(t)減去m(t)得h(t)，即：h(t)=x(t)-m(t)。

(2) 判斷h1(t)是否滿足IMF的兩個條件，如果滿足，h1(t)就是篩選出來的第一個本征模態函數分量，如果不滿足，把h1(t)作為原始數據繼續重復進行上述步驟一操作，直到h1(t)滿足IMF分量的條件。這樣經過多次篩選我們將得到第一個IMF分量C(t)。

(3) 將C(t)從原始信號x(t)分離出來得到剩余信號：r(t)=x(t)-C(t)

(4) 我們將剩余信號繼續進行上述一至三步驟操作，這樣得到第二個本征模態函數分量。如此重復K次，就會得到滿足IMF條件的K個IMF分量，以及一個最終殘余分量r(t)。殘余分量是一個不為零的分量，其代表信號的趨勢。此時滿足EMD終止準則，篩選結束。我們將信號分解結果可表示如下：

(1)

式中x(t)為原始信號，n表示IMF分量的個數，cj(t)為第j個IMF分量rn(t)為殘余分量。

經典EMD通過采用三次樣條插值擬合一階包絡，通過多次篩分的方式將復雜信號分解成各個本征模態函數。該方法具有普遍適用性，不管是復雜實際信號，還是多分量信號，都具有較好的分解效果，但是該方法對高階包絡沒有一定預測能力，會產生較大的估計誤差，在分解過程中，不能有效避免多余的IMF分量。

2 逆EMD算法的基本原理

逆EMD分解中的“逆”并非完全意義上的逆，而是因為它是一種基于高階極值點的EMD方法。它利用n階極值點信息構造成n階最優包絡通過逆向篩分的方式獲得各個本征模態函數。

原逆EMD算法[8]流程：

(1) 輸入給定信號s(t)，確定極值點階數n。

(2) 分別確定n階極大值點un(tk)，極小值點ln(tk)，極小極大值點lun(tk)，極大極小值點uln(tk)。

(3) 分別通過三次樣條插值構造上包絡線un(t)、下包絡線ln(t)。以及極小上包絡線lun(t)、極大下包絡線uln(t)。

(4) 計算出三個不同的包絡均值m1(t)=(un(t)+ln(t))/2、m2(t)=(uln(t)+lun(t))/2、m3(t)=(un(t)+ln(t)+uln(t)+lun(t))/4判斷包絡均值是否滿足給定條件OC，滿足的話，從三個均值中選擇最優均值執行步聚(5)，否則，跳轉至步驟(6)。

(5) 讓c1(t)+c1(t)+m(t)和s(t)=s(t)-m(t)重復以上操作直到m(t)不再滿足給定條件。這時剩余信號是最高頻率IMF1。

(6) 讓n=n-1，重復步驟(1)～(5)，直至n=1，得到剩余信號r(t)即為最高頻IMF1。

(7) 讓s(t)=s(t)-r(t)，重復步驟(1)～(6)，得到其他的IMF，直至所有IMF之和滿足給定條件OC1，則EMD循環終止，得到一個非本征模態函數的殘余信號。

3 逆EMD存在的問題以及改進的方法

逆EMD與標準EMD方法的最大不同之處在于最優包絡均值的擬合過程。通常，最優包絡均值是復雜的多分量、非平穩信號。標準EMD方法通過多次三次樣條插值來擬合，而逆EMD方法考慮到信號多分量的特點，采用不同階極值點來擬合該最優包絡均值。該方法的優勢是：可以有效的減小估計誤差，同時有效的抑制多余IMF分量。但是逆EMD方法中，由于不同階極值點是根據多分量信號特點提出的，所以不適用于復雜實際信號尤其是信號分量頻率差別較低的情況。逆EMD方法采用n階包絡深入式的逆向篩選，過分依賴不同階極值點擬合最優包絡均值，導致低階包絡均值篩分的強度不夠，從而得到的高頻分量誤差較大。所以本文提出了一種二階包絡逆向篩分與多次使用三次樣條插值擬合相結合的逆EMD改進方法。既保證了對信號分解時有效的估計，同時不會因高階極值點問題局限逆EMD僅適用于分解一些高頻率差的多分量信號，同時減小了高頻分量誤差。

一個IMF必須滿足以下兩個條件[9]：①在整個數據長度，極值和過零點的數目必須相等或最多相差一個；②在任意數據點，局部最大值的包絡和局部最小值的包絡的平均必須為0。在逆EMD分解過程中發現，當篩選進行一定程度時，其分解的信號已經滿足本征模態函數的基本條件，但由于上、下擬合的包絡有錯差問題致使其上、下包絡不對稱，導致篩選過程不能很快的收斂并且結束。這樣的情況如圖1所示，該測試信號完全滿足本征模態函數條件，其過零點次數與極值點個數相差為零，可是由于原始信號極小值點位置情況，導致將下包絡錯誤的擬合成了一條近似直線，使得包絡均值不能為零，如圖1的淺灰色曲線(原始的上、下包絡線)所示。所以為了避免類似的情況發生，我們采用將極大值點翻轉變為極小值點與信號極小值點一起擬合成下包絡，如圖1的深灰色曲線(新的下包絡線)所示。這樣它們的包絡均值為零。圖1中淺灰色曲線為原始包絡線，深灰色為修正后包絡線。

圖1 測試信號Fig.1 Test signal

為了避免類似情況的發生，在逆EMD算法中做了如下改進：將極大值點翻轉作為信號的極小值點來重新擬合新的下包絡，同樣方法將極小值點翻轉作為信號的極大值點來重新擬合成新的上包絡。然后通過條件判斷哪種上下包絡更符合要求(即包絡均值更接近零均值)而考慮選擇哪種包絡。該方法通過改變包絡擬合方式，從而改變了包絡均值，而IMF的正交性停止準則主要受包絡均值m(t)與輸入信號s(t)的相關性(m(t)=b·s(t)，b為相關系數)所影響，所以篩選過程中，采用該方法包絡均值快速收斂接近零均值，m(t)與s(t)相關性減小，同時正交性OC逼近停止條件，進而不致于過分的篩分出無物理意義的其他分量。由于整個篩分過程采用的是由OC選擇最優包絡均值的方式，選擇使得OC最大的m(t)作為最優包絡，也就是在篩分時篩分掉的是與信號s(t)相關性最大的包絡均值，所以一定程度上保證了分解得到IMF分量的正交性。

4 逆EMD的改進算法

改進逆EMD的算法流程如下：

(1) 設置初始值n=1，輸入給定信號s(t)。

(2) 分別確定信號極大值點U(t)，信號極小值點L(t)，將極小值點翻折構成新的極大值點集U1(t)，將極大值點翻折構成新的極小值點集L1(t)，利用三次樣條插值構造原始(上下)包絡以及翻折后的(上下)包絡，通過條件判斷原始包絡和翻折包絡來選擇作為一階包絡m1(t)。

(3) 分別確定信號中極小的極大值點LU(t)，極大的極小值點UL(t)，利用三次樣條插值構造上、下包絡，獲得二階包絡m2(t)。

(4) 分別獲得三種包絡均值：一階包絡均值m1(t)、二階包絡均值m2(t)、一階包絡與二階包絡的平均值：m3(t)=(m1(t)+m2(t))/2，判斷三種包絡均值是否滿足給定條件OC>1.05，滿足的話，從三種包絡均值中選擇最優均值(也就是從三種包絡均值中選擇使得OC最大的m(t)作為最優均值)執行步驟(5)，否則，跳轉至步驟(6)。

(5) 讓c1(t)=c1(t)+m(t) 和s(t)=s(t)-m(t)重復以上操作直到m(t)不再滿足給定條件。這時剩余信號是最高頻率IMF1。

(6) 讓n=n+1，s(t)=s(t)-IMF(t)，重復步驟(1)～(5)，得到其他的IMF，直至所有IMF之和滿足給定條件OC1，則EMD循環終止，得到一個非本征模態函數的殘余信號。

信號經逆EMD分解后可用如下式子表達：

S(t)=IMF1+IMF2+IMFn+rn(t)

(2)

n表示EMD分解所得的IMF個數。

改進前逆EMD分解效果如圖2、3所示；改進后逆EMD分解效果如圖4、5所示。

圖2 改進前逆EMD分解所得IMF1Fig.2 IMF1 of inverse EMD

圖3 改進前逆EMD分解所得IMF2Fig.3 IMF2 of inverse EMD

圖4 改進后逆EMD分解所得IMF1Fig.4 IMF1 of improved inverse EMD

圖5 改進后逆EMD分解所得IMF2Fig.5 IMF2 of improved inverse EMD

根據本征模態函數定義，其實上述測試信號已經滿足本征模態函數條件。可是當用改進前逆EMD算法來進行分解時，會得到兩個信號成分IMF1和IMF2，通過圖2、3與圖4、5比較可知，原逆EMD比改進后逆EMD篩分的程度更強，IMF信號的帶寬更窄。這表明，在原逆EMD的篩選過程，容易導致過分篩選現象。也就是在本來滿足條件的本征模態函數信號中，如果繼續不斷的篩選，從主幅度的IMF信號中會篩分出一個小幅度的IMF信號，即過分的篩分出無物理意義的其他分量。而當經改進后逆EMD算法分解時，其幾乎是不能再次分解的，且其得到的IMF2分量幾乎不含任何信號成分，信號幅度一直在零值附近。如圖4、5所示。由于包絡沒有嚴格的數學定義，采用這樣的一種方式構造新包絡，雖然一定程度上增加了IMF的帶寬，但可以有效的避免EMD分解過程中的IMF過篩選現象，可以保證IMF的信號完整性，而不會過分的篩分出其他分量。

實驗1 正交性測試

為了驗證本文方法的有效性，本次測試采用式(1)的正交性( Index of Orthogonality, IO) 指標來比較經典EMD和逆EMD分解的性能。

(3)

IO指標越接近于0，說明EMD分解得到的各個內蘊模式函數之間越接近于正交。

基于P.Flandrin數值實驗結論[10]：

當兩信號頻率比a<0.5(a=f2/f1，f20.7，兩信號成分完全不能分離到兩個IMF中。

我們構造信號針對于完全分離、部分分離、完全不能分離三種情況，分別測試正交性指標來比較經典EMD與逆EMD分解的性能。

構造完全分離時的仿真信號：

y1=cos((α·f0(i))·2·pi·t)

y2=cos((β·f0(i))·2·pi·t)

t∈(1,1 024)

y=y1+y2

其中f0在頻率0至0.5范圍內，頻域等距的500個仿真信號，完全分離時設置α=1，β=0.4；部分分離時設置α=1，β=0.6；完全不能分離時設置α=1，β=0.7。

從圖6～圖8可以看出，當兩種信號成分無論是完全分離、部分分離、還是完全不能分離，通過逆EMD分解方法得到兩成分IMF1和IMF2，其正交性都比標準EMD效果要好。而且隨著頻率差值的逐漸接近，二階包絡逆EMD體現出更好的分解效果。換言之，該方法比經典EMD、逆EMD更適用于分解一些頻率差值較低的多分量信號。

圖6 完全分離Fig.6 Completely separate

圖7 部分分離Fig.7 Partly separate

圖8 完全不能分離Fig.8 Completely unable to separate

實驗2 能量泄漏測試

在EMD算法中，每個IMF分量都是經過多次篩選的過程來獲得，而每次篩選的過程都是通過信號的局部極大、極小值點用三次樣條插值的方法擬合而成的上、下包絡線的均值。但是由于采樣間隔的存在，信號的局部極值點很有可能定位失配，從而會產生極值點定位的偏離誤差，最終影響會體現在EMD分解中，我們稱之為能量泄漏。在理想的情況下，如果原始信號是單頻率的正余弦信號，對其進行EMD分解，由于信號的上下包絡線是完全對稱的，所以分解的結果就只會有一個IMF分量即為此正余弦信號。而實際上，將信號進行EMD分解時，由于采樣頻率的緣故，離散的單頻率信號其上下包絡線并不能完全的對稱，最終導致會有多個IMF的結果，我們就說，采樣引起的非對準使得EMD分解的部分能量泄漏到其他的IMF中了。

對不同采樣頻率的單位幅度單頻率正弦信號進行EMD分解，定義由于采樣的非對準而造成其能量泄漏公式為：

(4)

式中N為數據段長度，r為EMD分解的層數，求和到r+1意味著包括余量。End_error指標越接近于0，說明EMD分解越接近無損。

測試信號：

y=cos((1.0·f0(i)·2·pi·t)

t∈(1,1 024)

其中f0在頻域0至0.5范圍內，頻率等距選擇500個頻率點用來構造仿真信號。

圖9 EMD分解與逆EMD分解能量泄漏的比較Fig.9 Energy leakage contrast of EMD and inverse EMD

由圖9可知三種EMD方法的能量泄漏均滿足非對準誤差曲線。從整個頻域角度來看，在一些頻率點位置，逆EMD分解的能量泄漏明顯低于標準的EMD分解能量泄漏值。當用SUM_EMD來求三種方法在全頻段的能量泄漏占用總面積時，逆EMD能量泄漏總值小于標準EMD，即在全頻段，逆EMD效果更好。改進前、后逆EMD方法的能量泄漏曲線很接近，但是改進后的逆EMD能量泄漏普遍小于改進前逆EMD。通過對三種EMD方法能量泄漏測試，實驗證明，改進的逆EMD方法效果最好。

前面進行了三種EMD方法正交性以及能量泄漏的測試與比較，下面探討在逆EMD分解過程中分別采用三種包絡曲線去實現信號的篩選，并測試其分解后的能量泄漏。第一種：在逆EMD分解過程中，整個過程使用的是一階包絡曲線求取包絡均值(m1)去完成信號的篩選；第二種：整個過程采用的是一階與二階包絡曲線共同求取包絡均值(m3)去完成信號的篩選；第三種：整個過程將一階包絡均值、二階包絡均值以及一階與二階兩個包絡的平均值通過一個準則從中選擇一個作為最優包絡均值(option)去完成信號的篩選。

實驗3 針對這三種方式，測試其在逆EMD分解過程中能量泄漏指標

圖10 三種包絡方式的能量泄漏Fig.10 Energy leakage of three kinds of envelope

圖10中深灰色為最優包絡方法，淺灰色為一階包絡方法，處于深灰色與淺灰色曲線之間的曲線為二階包絡方法。這三種包絡曲線方法在逆EMD分解過程中的能量泄漏均滿足非對準誤差曲線。三種方法的能量泄漏值比較接近，但是最優包絡在逆EMD中效果更佳。

實驗4 改進的逆EMD在信號處理中的應用

心電圖信號是一種微弱的生物電信號，也是一種典型的非平穩隨機信號。該方法可以有效地減小高頻分量誤差，抑制多余IMF分量等諸多優勢。在此，將改進后的逆EMD方法應用到心電圖信號的分解中，并把處理結果同改進前的EMD分解方法作比較，以對比處理效果。

本次測試采用的心電圖信號來源于近年來應用比較廣泛的國際上公認的美國麻省理工學院提供的MIT-BIH數據庫。圖11和圖12中最上面的信號s(t)為一段心電圖信號，含有高頻噪聲干擾和電平整體偏置；最下面的信號res.為剩余信號，主要為整體偏置電平；c1(t)含有較多的高頻噪聲；c2(t)主要包含了心電信號。 改進后EMD得到的c2(t)與改進前EMD的相比高頻噪聲幅度明顯較小，可以忽略不計，說明改進后EMD方法可以有效地減小高頻分量誤差。改進前分解得6個IMF分量，改進后分解得4個IMF分量，說明改進后EMD方法可以有效抑制多余IMF分量。

圖11 改進前逆EMD分解結果Fig.11 Result of inverse EMD

圖12 改進后逆EMD分解結果Fig.12 Result of improved inverse EMD

4 結 論

本文提出基于二階包絡逆向篩分的EMD方法，該方法有效的結合經典EMD與高階逆EMD的優勢。數值實驗證明該方法不僅適用于多分量信號，而且適用于復雜實際信號。對于分解一些頻率差值較低的多分量信號，可以減小高頻分量誤差。可以有效的減小估計誤差，抑制多余IMF分量。

[ 1 ] HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-s tationary time series analysis[J]. Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences. The Royal Society, 1998, 454(1971):

903-995.

[ 2 ] 徐冠雷, 王孝通, 徐曉剛, 等. 多分量到單分量可用 EMD 分解的條件及判據[J]. 自然科學進展, 2006, 16(10): 1356-1360. XU Guanlei, WANG Xiaotong, XU Xiaogang, et al. Conditions and criteria of EMD for multi-component to single component [J]. Progress in Natural Science, 2006, 16(10): 1356-1360.

[ 3 ] KOPSINIS Y, MCLAUGHLIN S. Enhanced empirical mode decomposition using a novel sifting-based interpolation points detection[C]. IEEE/SP Workshop on Statistical Signal Processing, 2007,53(3): 315-318.

[ 4 ] KOPSINIS Y, MCLAUGHLIN S. Investigation and performance enhancement of the empirical mode decomposition method based on a heuristic search optimization approach[J]. Signal Processing, IEEE Transactions on, 2008, 56(1): 1-13.

[ 5 ] XU Guanlei, WANG Xiaotong, XU Xiaogang. Neighborhood limited empirical mode decomposition and application in image processing[C]∥ Proceeding of the fourth International Conference on Image and Graphics. Washingdon, DC:ICIG, 2007.

[ 6 ] LI Lin, JI Hongbing. Signal feature extraction based on an improved EMD method[J]. Measurement, 2009, 42(5): 796-803.

[ 7 ] HUANG N E, SHEN Z, LONG S R. A new view of nonlinear water waves: The Hilbert Spectrum 1[J]. Annual review of fluid mechanics, 1999, 31(1): 417-457.

[ 8 ] 江莉, 李林, 董惠. 基于改進 EMD 方法的多分量信號分析[J]. 振動與沖擊, 2009, 28(4): 51-53. JIANG Li, LI Lin, DONG Hui. An improved empirical mode decomposition method for multicomponent signal represntation [J]. Journal of Vibration and Shock, 2009, 28(4): 51-53.

[ 9 ] 鄭天翔, 楊力華. 經驗模式分解算法的探討和改進[J]. 中山大學學報 (自然科學版), 2007, 46(1): 1-6. ZHENG Tianxiang, YANG Lihua. Discussion and improvement on empirical mode decomposition algorithm[J]. Journal of Zhongshan University (Natural Science Edition) , 2007, 46(1): 1-6.

[10] RILLING G, FLANDRIN P. One or two frequencies?The empirical mode decomposition answers[J]. Signal Processing, IEEE Transactions on, 2008, 56(1): 85-95.

Improvement and implementation of inverse EMD algorithm based on second order envelope

HE Jingwei， HU Weiping， MO Jialing

(Guangxi Key Lab of Multi-source Information Ming & Security, College of Electronic Engineering, Guangxi Normal University, Guilin 541004，China)

On the basis of analysing the classic empirical mode decomposition (EMD) method and the inverse EMD method in which the high-order extreme value point is used in reverse screening, a method based on second order envelope inverse EMD was proposed. The method integrates the good properties of the classic EMD and inverse EMD. Making use of the first and second order envelope and introducing certain conditions in the screening process, the mean value of an optimal envelope was chosen for screening.The experimental results show that the method is suitable for being applied to decompose multi-component signals with low frequency difference between components. It can reduce the error of high frequency components, effectively decrease the estimation error and restrain the excess IMF components.The EMD decomposition effect will be better when the optimal envelope is adopted in screening.

classic EMD; inverse EMD; inverse EMD based on second order envelope

國家自然科學基金資助項目(61062011；61362003)；廣西多源信息挖掘與安全重點實驗室基金項目(13-A-03-02)；廣西壯族自治區高等學校科學研究項目(KY2015YB034)

2015-07-14 修改稿收到日期： 2016-01-21

何經偉 男，碩士，1988年生

胡維平 男，博士，教授，1963年生

TN911.7

A

10.13465/j.cnki.jvs.2017.07.019