基于延拓與可變余弦窗的經驗模態分解改進算法研究

丁雪娟，時培明

（燕山大學電氣工程學院，秦皇島 066004）

基于延拓與可變余弦窗的經驗模態分解改進算法研究

丁雪娟，時培明

（燕山大學電氣工程學院，秦皇島 066004）

針對經驗模態分解（Empirical Mode Decomposition，EMD）中存在的邊界效應及邊界發散現象隨著篩選層次的增加而增加的問題，提出一種利用延拓與可變余弦窗相結合的改進新方法。首先對信號進行延拓處理，增加一定長度的數據，實現延拓數據與原始信號交界處的光滑過度。其次，根據信號邊界的發散程度，在逐層提取各階本征模函數（Intrinsic Model Function，IMF）之前，在信號兩端加上寬度可變的余弦窗函數，使得每一個IMF分量邊界發散問題最小化，保證信號有效數據的正確分解，實現EMD邊界處理算法的改進。仿真和實例信號分析表明，該方法能較好地抑制EMD邊界效應，有效地提取故障信號中的特征信息。

信號處理；經驗模態分解；邊界效應；可變余弦窗

經驗模態分解是Huang首先提出的一種新的自適應的時頻分析方法［1］。它根據信號的局部時變特征進行自適應的分解，避免了人為因素，克服了傳統方法中的非局域化及非自適應性，具有良好的時頻聚集性，適合處理瞬變、調幅或調頻等非線性、非平穩信號。自EMD方法問世以來，它就引起了眾多學者的廣泛關注，已被應用于設備診斷、醫學信號分析、語音信號的檢測、圖像信息分離等領域［2－4］。

邊界效應一直是制約EMD發展的一個關鍵問題，針對這一問題，研究人員已經提出了一些改進方法，如鏡像延拓［5］、神經網絡延拓［6］、相似極值延拓［7］、波形特征匹配延拓［8］、支持矢量回歸機［9］等。但是，這些延拓方法仍然存在一定的誤差，導致其邊界發散問題依然存在。窗函數可以將延拓誤差控制在信號的兩端，使其無法（或以較慢速度）向數據內部發展。任達千等［10］提出了窗函數法對抑制邊界效應有明顯的效果，但指出加窗會改變原始信號，影響分解精度。Qi等［11］提出了在原始信號上加余弦窗函數，雖然對提取低頻微弱信號具有一定作用，但是固定窗函數只能抑制特定信號的邊界效應，對最后的分解結果有一定的影響。Parey等［12］提出了在IMF分量上加可變余弦窗函數，該方法必須保證每個IMF分量具有明確的物理意義，才能使得方法有效，但是在EMD分解過程中有些微弱低頻分量難以篩選出來，無法反映真實的物理過程，使得該方法具有局限性。

本文提出了一種基于延拓與可變余弦窗相結合的改進方法。該方法首先將信號延拓處理，增加一定長度的數據，避免在加窗函數時造成數據的改變。其次，每次提取IMF前，在信號兩端加上寬度可調的余弦窗函數，使得包絡線變得比較平滑，沿著正確的方向延展。為了降低窗函數帶來的誤差，該方法設計窗函數的寬度可調，提高了IMF的分解精度。最后對EMD分解得到的IMF進行Hilbert變換，按原始信號長度及位置提取分析結果，得到有效數據的分析結果。仿真和應用實例分析驗證了該方法的有效性。

1 EMD邊界效應及改進算法

1.1 邊界效應

EMD時頻分析方法處在發展階段，其理論還需要不斷完善。邊界效應是其應用過程中的一個重要問題。EMD分解中由于無法保證數據端點處的極值點，導致求包絡平均過程中，會在樣條插值時產生數據的擬合誤差。并且隨著分解的進行，誤差不停積累，由端點處向內逐漸傳播，嚴重時會使分解的數據失去意義。同時，用數字方法實現的Hilbert變換是基于Fourier變換的，在對本征模分量進行Hilbert變換時，由于首尾兩點的數值不同，導致頻譜的泄漏，表現在時頻譜圖上便是邊界效應，影響了時頻分析的精度。

1.2 基于延拓和可變余弦窗的邊界效應改進算法

在EMD分解過程中，邊界效應是在樣條差值時產生的擬合誤差，從端點附近開始出現，并且逐漸向內“污染”。隨著分解的進行，誤差不停積累，導致IMF從高頻到低頻的發散程度越來越嚴重，甚至可能無法正確的提取出低頻信號。如果信號的端點也是上、下包絡線的端點，邊界效應即可大大減輕。余弦窗窗函數在保證中心點附近有用信號的特征時，使得端點趨近于0，進而使得上下包絡線收斂于端點，有效地抑制了邊界發散現象。

為防止加窗改變原始信號，首先，采用波形特征匹配延拓方法，將信號延拓一定的數據長度。由于波形特征匹配延拓方法兼顧了原始極值點和非極值點數據，可以使得延拓數據與原始信號交界處的比較光滑，進而避免了加窗引起的原始頻譜泄漏以及邊界瞬時頻率的跳躍。波形特征匹配延拓通過采用原始信號內部和邊界處變化趨勢最為相似的子波來對端點處數據進行延拓，是一種自適應的方法。在具體實現中，通過計算波形匹配來量化兩端波形的變化趨勢。

如圖1所示，以左邊界第一個極值點為極大值為例，mi、ni（n＝1，2，3，…）分別為曲線的極大值、極小值點，分別對應時間tmi、tni，a1為左端點，波形特征匹配延拓法以a1－m1－n1為邊界特征波形，在全部數據中找到與a1－m1－n1構成的三角形最接近的波形為匹配波形，如ai－mi－ni，從ai的前一點（右邊界為后一點）數據開始，向前（右邊界向后）延拓波形數據，使延拓數據符合信號的自然走向。

圖1 波形特征匹配原理圖Fig.1 Diagram ofwave characteristicmatching

其次，為降低窗函數帶來的誤差，改進方法設計窗函數的寬度為可調。逐層提取IMF前，根據邊界發散的程度，在信號兩端加上可調的余弦窗函數，使得包絡線沿著正確的方向延展，提高了IMF的分解精度。

可變余弦窗函數定義為：

式中：L為信號延拓后的長度，ΔT為控制窗體中間部分時間范圍的可變參數，它根據信號的改變而改變。由于邊界發散程度是隨著IMF的增加而不斷增加的，根據篩選不同IMF前信號的發散程度而選取不同的值選取ΔT的大小。由圖2可見，在窗函數的中部，其幅值等于1，而在窗函數的兩端幅值逐漸衰減，直至窗函數的2個端點處幅值減至為0。寬度可調的余弦窗函數能夠較準確地將延拓誤差控制在信號兩端，抑制邊界的發散，實現算法的改進。

圖2 可變余弦窗函數示意圖Fig.2 Variable cosine window function diagram

基于延拓和可變余弦窗的EMD改進方法的步驟如下：

（1）對原始信號x（t）進行延拓，加一定ΔT的余弦窗ω（t）。將延拓信號x′（t）與余弦窗函數進行內積運算，得到信號y（t）＝［x′（t），ω（t）］，確定處理信號y（t）的所有極值點。

（2）求由局部極大值點和局部極小值點確定的上包絡線和下包絡線，Lmax（t）和Lmin（t）。

（3）計算包絡線均值

a（t）＝（Lmax（t）＋Lmin（t））／2（2）

（4）求出h（t）＝y（t）－a（t）。

（5）如果h（t）不滿足IMF條件，重復上面循環；如果h（t）是一個IMF，將信號f（t）＝x（t）－a（t）作為原始數據重復上面的循環，此時窗函數ΔT根據發散程度而改變。

當f（t）成為一個單調函數不能再從中提前出滿足IMF條件的分量時，循環結束。這樣就得到：

（6）對式（3）中的每一個本征模函數ci（t）作Hilbert變換，得到Hilbert譜及其邊際譜。

2 端點效應的評價指標

（1）主觀評價通過觀察信號分解結果的發散程度，主觀判斷邊界效應對于EMD結果的影響程度。

（2）相關系數對于仿真信號，一般有理論可以參照，邊界效應的影響也容易評估。利用相關系數作為參數指標，可表征EMD分解結果與理想分解結果的接近程度，其基本原理如下［10］：

設a（n）、b（n）是兩個能量有限的確定性信號，并假定它們是因果的，則定義

由上式可知，當a（n）＝b（n）時，rab＝1，則表明兩個信號完全相關；當a（n）與b（n）完全不相關時，rab，當a（n）與b（n）有某種程度相似時，0＜｜rab｜＜1。因此，rab可以用來描述a（n）與b（n）之間的相似程度。

（3）能量角度邊界效應會使IMF總能量會相應增加，因此可以比較EMD分解前后的能量來評估邊界效應的影響程度。

原始信號和EMD得到的IMF的均方根有效值計算公式為：

式中：RMS為信號有效值，u（n）為信號序列，n為信號的采樣點數。

定義評價指標θ

式中：RMSorginal原始信號有效值，RMSi為第i個IMF的有效值，n為IMF總個數。如果θ＝0，說明沒有邊界效應，θ值越大說明邊界效應的影響越大。

3 仿真實驗

仿真信號的表達式：

x（t）＝sin（2π×f1t）＋sin（2π×f2t）＋0.02sin（2π×f3t）（7）

式中：f1＝200 Hz，f2＝100 Hz，f3＝50 Hz，采樣頻率為2 000Hz。信號由三個正弦信號組成，低頻信號比較微弱。仿真信號中由幾個正弦諧波分量組成，呈現出一定的周期性。EMD分解過程實質是將信號不斷平穩化的過程，最終將信號分解為滿足一定條件的有限個IMF組，每一個IMF在局部吻合標準的正弦曲線，而且不同的IMF分量包含了不同的特征時間尺度。

如圖3所示，分解時沒有對原始數據進行任何處理。信號x（t）分解出來的3階IMF都有一定程度的邊界效應，而且邊界發散程度逐漸增加，尤其是IMF3分量，已經無法表達出原信號中的低頻小能量成分。

對原始信號進行延拓，如圖4所示，對延拓后的信號進行EMD處理，得到圖5，由圖可以看出延拓后的數據仍然存在邊界發散現象，尤其對于低頻小能量的微弱信號尤為嚴重。將延拓信號x（t）與余弦窗函數進行內積運算，得到信號y（t）＝［x（t），ω（t）］。

圖3 未處理的EMD分解結果圖Fig.3 Untreated EMD decomposition resultsmap

圖4 延拓后的信號Fig.4 Continuation signal

圖5 延拓后的EMD分解Fig.5 The continuation EMD decomposition

圖6為不同窗函數對應篩選不同IMF前的加窗結果，圖6（a）中ΔT分別為43，80，110，逐層分解IMF前的加窗結果如圖6（b）－（d）。從圖中可以看出，信號兩端的發散現象都得到了較好的抑制。圖7為改進算法對信號x（t）分解截取有效數據與理想值的比較結果，由圖中可以看出，3個IMF分量包括低頻微弱信號被有效地分離出來，并且邊界效應得到了很好的抑制。

圖6 窗函數和加窗結果Fig.6 Window function and window results

圖7 改進算法的分解結果Fig.7 Decomposition results of improved algorithm

圖8、9分別是信號x（t）未經過處理和用改進方法得到的Hilbert譜及其邊際譜圖。由圖8可以看出，未處理的Hilbert譜在信號兩端有比較嚴重的發散現象，而且發散程度從高頻到低頻逐漸增加，而通過本文方法處理后的Hilbert譜效果有明顯的改善，同時低頻分量的頻率和幅值也顯示出來。從圖9中也能清楚地看到未經過處理的圖像兩端有很多離散的能量，而處理后的圖像兩端能量就相對集中，信號失真情況有明顯改善。圖9（a）中的幅值顯示明顯有很多微弱振蕩，尤其是在低頻部分出現了虛假成分。通過邊界效應處理后得到的邊界譜（b）中，幅值振蕩基本消除，信號所包含的3個頻率成分準確地顯示出來。

圖8 三維Hilbert譜對比圖Fig.8 Comparison of three dimensional Hilbert spectroscopy

圖9 邊際譜對比圖Fig.9 Comparison ofmarginal spectrum

觀察信號分解結果的的發散程度，可以看出改進的方法有效地抑制了邊界發散現象。同時表1給出了未處理方法和改進算法的邊界效應的性能評價，包括IMF分量與理想分量的相關性，以及EMD分解前后的能量對比參數θ值。相關系數r越接近于1，表明得到的結果越接近理想的分解結果。θ越小，表明虛假成分越少以及IMF分量的發散程度越小，邊界效應的影響越小。從表1中可以看出未處理的EMD方法篩選出的IMF3與理想值基本上不相關，無法正確表達原信號中的低頻小能量成分。而用改進算法的分解后，IMF1和IMF2與理想值的相關度提高，IMF3的rab也接近1，證明分解結果非常接近理想的分解結果，同時改進的算法使得θ值大大降低，說明邊界效應減小，驗證了改進算法的有效性。

表1 未處理的EMD算法以及改進算法的邊界效應評價表Tab.1 Boundary effect evaluation of untreated EMD algorithm and an im proved algorithm

4 實例分析

將本文提出的EMD改進算法應用于橫向裂紋故障信號的特征提取及診斷。由于裂紋故障信號包含信號的特征頻率及其倍頻成分，具有周期性，進而使得每個IMF分量趨近于正弦曲線，便于觀察IMF的發散程度。另外，信號內部的周期性，使得本文提出的方法更具優勢，這種規律使我們在尋找與邊緣處波形相似的子波更加容易，減小了因為延拓產生的誤差。圖10是一個含橫向裂紋故障轉子的振動信號，轉速為480 r／min。

圖11是未處理算法分解結果，前3個分量中IMF1是噪聲信號，IMF2、IMF3是含有裂紋故障特征的轉子振動信號。從分解圖中可以看出，分量兩端存在發散現象，破壞了信號的規律性，并且隨著分解層次的增加，發散現象越明顯。

圖10 轉子裂紋的實際信號Fig.10 The actual signal of rotor crack

圖11 未處理的EMD分解結果Fig.11 Untreated EMD decomposition results

根據觀察IMF兩端發散的程度，確定窗函數的寬度，使得信號的兩端幅值逐漸衰減，直至歸為0，進而抑制信號邊界發散現象。圖12是對原始信號延拓后，逐層提取IMF分量之前，加可調余弦窗的結果。圖12（a）中窗函數的ΔT分別為43，53，63。圖12（b）～（d）分別是所加的窗函數及加窗后效果。

圖13為改進算法的分解結果，對比圖11、圖13可以看出，IMF2與IMF3的邊界發散現象得到了很好的改善，各個分量的規律性也更明顯。同時，從能量角度判定邊界效應度的θ值由0.004 3減小為0.000 724，說明了改進算法的有效性。

圖14為改進前后信號的Hilbert譜，圖15為信號的Hilbert邊際譜。比較圖14，改進后的Hilbert譜的邊界發散現象得到了抑制，信號的的能量也相對集中。比較圖15，改進后的邊際譜幅值微弱振蕩基本消除，且2倍頻分量得到了突出，其調頻現象更加的突出。特征頻率和2倍頻的存在，是橫向裂紋故障的特征，但是與不對中故障區分不明顯。

圖12 加可調余弦窗的結果Fig.12 Adjustable cosine window results

圖13 改進算法的分解結果Fig.13 decomposition results of improved algorithm

圖14 改進前后的Hilbert譜對比圖Fig.14 Comparison of Hilbert spectrum before and after improvement

為進一步說明設備存在的是橫向裂紋故障，采集一組高轉速下的數據，轉速為900 r／min。如圖16所示，該組數據為單模態分量。圖17為原始數據的Hilbert譜。可以看出較高轉速下的軸頻信號產生了頻率調制，其調制頻率與軸頻相同，由此現象可作為轉子橫向裂紋故障的評判依據。根據以上分析，診斷該轉子的故障為裂紋故障，同時驗證了改進算法的有效性。

圖15 改進前后的邊際譜對比圖Fig.15 Comparison ofmarginal spectrum before and after improvement

圖16 高轉速下的振動信號Fig.16 High speed vibration signal

圖17 振動信號的Hilbert譜Fig.17 Hilbert spectrum of vibration signal

5 結 論

本文提出了一種延拓與可變余弦窗相結合的抑制經驗模態分解邊界效應的新方法。該方法先將信號延拓處理，根據信號的發散程度，逐層提取IMF之前，在信號兩端加上寬度可變的余弦窗函數，使得包絡線沿著正確的方向延展。由于窗函數的寬度設計為可以調節，使得每一個IMF分量邊界發散問題最小化，降低了窗函數帶來的誤差，提高了IMF的分解精度。通過仿真實驗改進算法與未經處理EMD算法的分解結果對比，可以看出改進算法使得邊界效應得到了很好的抑制。另外利用邊界效應性能評價指標分析了改進算法的分解結果，同樣得出了改進算法抑制邊界效應的有效性。EMD改進算法在故障信號特征提取中的應用實例，說明了改進方法能夠有效地提取出故障信號的特征信息。

［1］Huang N E.The empirical mode decomposition and the hilbert spectrum for nonlinear and non-stationary time series analysis［J］.Proc Rsoc Lond A，1998，454：903－955.

［2］Shen Z J，Chen X F，Zhang X L，et al.A novel intelligent gear fault diagnosis model based on EMD and multi－class TSVM［J］.Measurement，2012，45：30－40.

［3］沈 路，周曉軍，張志剛，等.Hilbert-Huang變換中的一種端點延拓方法［J］.振動與沖擊，2009，28（8）：168－174.

SHEN Lu，ZHOU Xiao-jun，ZHANG Zhi-gang，et al.Hilbert-Huang an endpointextensionmethod［J］.Journal of Vibration and Shock，2009，28（8）：168－174.

［4］于慧君，陳章位，王慶豐.一種加窗重疊信號平滑連接方法及其在振動信號預處理中的應用［J］.振動與沖擊，2007，26（8）：39－40

YU Hui-jun，CHEN Zhang-wei，WANG Qing-feng.A windowed overlapping signals smooth connection method and its application invibration signal preprocessing［J］.Journal of Vibration and Shock，2007，26（8）：39－40.

［5］韓建平，錢 炯，董小軍.采用鏡像延拓和RBF神經網絡處理EMD中端點效應［J］.振動、測試與診斷，2010，30（4）：414－417.

HAN Jian-ping，QIAN Jiong，DONG Xiao-jun.Mirror extension and RBF neural network processing end effect of EMD［J］.Journal of Vibration，Measurement＆Diagnosis，2010，30（4）：414－417.

［6］胡勁松，楊世錫.EMD方法基于徑向基神經網絡預測的數據延拓與應用［J］.機械強度，2007，29（6）：894－899.

HU Jin-song，YANG Shi-xi.EMD method is based on the extension and application of radial basis function neura network forecast data［J］.Journal of Mechanical Strength，2007，29（6）：894－899.

［7］黃先祥，李勝朝，謝 建.新型經驗模式分解端點效應消除方法［J］.機械工程學報，2008，44（9）：1－5.

HUANG Xian-xiang，LI Sheng-chao，XIE Jian.The new empirical mode decomposition endpoint effect elimination［J］.Chinese Journal of Mechanical Engineering，2008，44（9）：1－5.

［8］胡愛軍，安連鎖，唐貴基.Hilbert-Huang變換端點效應處理新方法［J］.機械工程學報，2008，44（4）：154－158.

HU Ai-jun，AN Lian-suo，TANG Gui-ji.Hilbert-Huang Transform endpoint effects processing new method［J］.Chinese Journal of Mechanical Engineering，2008，44（4）：154－158.

［9］程軍圣，于德介，楊 宇.基于支持矢量回歸機的Hilbert-Huang變換端點效應問題的處理方法［J］.機械工程學報，2006，42（4）：23－31

CHENG Jun-sheng YU De-jie YANG Yu.Processingmethod based on the Hilbert-Huang Transform，support vector regression machine end effects［J］.Chinese Journal of Mechanical Engineering，2006，42（4）：23－31.

［10］任達千，吳昭同，嚴拱標.EMD端點效應的評價指標及抑制邊界效應的窗函數法［J］.制造業自動化，2007，29（1）：21－24.

REN Da-qian，WU Zhao-tong，YAN Gong-biao.Window functionmethod of evaluation of the end effect of EMD and suppress boundary effects［J］.Manufacturing Automation，2007，29（1）：21－24.

［11］Qi K Y，He Z J，Zi Y Y.Cosine window-based boundary processing method for EMD and its application in rubbing fault diagnosis［J］.Mechanical Systems and Signal Processing，2007，21：2750－2760.

［12］Parey A，PachoriR B.Variable cosinewindowing of intrinsic mode functions：Application to gear fault diagnosis［J］.Measurement，2012，45（3）：415－426.

Improved algorithm against end effect of EMD based on extension and variable cosine w indow

DING Xue-juan，SHIPei-ming
（College of Electrical Engineering，Yanshan University，Qinhuangdao 066004，China）

For the end effect of empirical mode decomposition（EMD）keeping pace with increase in filtering level，a new improved method using a combination of extension and variable cosine window was proposed.Firstly，the signalwas extended with a certain length of data.Then，the extended signal was processed at both ends with variable cosine window before extracting every intrinsicmodel function（IMF）tomake the end effect of each IMFminimized and to ensure the correct decomposition and the improvement of the EMD algorithm.Simulation analysis and crack fault diagnosis examples showed that the proposedmethod can inhibit end effects effectively.

signal processing；EMD；end effect；variable cosine window

TH165；TN911

A

國家自然科學基金（51005196）；高等學校博士學科點專項科研基金（20101333120004）

2013－01－17 修改稿收到日期：2013－03－23

丁雪娟女，碩士生，1986年生