基于鏡像延拓和窗函數的端點效應抑制方法

徐力彬，宋余慶，劉 毅

(江蘇大學計算機科學與通信工程學院，江蘇鎮江212013)

1 概述

Hilbert-Huang變換(Hilbert-Huang Transform，HHT)是由Huang N E等人通過深入分析和總結后提出的一種新的非平穩非線性信號分析方法［1］。該方法主要包含2個部分:經驗模態分解(Empirical Mode Decomposition，EMD)和 Hilbert譜 分 析(Hilbert Spectral Analysis，HAS)，其中 EMD 為核心部分［2］。與傳統方法相比，EMD無須預先設定任何基函數，而是根據信號特點自適應地將原始信號分解為一系列固有模態函數(Intrinsic Mode Function，IMF)，去除了快速傅里葉變換方法中引入的諧波基函數，也沒有小波分析復雜信號過程中產生信號模糊與失真的問題［3-6］。然而，HHT變換也有一定的局限性，例如端點效應問題［7］和模態混疊問題［8］等。尤其是端點效應問題會嚴重妨礙信號的特征提取，從而導致信號分析的誤差。

EMD中使用三次樣條函數對包絡進行擬合時，由于難以確定端點處是否為極值點，導致出現了端點效應問題。國內外專家針對端點效應已經提出了多種有效的抑制方法，包括鏡像延拓［9］、極值延拓［10］、基于神經網絡預測的數據延拓法［11］以及多項式擬合延拓法［12］等。這些方法在抑制端點效應方面都有著良好的效果，但同時也存在弊端，如多種延拓方法存在同一種缺陷，即在對原始信號進行延拓處理以后，其端點依然是不固定的，擬合的上下包絡線兩端還是會出現發散現象，只是在EMD完成后截去延拓部分，端點效應相對減輕。本文在鏡像延拓的基礎上，利用余弦窗函數法對延拓后的信號進行處理，使樣條函數可以更好地擬合包絡線，對端點效應起到抑制作用。

2 EMD算法

2．1 算法簡介

EMD算法是HHT的核心算法之一，用來將一個信號的能量按照時域各種固有尺度的波動進行分解，得到一系列頻率從小到大的IMF。這里IMF須滿足如下2個性質［13］:

(1)信號的極值點數目和過零點數目相等或最多相差一個;

(2)由局部極大值構成的上包絡線和由局部極小值構成的下包絡線的平均值為0。

EMD的分解過程如下:

Step1設原始信號曲線為X(t)，求出所有的極大值、極小值。

Step2對極值使用三次樣條函數進行處理，描繪出原始信號的上、下包絡線，分別記為 u(t)和v(t)。

Step3將上、下包絡線的平均值記為:

Step4將均值曲線m(t)從原信號X(t)中分離出來，得到的剩余部分記為:

Step5判斷h1(t)是否滿足IMF的2個性質。若滿足，則h1(t)為第一個IMF分量;否則，記h1(t)為X(t)，重復Step1～Step4，直到第k次迭代之后得到一個IMF，將其記為c1(t)=h1k(t)。這里，將標準差準則作為第k次迭代的終止原則，使其位于0．2～0．3 之間，停止篩選［14］。

Step6記r1(t)=X(t)-c1(t)為新的待處理信號，再重復Step1～Step5，得到第2個IMF分量，將其記為c2(t)。繼續重復Step1～Step6，當所得余項rn(t)是一個單調信號或者rn(t)的值小于事先給定的閾值時，分解結束［15］。

綜上所述，原始信號X(t)為所有IMF分量與余項rn(t)的和:

2．2 端點效應

在EMD中，為得到信號的瞬時平均包絡線，需要對原信號的局部極大值和極小值分別進行三次樣條插值算法擬合出上、下包絡線后，再計算出信號的局部平均值。由于信號兩端無法同時處于極大值和極小值點，因此上、下包絡線在數據信號兩端不可避免地會出現偏移，并且這種偏移現象會隨EMD的進行逐漸傳播到信號內部，從而嚴重影響信號分解的質量［16］，這種現象就是端點效應，如圖1所示，其中，實心點為延拓后的點，實線為真實的包絡線，虛線為無端點約束產生的錯誤包絡線。

圖1 端點效應示意圖

2．3 抑制端點效應的常用方法

2．3．1 常用延拓方法

延長數據序列或者在數據兩端增加極值點來抑制端點效應是得到普遍認可的一種方法［17］。對于復雜的信號，準確地對數據序列進行延拓是不現實的，但可以使延拓后得到的平均包絡線與真實平均包絡線相對較為接近，這是到目前為止普遍認為的抑制端點效應的有效途徑。常用的延拓方法如下:

(1)鏡像延拓法

鏡像延拓法的實質是將原始信號對稱地延拓成一個閉合的環形信號，這種環形信號沒有端點，其上、下包絡線完全由信號內部確定，從而在根本上避免了端點效應的產生。但是，當無法確定端點是否為極值點時，可能會使延拓部分和原始信號的均值有明顯差異，從而影響EMD產生的最終結果。

(2)極值延拓法

極值延拓法不需要對信號本身延拓，可以只對信號的極值進行延拓。該方法在鏡像延拓的基礎上，對延拓后的對稱點是否為極值點進行了判斷，提高了處理效果［18］。但該方法也存在缺陷，即若將端點都作為極值點進行處理，端點處的包絡線會產生收縮，從而導致包絡線失真。

(3)多項式擬合延拓法

通過多項式擬合法算出信號端點處對應的函數值，將其作為端點處極值的近似值。用正交多項式擬合的方法會比一般的多項式擬合更加準確，算法復雜度也不搞。該方法對準周期信號延拓效果比較好，但對隨機信號等變化規律不明顯的信號在端點效應抑制方面的效果一般。

2．3．2 窗函數法

窗函數法是一種能夠有效抑制端點效應的方法，對信號進行窗函數處理可以增強函數中心處的信號，同時抑制遠離中心處的信號［19］。在EDM 分解過程中，端點效應從端點逐漸向信號內部發散，誤差不斷積累，導致發散程度越來越嚴重。

在對信號進行窗函數處理后，信號端點處的值變為0，包絡線也會相對平滑。因此，三次樣條函數能夠更真實地擬合包絡線，最終上、下包絡線收逐漸收斂于端點，不再發散。

在窗函數的選取問題上，筆者不需要自己設計窗函數，國內外數學家已經提出了很多能適合不同需要的窗函數。在處理端點效應的問題上，應用較多的有漢寧窗、海明窗、矩形窗和余弦窗等。

3 本文算法

當原信號進行窗函數處理后，原信號的包絡線將被改變，可能會對接下來的過程帶來誤差，這時可以考慮先將信號進行延拓處理，再進行加窗，最后在EMD過程結束后舍去延拓部分。

仿真信號表達式如下:

其中，f1=200 Hz;f2=100 Hz;f3=50 Hz。采樣頻率為1 000 Hz。信號由3個正弦信號組成。圖2為原信號x(t)分解出來的3階IMF，分解時沒有對原始信號進行任何處理，可見IMF1，IMF2，IMF3都產生了明顯的端點效應，而且隨著分解過程的進行，端點的發散程度越來越嚴重，在IMF3中甚至已經無法表達出原信號中的低頻部分。

圖2 未經處理的EMD分解結果

本文方法的具體步驟如下:

(1)將原信號x(t)進行鏡像延拓處理，得到延拓后的信號y(t)。

(2)將延拓后的信號y(t)用一定周期ΔT的余弦窗函數α(t)進行處理。即將y(t)與α(t)進行內積運算，得到信號 σ(t)=［x(t)，α(t)］，確定處理信號σ(t)的所有極值點。

(3)計算由局部極大值和極小值點確定的上下包絡線Mmax(t)和Mmin(t)。

(4)計算包絡線均值:

m(t)=(Mmax(t)+Mmin(t))/2

(5)求出h(t)=y(t)-m(t)。

(6)如果h(t)不滿足IMF的條件，則重復以上循環;如果h(t)是一個IMF，將信號y1(t)=y(t)-h(t)作為原始信號重復以上循環。當yk(t)成為一個單調信號時，循環結束。

4 仿真實驗及結果分析

4．1 實驗分析

4．1．1 延拓處理

采用余弦窗函數對原始信號進行處理后，有可能會改變信號的特征，因此，改進算法先對信號進行鏡像延拓，使延拓信號和原始信號的交界處變得光滑。

本文將仿真信號x(t)作為延拓對象，在從左向右的第L個極值點處和從右向左的第R個極值點處分別設置2個光滑的反射面，將信號向外反射，得到數據序列長度兩倍于反射信號長度的周期性信號。這時，信號作為一個環形的閉合曲線，滿足了原信號曲線在端點處一階、二階倒數存在的要求。然后將閉合曲線作為一個整體進行EMD分解，最終取上半部分作為分解結果。圖3所示結果為延拓后的信號，可以看出，此時信號兩端點處已呈一定的周期性。

圖3 鏡像延拓后的EMD分解結果

4．1．2 加窗處理

余弦窗是信號處理中應用較為廣泛的一種窗函數，在本文中，余弦窗函數可被定義為:

式(6)的曲線圖如圖4所示，可見，在窗函數中間部分的值為1，窗函數兩端的值逐漸減少至0。這樣，信號端點處的值變為0，包絡線變得相對平滑，擬合出的包絡線也變得更為真實。隨著EMD分解的進行，上下包絡線逐漸收斂于端點，發散現象就不會再發生。圖5為延拓后信號進行余弦窗函數處理后的情況，可見IMF1，IMF2，IMF3的端點效應均得到了明顯的抑制，而且IMF4分量中的低頻部分仍能表達出來。

圖4 余弦窗函數α(t)

圖5 加窗處理后的EMD分解結果

4．2 評價指標

目前，所有提出的端點效應抑制方法都需要一個定性與定量相結合的評價指標來綜合考慮方法的優劣。本文主要采用均方根有效值評價法來定性分析端點效應抑制方法的效果:可以先計算出原始信號和通過篩選得到的各個IMF分量的均方根有效值。如式(7)所示。

其中，RMS表示信號有效值;s(i)為信號序列;n為信號的采樣點數。

按式(7)比較各IMF分量有效值的總和與原信號有效值，得到一個評價指標θ。

其中，RMSoriginal表示原始信號的有效值;RMSi表示第i個IMF分量的有效值;n為IMF分量的個數，其中包括EMD的殘留項。

根據定義，θ≥0，當端點效應對EMD分解過程沒有影響時，θ=0。θ值越大，說明端點效應對EMD分解的影響越大。θ值越小，表明IMF分量的發散程度越小，端點效應的影響也就越小。

4．3 實驗結果

實驗給出了未加任何處理的EMD算法、極值對稱延拓法、鏡像延拓法、窗函數法和改進算法對處理端點效應的效果評價，給出了EMD分解前后的效果對比參數值 θ。從表1中可以看出，未加處理的EMD算法的θ值明顯偏大，極值對稱延拓法、鏡像延拓法、窗函數法的θ值都有不同程度的減小，而本文改進算法的θ值較其他4種算法都小，說明對端點效應的抑制效果更好。

表1 多種算法的端點效應抑制評價結果

5 結束語

本文提出了一種將鏡像延拓和余弦窗函數結合的端點效應抑制方法。首先將原信號進行鏡像延拓處理，然后在對延拓后信號進行IMF提取之前，在其兩端加上余弦窗函數，使得上下包絡線能夠更好地擬合。該方法彌補了鏡像延拓可能使延拓部分信號與原始信號的均值有明顯差異的不足，同時也解決了單純窗函數法因為改變原信號而可能帶來其他問題的缺陷。另外，本文利用θ值作為評價指標綜合評價了多種抑制算法的效果，結果顯示出改進方法較傳統方法略有提升。但是，在EMD分解過程中，由于每次篩選得出的新信號發散程度不同，如果使用單一窗函數對篩選后的新信號進行處理可能還會帶來信號失真等問題，隨著EMD的進行，這種信號失真現象可能越來越嚴重。因此，對于EMD分解中的端點效應問題還有待進一步研究。

［1］ Huang N E，Shen Zheng，Long S R，et al．The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear and Nonstationary Time Series Analysis［J］．Proceedings ofRoyalSociety London，1998，454:903-995．

［2］ 王黎黎．基于希爾伯特-黃變換的時頻分析算法研究［D］．西安:西安電子科技大學，2009．

［3］ 秦 毅，秦樹人，毛永芳．改進的Hilbert-Huang變換在信號瞬態特征提取中的應用［J］．沖擊與振動，2008，27(11):129-130．

［4］ Guo D，Peng Z P．Vibration Analysis of a Cracked Rotor Using Hilbert-Huang Transform［J］．Mechanical Systems and Signal Processing，2007，21(8):3030-3041．

［5］ Li X L，Bassiunya A M．Flute Breakage Detection During End Milling Using Hilbert-Huang Transform and Smoothed Nonlinear Energy Operator［J］．International Journal ofMachine Toolsand Manufacture，2007，47(6):1011-1020．

［6］ Feldman M．Considering High Harmonics for Identification of Non-linear System by Hilbert Transform［J］．MechanicalSystems and SignalProcessing，2007，21(2):943-958．

［7］ Yan Ruqiang，Gao R X．A Tour of the Hilbert-Huang Transform:An Empirical Tool for Signal Analysis［J］．IEEE Instrumentation ＆ Measurement Magazine，2007，10(5):40-45．

［8］ 張 超，陳建軍．EEMD方法和EMD方法抗模態混疊對比研究［J］．振動與沖擊，2010，29(8):87-90．

［9］ 蘇玉香，劉志剛，李科亮，等．一種改善EMD端點效應的新方法及其在諧波分析中的應用［J］．電工電能新技術，2008，27(2):33-36．

［10］ Zhao Jinping，Huang Daji．Mirror Extending and Circular Spline Function for Empirical Mode Decomposition Method［J］．Journal of Zhejiang University，2001，2(3):247-252．

［11］ 顧小麗，李培良，譚海濤，等．基于 RBF神經網絡的EMD方法在海平面分析中的應用［J］．海洋與湖沼，2009，40(5):532-534．

［12］ 朱金龍，邱曉暉．正交多項式擬合在EMD算法端點問題中的應用［J］．計算機工程與應用，2006，42(23):72-74．

［13］ 甘錫林，黃韋艮，楊勁松，等．基于希爾伯特-黃變換的合成孔徑雷達內波參數提取新方法［J］．遙感學報，2007，11(1):39-40．

［14］ 張 博．基于Hilbert-Huang變換的滾動軸承故障診斷方法研究［D］．西安:西安理工大學，2009．

［15］ 戴桂平．基于EMD的時頻分析方法研究［D］．秦皇島:燕山大學，2005．

［16］ Zeng Kan，HeMingxia．A SimpleBoundaryProcess Technique for Empirical Decomposition［C］//Proceedings of International Geoscience and Remote Sensing Symposium，September 20-24，2004，Anchorage，USA．Washington D．C．，USA:IEEE Press，2004:4258-4261．

［17］ 許寶杰，張建民，徐小力，等．抑制EMD端點效應方法的研究［J］．北京理工大學學報，2006，26(3):196-200．

［18］ 黃大吉，趙進平，蘇紀蘭．希爾波特-黃變換的端點延拓［J］．海洋學報，2003，25(1):1-9．

［19］ 任達千，吳昭同，嚴拱標．EMD端點效應的評價指標及抑制端點效應的窗函數法［J］．制造業自動化，2007，29(1):21-24．