郭 迪，沈洋洋，尹 兵

(南京郵電大學(xué) 電子科學(xué)與工程學(xué)院，江蘇 南京 210003)

經(jīng)驗?zāi)Ｊ椒纸舛它c效應(yīng)抑制方法的研究

郭 迪，沈洋洋，尹 兵

(南京郵電大學(xué) 電子科學(xué)與工程學(xué)院，江蘇 南京 210003)

針對經(jīng)驗?zāi)Ｊ椒纸夥?Empirical Mode Decomposition，EMD)中的端點效應(yīng)嚴(yán)重影響算法精度的情況，為了減小端點效應(yīng)，文中提出一種新的先匹配后鏡像延拓方法。該方法借鑒鏡像對稱延拓與波形特征匹配延拓，將信號先進(jìn)行波形特征匹配延拓再進(jìn)行鏡像對稱延拓。利用仿真數(shù)據(jù)對先匹配后鏡像延拓的方法進(jìn)行了驗證，并與原有幾種延拓方法在不同評價指標(biāo)下進(jìn)行了對比。結(jié)果表明，先匹配再鏡像延拓后，信號的包絡(luò)線發(fā)生畸變最小，同時新方法分解得到的IMF精度較高，正交性好。先匹配后鏡像延拓方法提高了經(jīng)驗?zāi)Ｊ椒纸獾木龋芨行У匾种平?jīng)驗?zāi)Ｊ椒纸夥ㄖ械亩它c效應(yīng)。

經(jīng)驗?zāi)Ｊ椒纸?；延拓；端點效應(yīng)；鏡像法；波形特征匹配法

0 引 言

經(jīng)驗?zāi)Ｊ椒纸夥?Empirical Mode Decomposition，EMD)是由N.E.Huang[1]等提出的新型自適應(yīng)信號時頻處理方法，它把復(fù)雜信號分解為模態(tài)函數(shù)(Intrinsic Mode Function,IMF)之和，同時通過希爾伯特變換對相位微分求解得到瞬時頻率，賦予了瞬時頻率實際的物理意義。由于分解基于信號的局部特征時間尺度，使得分解具有自適應(yīng)性與高效性，特別適用于分析非線性、非平穩(wěn)信號。在EMD分解過程中，需要根據(jù)信號提供的極值點信息進(jìn)行三次樣條插值，以便畫出信號的上、下包絡(luò)線，信號的左右兩個端點可能不是極值點，這樣會使擬合時產(chǎn)生誤差。所得到的上、下包絡(luò)線在信號的邊界附近會產(chǎn)生大幅擺動，并且這種擬合誤差會隨著分解逐漸向信號內(nèi)部傳播，導(dǎo)致得到的IMF分量的波形嚴(yán)重失真。

文中在對EMD分解端點效應(yīng)深入分析的基礎(chǔ)上，提出一種波形匹配法和鏡像法相結(jié)合的延拓算法，并通過實驗仿真證明了其有效性。

1 EMD基本原理與端點效應(yīng)

1.1 EMD基本原理

EMD分解最基本的步驟是將非線性、非平穩(wěn)的信號分解為表征信號特征時間尺度的IMF，它滿足單分量信號的物理解釋。任何IMF必須滿足兩個條件：

(1)整個數(shù)據(jù)的極值點的個數(shù)和零點的個數(shù)相等或是相差一個。

(2)在任意時刻，由局部極值點形成的上、下包絡(luò)線的平均值為零，即上、下包絡(luò)線關(guān)于時間軸局部對稱。

EMD分解步驟如下：找出信號全部局部極大值點與局部極小值點的位置坐標(biāo)，然后用三次樣條插值法分別得到局部極大值包絡(luò)線與局部極小值包絡(luò)線，求上、下包絡(luò)線坐標(biāo)的平均值，得到一條新曲線(均值線)并記為m1,記其與原始信號x(t)的差值為h1。

h1=x(t)-m1

(1)

判斷h1是否滿足IMF的條件，如果滿足，那么h1就是x(t)的第一個IMF分量；如果不滿足，用h1代替原始數(shù)據(jù)，如果不滿足，則重復(fù)前述步驟直到h1k=h1(k-1)-m1k滿足IMF的條件。h1k為信號x(t)的第一個本征模態(tài)函數(shù)。記c1=h1k，c1為輸入信號中頻率最高的分量。再對余量r1=x(t)-c1進(jìn)行分解，得到c2，重復(fù)上述篩選過程直到：

(1)當(dāng)分量cn或殘余分量rn比預(yù)定值小時停止循環(huán)；

(2)當(dāng)rn變?yōu)閱握{(diào)函數(shù)，從中不能再篩選出IMF分量為止。

分解停止后，得到信號x(t)的全部本征模態(tài)函數(shù)。這樣就有：

(2)

以式(3)信號為例，采樣頻率取fs=2 048Hz。

(3)

按照如上步驟，對信號進(jìn)行經(jīng)驗?zāi)B(tài)分解，得到3個IMF分量和一個res分量。如圖1所示，IMF1至IMF3分別代表從原始信號分解出來的第1至第3個本征模態(tài)函數(shù)，其頻率由高到低排列；res代表殘余分量。

圖1 原始信號及各個IMF分量

然后求每個本征模態(tài)函數(shù)的希爾伯特譜，并將它們相疊加，得到原始信號的希爾伯特譜。如圖2所示，圖中的3條曲線從上到下分別為IMF1至IMF3的希爾伯特譜。

圖2 原始信號的希爾伯特譜

1.2 端點效應(yīng)

經(jīng)驗?zāi)B(tài)分解過程中，需要根據(jù)信號提供的極值點信息進(jìn)行三次樣條插值以便畫出信號的上、下包絡(luò)線，由于不確定信號的左、右兩個端點是否為極值點，會使擬合時產(chǎn)生誤差。所得到的上包絡(luò)線和下包絡(luò)線在信號的邊界附近產(chǎn)生大幅擺動，并且這種擬合誤差會隨著分解逐漸向信號內(nèi)部傳播，導(dǎo)致得到的IMF分量的波形嚴(yán)重失真，尤其是原始信號長度比較短時，會對經(jīng)驗?zāi)B(tài)分解的效果造成嚴(yán)重影響。

圖1為原始信號經(jīng)過分解后得到的各個IMF分量，可見IMF2的兩端有明顯的震蕩。這是HHT的第一個端點效應(yīng)[2-5]。而且在對各個IMF分量進(jìn)行Hilbert譜變換時，由端點效應(yīng)引起的誤差從數(shù)據(jù)的中心向外逐漸擴散，希爾伯特譜的兩端也產(chǎn)生突變和振蕩，如圖2所示。這是HHT中的第二個端點效應(yīng)。端點效應(yīng)的存在不利于信號的分析和處理。

2 抑制端點效應(yīng)的常用方法

2.1 鏡像對稱延拓

鏡像延拓的思路是在原始信號的兩端某位置各放一面鏡子，得到原始信號在鏡中的像，鏡子中原始信號的像關(guān)于鏡子與原數(shù)據(jù)對稱。

鏡像對稱延拓法關(guān)于鏡子位置的選擇，分為以下幾種情況：

對采樣起始點：

(1)第一個極值點是極大值時，如果第一個采樣點大于第一個極小值點，以第一個極大值為對稱中心；如果第一個采樣點小于第一個極小值點，則將其作為極小值點，并將其作為對稱中心。

(2)第一個極值點是極小值時，如果第一個采樣點小于第一個極大值點，以第一個極小值為對稱中心；如果第一個采樣點的值大于第一個極大值點，則將其作為極大值點，并將其作為對稱中心。

對采樣結(jié)束點：

(1)最后一個極值點是極小值時，如果最后一個采樣點的值小于最后一個極大值，將其作為對稱中心；如果最后一個采樣點的值大于最后一個極大值，則將其作為極大值點，并將其作為對稱中心。

(2)最后一個極值是極大值時，如果最后一個采樣點的值大于最后一個極小值，則把最后一個極大值點作為對稱中心；如果最后一個采樣點的值小于最后一個極小值，則將其作為極小值點，并將其作為對稱中心。

2.2 波形特征匹配延拓法

波形特征匹配延拓法基于波形匹配的思想，其基本思路是從信號內(nèi)部尋找一個與邊界處特征最相近的波形作為延拓數(shù)據(jù)的起點，按照所需延拓長度從信號內(nèi)部取出匹配波形數(shù)據(jù)前面(或后面)相應(yīng)長度的數(shù)據(jù)連接到邊界處，完成延拓。通過波形特征匹配延拓得到的數(shù)據(jù)符合信號的自然趨勢，能有效抑制Hilbert變換的端點效應(yīng)。

波形匹配延拓原理如圖3所示。

圖3 波形特征匹配延拓示意圖

以信號左邊界第一個極值點為極大值為例，Mi,Ni(i=1,2,…)分別為波形的極大值、極小值點，其時間點為tmi,tni,S1為第一個采樣點，波形特征匹配延拓法將S1-M1-N1作為特征波形，在全部波形中找到與S1-M1-N1構(gòu)成的三角形最接近的波形為匹配波形Sk-Mk-Nk，從匹配波形Sk-Mk-Nk開始，向前對波形數(shù)據(jù)進(jìn)行延拓。

具體步驟如下：

(1)根據(jù)S1與邊界處極值點的時間相對位置尋找對應(yīng)點Si，其時間坐標(biāo)為：

(4)

如果tsi不在采樣點上，可以用線性插值法計算其精確值si。時間點和采樣點可以互相轉(zhuǎn)化。

(2)計算匹配誤差。

Ei=|Si-S1|+|Ni-N1|+|Mi-M1|+|Mi+1-M2|

(5)

式中，|Mi+1-M2|為匹配波形的趨勢項，它反映了特征波形相對后一個極值點的位置。

(3)選取Ei值最小的波形為匹配波形，如果有多個Ei相等且為最小值，則取與左端點距離最遠(yuǎn)的波形為匹配波形。

(4)從匹配波形的右端的前一點開始，將實際波形復(fù)制到S1前，根據(jù)需要選擇延拓點數(shù)，如果信號中Si前數(shù)據(jù)點個數(shù)少于需要延拓的點數(shù)，可反復(fù)延拓此段波形。

(5)采用同樣的原理對右邊界進(jìn)行延拓。

3 波形匹配法和鏡像法相結(jié)合的新算法

3.1 算法步驟

對以上兩種延拓方法進(jìn)行理論研究，可發(fā)現(xiàn)：鏡像法雖然對端點是否為極值點進(jìn)行了判斷，但在以內(nèi)部極值點為對稱中心對信號進(jìn)行延拓時，需要截斷原始數(shù)據(jù)，這對短數(shù)據(jù)序列是不適用的；波形匹配法在出現(xiàn)信號中Si前數(shù)據(jù)點個數(shù)小于需要延拓的點數(shù)這種情況時，給出的解決方案是反復(fù)延拓此段波形。這會造成信號幅值的突變。為此，文中提出一種新的波形延拓方法：波形特征匹配法和鏡像法相結(jié)合的新算法。

算法步驟如下：

(1)先對原始信號x1的兩端進(jìn)行波形匹配延拓，得到信號x2。延拓的長度為端點到第2個極值點的長度。

(2)對x2進(jìn)行鏡像延拓，得到信號x3。鏡子的位置選在x2最外端的極值點處。

下文中將波形匹配法和鏡像法相結(jié)合的算法稱為先匹配后鏡像延拓算法。

將式(3)中的原始信號進(jìn)行先匹配后鏡像延拓，得到的延拓信號如圖4所示。

圖4 先匹配后鏡像延拓效果圖

3.2 端點效應(yīng)評價

一些學(xué)者對于端點效應(yīng)的評價標(biāo)準(zhǔn)做了研究。文獻(xiàn)[6]提出通過計算EMD后信號的總能量來評估端點效應(yīng)的影響程度;文獻(xiàn)[7]提出利用IMF與對應(yīng)信號分量的相關(guān)系數(shù)來評價EMD得到的IMF精度;文獻(xiàn)[8]提出考慮比較不同算法的運行時間。

文中采用下述4個指標(biāo)來評價多種端點效應(yīng)處理方法的效果:

(1)從能量的角度，端點效應(yīng)的實質(zhì)是信號的包絡(luò)線發(fā)生畸變，而產(chǎn)生了一些虛假成分，各IMF總能量會相應(yīng)增加，因此，可以比較EMD分解前后的能量來評估端點效應(yīng)的影響程度[9-14]。

(6)

其中：RMSorignal為原信號有效值;RMSi為第i個IMF的有效值;n為IMF總個數(shù)(包括EMD的殘余分量);θ的值越小，表示端點效應(yīng)的影響越小。

(2)比較EMD分解得到的模態(tài)分量和原信號相應(yīng)分量的相關(guān)系數(shù)ρ來評價IMF的精度。

(3)比較各個方法的運算時間t(單位為s)，在保證處理效果的同時，要避免算法過于復(fù)雜。

(4)EMD分解的正交性要滿足實際意義，即：

(7)

如果分解的信號滿足正交性，則式(7)右邊第二項為零，因此可用正交性系數(shù)ORT來衡量分解精度，ORT值越小表示時間序列的IMF分量正交性越好。

(8)

按照上述評價指標(biāo)，求得的幾種端點處理方法的評價指標(biāo)如表1所示。

從表中可以看出，從能量角度，先匹配再鏡像延拓后，信號的包絡(luò)線發(fā)生畸變最小，同時新方法分解得到的IMF精度較高，正交性好。

表1 端點處理方法性能比較

對式(3)的原始信號使用先匹配再鏡像的延拓方法，然后進(jìn)行EMD分解，并根據(jù)采樣點數(shù)截取延拓組IMF分量中對應(yīng)著原始信號的部分，得到的結(jié)果如圖5所示。

從圖中可看出，先匹配再鏡像延拓法分解信號的效果非常好，很好地抑制了EMD分解中出現(xiàn)的端點效應(yīng)，把端點效應(yīng)造成的誤差限制在了信號端點的外側(cè)，得到的IMF分量比較理想。

圖5 先匹配再鏡像延拓信號的IMF

由此可見，先匹配再鏡像延拓方法是一種有效的延拓方法，當(dāng)原始信號長度較短時，使用此方法對其進(jìn)行延拓，然后再使用EMD方法對其進(jìn)行分解，可有效抑制端點效應(yīng)的產(chǎn)生。

4 結(jié)束語

介紹了兩種常用的延拓方法，并提出了一種先匹配再鏡像延拓新方法，給出了具體操作步驟，并進(jìn)行了實例仿真。

給出了端點效應(yīng)評價的4種指標(biāo)，對先匹配再鏡像延拓進(jìn)行了端點效應(yīng)評價，并與原有的鏡像匹配延拓與波形匹配延拓在不同評價指標(biāo)下進(jìn)行了詳細(xì)比較。

通過比較可以證明，使用先匹配再鏡像延拓的方法，增加了EMD分解得到的IMF的正交性，提高了EMD分解精度，改進(jìn)了原有延拓方法的不足。

[1]HuangNE,WuMC,LongSR,etal.AconfidencelimitfortheempiricalmodedecompositionandtheHilbertspectralanalysis[J].ProcofRocSocLondon,SerA,2003,459:2317-2345.

[2] 高 強,段晨東,趙艷青,等.基于最大相關(guān)波形延拓的經(jīng)驗?zāi)Ｊ椒纸舛它c效應(yīng)抑制方法[J].振動與沖擊,2013,32(2):62-66.

[3] 楊 斌,陳桂明,劉建友.基于標(biāo)準(zhǔn)互相關(guān)函數(shù)的經(jīng)驗?zāi)B(tài)分解端點效應(yīng)處理方法[J].機械工程學(xué)報,2013,49(5):63-68.

[4] 白春華,周宣赤,林大超,等.消除EMD端點效應(yīng)的PSO-SVM方法研究[J].系統(tǒng)工程理論與實踐,2013,33(5):1298-1306.

[5] 任達(dá)千,楊世錫,吳昭同,等.LMD時頻分析方法的端點效應(yīng)在旋轉(zhuǎn)機械故障診斷中的影響[J].中國機械工程,2012,23(8):951-956.

[6]BenkherratM,BauguerraR,ChoufaT.RemovalofocularartifactsfromrelatedevokedpotentialsusingVSSLMSadaptivefilter[C]//Procofinternationalconferenceon"computerasatool".Warsaw,Poland:[s.n.],2007:349-352.

[7]VenkataramananS,KalpakamNV.Anovelwaveletbasedtechniquefordetectionandde-noisingofqcularartifactinnormalandepilepticelectroencephalogram[C]//ProcofICCCAS.[s.l.]:[s.n.],2004:27-29.

[8]KumarPS.RemovalofocularartifactsintheEEGthroughwavelettransformwithoutusinganEOGreferencechannel[J].IntJOpenProblemsComptMath,2008,1(2):188-200.

[9] 孟 宗,顧海燕,李姍姍.基于神經(jīng)網(wǎng)絡(luò)集成的B樣條經(jīng)驗?zāi)B(tài)分解端點效應(yīng)抑制方法[J].機械工程學(xué)報,2013,49(9):106-112.

[10] 曹沖鋒.基于EMD的機械振動分析與診斷方法研究[D].杭州:浙江大學(xué),2009.

[11] 黃先祥,李勝朝,謝 建.新型經(jīng)驗?zāi)Ｊ椒纸舛它c效應(yīng)消除方法[J].機械工程學(xué)報,2008,44(9):1-5.

[12] 朱曉軍.HHT變換及其在腦電信號處理中的應(yīng)用研究[D].太原:太原理工大學(xué),2012.

[13] 楊賢昭.基于經(jīng)驗?zāi)B(tài)分解的故障診斷方法研究[D].武漢:武漢科技大學(xué),2012.

[14] 張進(jìn)林,張榆鋒,張 燕,等.經(jīng)驗?zāi)B(tài)分解端點效應(yīng)抑制的常用方法比較研究[J].云南大學(xué)學(xué)報:自然科學(xué)版,2010,32(4):406-412.

Research on Method for End Effects Reduction of Empirical Mode Decomposition

GUO Di,SHEN Yang-yang,YIN Bing

(College of Electronic Science and Engineering,Nanjing University of Posts and Telecommunications,Nanjing 210003,China)

The precision of Empirical Mode Decomposition (EMD) is reduced greatly by its end effects,so a new end extending method combining the waveform feature matching extending method and mirror extending method is presented.The waveform feature matching extending method and mirror extending are used for reference in this new method,the signal is extended by waveform feature matching method at first and then extended by mirror method.A simulation signal is applied to test the performance of the new method,and a comparison under different evaluating indicators between the new method and old methods is made and analyzed.The result shows that the signal envelope has minimum distortion,at the same time,the IMF has high precision and good orthogonality decomposed by new method.The proposed improved extending method can improve the precision of EMD and restrain the end effect effectively.

EMD;extending;end effect;mirror method;waveform feature matching method

2015-06-04

2015-09-08

時間：2016-02-18

國家自然科學(xué)基金資助項目(61003237)；江蘇省高校自然科學(xué)研究(10KJB510018)

郭 迪(1991-)，男，碩士，研究方向為智能信息處理。

http://www.cnki.net/kcms/detail/61.1450.TP.20160218.1630.022.html

TP301

A

1673-629X(2016)03-0089-04

10.3969/j.issn.1673-629X.2016.03.021