周 義, 李鴻光
(上海交通大學 機械系統與振動國家重點實驗室,上海 200240)
?
快速自適應經驗模態分解方法的基本原理及其性能評估
周義, 李鴻光
(上海交通大學 機械系統與振動國家重點實驗室,上海200240)
摘要:經驗模態分解是一種有效的信號分解方法,尤其是針對非平穩非線性信號。然而,隨著研究的深入,學者們發現該方法中存在著諸多弊端。根據Bhuiyan的研究,提出了一種針對一維信號的快速自適應經驗模態分解方法。通過大量的數值仿真,證明這種方法不但能克服傳統方法的弊端、得到高質量的分解結果,還能大幅度地提高計算效率。
關鍵詞:經驗模態分解;快速自適應經驗模態分解;數值仿真
傳統的信號處理方法大致可以分為兩類:時域統計分析和傅里葉變換。然而,這些方法都基于一個共同的假設,處理對象是平穩線性信號。然而,工程實際中遇到的都是非平穩、非線性的瞬態信號,因此如果仍使用傳統方法,往往會得到錯誤的結果。為此,學者們提出了許多針對非平穩非線性信號的時頻分析方法,如小波變換。但是小波變換需要人為地選擇基函數,只有在基函數的形狀與信號性質時,才能產生高幅值小波系數,這種非自適應性極大地限制了這種方法的發展。
Huang[1]提出了一種自適應信號處理方法,經驗模態分解(Empirical Mode Decomposition,EMD)。該方法基于信號本身的局部特征時間尺度,能夠將一個非平穩非線性信號分解成為一系列完備的、近正交的分量,稱為“固有模態函數”(Intrinsic Mode Function,IMF)。IMF作為EMD方法的基函數,是嵌入在原信號中的固有振蕩模態。它不需要提前人為決定,而是由信號本身在分解過程中決定。因此,經驗模態分解是從根本上脫離傅里葉變換,而從信號本身進行分析處理,具有完全的自適應性、無監督性。
經過十幾年的發展,EMD的理論研究逐漸深入,在一維信號處理方面的應用更為廣泛。Huang等主要建立了EMD的基本框架,分析了基本原理,比較了EMD算法和小波以及其他信號分析方法的區別。其他學者在理論框架的基礎上,成功地將EMD應用在地震物理學[2]、生物醫學[3]、設備故障診斷[4-5]等領域。
然而,隨著研究的深入,學者發現EMD存在諸多弊端:①終止準則:Rilling等[6]將終止準則中的SD,變為雙閾值準則。②端點效應:郭明威等[7]使用基于時間尺度的LS-SVM端點延拓方法,蔡艷平等[8]使用Lyapunov指數預測模型進行邊界延拓,時培明等[4]提出一種波形特征匹配延拓與余弦窗函數相結合的方法。③模態混疊:Wu等[9]提出集合經驗模態分解,消除模態混疊現象。④運算效率:Rilling等[10]又把單變量EMD擴展到雙變量EMD,減少計算量;Damerval等[11]提出基于Delauney三角和分段三次多項式插值方法提高EMD的分解速度。
盡管如此,還沒有一種算法,可以綜合考慮上述四個問題,只是針對某一個問題展開討論。因此,本文根據學者Bhuiyan的研究[12],提出了一種旨在提高效率的新算法——快速自適應經驗模態分解(Fast and Adaptive Empirical Mode Decomposition,FAEMD),同時指出此算法能在不同程度上優于目前現有的改進方法,能較好的同時解決上述幾個問題。
1傳統經驗模態分解
經驗模態分解,作為時頻分析方法的基礎,可以將一個非平穩非線性信號分解成多個被稱為“固有模態函數”的振蕩信號和趨勢項,所有IMF按頻率從高到底排列,最后一項為趨勢項。
固有模態函數是滿足單一振蕩模式以及單分量信號物理解釋的一類信號。單分量信號是相對于多分量信號而言的,只有單分量信號的瞬時頻率才有物理意義。Huang定義的固有模態函數必須滿足以下兩個條件[1]:①在整個數據范圍內,極值點和過零點的數目必須相等或最多相差一個;②在任何時刻,由極大值點定義的上包絡線和極小值定義的下包絡線的平均值為零。
Huang認為,任何信號,無論是平穩還是非平穩信號,都能通過不斷迭代分解出滿足條件的固有模態函數,這個迭代過程即為“篩分(Sifting)”。篩分具有兩個目的:一是消除模態波形的疊加;二是使波形輪廓更加對稱。每一輪篩分過程都將分離出1個固有模態函數,然后將余量繼續執行篩分操作,當全部固有模態函數提取出來后,就剩下一個余量或趨勢項,原始信號x(t)可表示為:
(1)
式中:ci(t)表示第i個IMF分量,r(t)表示篩選到最后剩下的趨勢項。
2快速自適應經驗模態分解
2.1基本原理
Bhuiyan[12]提出了快速自適應二維經驗模態分解,完全突破了傳統經驗模態分解的實現模式。基于此,我們做了大量改進工作,提出了針對一維信號的快速自適應經驗模態分解方法。通過使用順序統計濾器替代樣條插值函數獲取上下包絡線,并在保證分解質量的前提下,將控制迭代次數,不僅降低了程序的計算時間,結果在許多指標上也比傳統EMD出色。
順序統計濾波器(Order Statistics Filter, OSF)是一種非線性濾波器,最早用于圖像處理,增強圖像的平滑度。其思路是:以極值點間距的最小值作為濾波器寬度對信號進行順序統計濾波,求上下包絡;再用均值濾波器做平滑處理。
假設原信號含有n個極大值點,OSF統計每個極大值點與相鄰極大值點的最小距離,并存入數組dadj-max中,顯然該數組有n個數據;同理,對于極小值點,也能得到一個數組dadj-min。可以通過如下四種準則確定OSF窗口寬度w:
w=d1=min{min{dadj-max},min{dadj-min}}
(2)
w=d2=max{min{dadj-max},min{dadj-min}}
(3)
w=d3=min{max{dadj-max},max{dadj-min}}
(4)
w=d4=max{max{dadj-max},max{dadj-min}}
(5)
式中:max代表取數組中所有元素的最大值;min代表取數組中所有元素的最小值。分別定義式中四種準則為1型、2型、3型和4型,四個所得值之間的關系為:d1 確定OSF窗口寬度后,用寬度為w、平行于時間坐標軸、相切于極值點的線段,初步繪制包絡線,形成上下包絡線: (6) (7) 式中:i,j代表提取第i階IMF中的第j次迭代;uij(t)為上包絡,lij(t)為下包絡;hmaxij(t)代表極大值譜,hminij(t)代表極小值譜;zt是窗口寬度為w的線段。 在源信號上(下)方繪制了若干個線段,但這些線段之間有的沒有接壤,有的接壤卻不是平滑過渡,因此還需要進行均值平滑濾波。這種平滑濾波在方程所得結果上直接進行,推薦使用與次序統計濾波濾波器同樣的窗口寬度,具體操作如下所示: (8) (9) 式中:w′是均值平滑濾波器的窗口寬度,令w′=w。這一步平滑濾波操作是一種算術平均濾波,使數據內的局部變量更加平滑。 此外,還要設置具體的單次篩分迭代次數,一般1~3次就可以保證精度。其余步驟與傳統EMD一致。 以仿真信號為例,圖1為仿真信號通過FAEMD在獲取包絡均值的過程。顯然通過順序統計濾波器得到的包絡均值存在許多一階不可導點。所以我們在獲得包絡均值后,對包絡均值用算數平均的方式做平滑處理,獲得圖2所示光滑曲線。 圖1 仿真信號在獲取包絡均值的過程Fig.1 The process obtaining mean envelop for simulated signal 圖2 均值平滑濾波處理過程Fig.2 The process of average smoothing filtering 2.2窗口寬度設定 根據大量的仿真計算可知,式(2)~(5)中:1型窗口尺寸最小,4型窗口尺寸最大;使用3型和4型提取第i+1階IMF,窗口尺寸會大于第i階IMF;使用1型和2型提取第i+1階IMF,窗口尺寸有時可能并不會大于第i階IMF,這種現象有悖于經驗模態分解的基本原則,即后階IMF一定是含有頻率較低的成分。因此,有必要額外添加一個操作,令后面的窗口尺寸大于前面的,例如可以令當前窗口尺寸為前一個窗口尺寸的1.5倍。盡管這個操作可能沒有必要,但它可以充分保證經驗模態分解的基本原則。 一般情況下,建議使用統一的OSF窗口寬度繪制上下包絡。如果分別確定上下OSF窗口寬度,可能會遇到如下問題:①源信號中極小值點足夠多而極大值點只有一個,會導致數組dadj-max為空,無法繪制上包絡面,卻沒有滿足終止準則,不能終止篩分迭代的過程。②對于同一階IMF的提取,使用相同寬度的OSF,有助于從源信號提取尺度相近的成分;而如果寬度不一致,得到的IMF分量中尺度信息可能相差很遠,這就是 “模式混疊”現象。 2.3小結 綜上所述,FAEMD與原始EMD的最大不同之處在于: (1)傳統EMD方法均使用插值技術繪制包絡曲線,插值法本身就是一個復雜的迭代過程,需要大量的計算資源;而FAEMD突破傳統的插值法,使用估計法繪制包絡曲線,這種方法尤其適用于含噪信號,可以快速繪制包絡線。 (2)EMD方法采用篩分的方法提取IMF分量,每一次篩分都需要幾十次甚至上百次的循環迭代,極度耗時;而FAEMD在保證分解質量的前提下,把單次篩分的迭代次數降低到5次以內,甚至是1次,大大提高了計算速度。 3方法的性能對比實驗 雖然EMD發展了10余年,但仍有許多問題需要解決:①終止準則;②端點效應;③模式混疊;④運算效率。 3.1終止準則 過多的篩分會導致IMF分量變成純粹的頻率調制信號,而幅度趨于恒定。因而存在篩分過程的終止準則問題。為了保證IMF在幅值和頻率上都有物理定義,必須定義一個篩分過程的停止準則。這個準則可以通過限制標準差的大小來實現。標準差SD通過兩個連續的處理結果來計算得出: (10) 式中:T代表信號長度,hij(t)代表提取第i階IMF過程中的第j次迭代所得信號,hi(j-1)(t)代表提取第i階IMF過程中的第(j-1)次迭代所得信號。SD稱為“篩分閥值”,一般取0.2。如果計算所得SD若小于此值,篩分過程停止,認為為第一階IMF;否則,還要繼續循環迭代,直到滿足該條件為止。 ①課程內容按照學科知識點的邏輯順序組織。為了促進學生主動思維能力的發展,教師設定每個知識單元學習目標和完成標準,指引學生跟隨任務完成學習。如在教學三相異步電動機的正反轉控制時,教師事先提出預習思考題:異步電動機的轉動原理是什么?如何實現異步電動機的正轉和反轉?使學生帶著學習目標學習并思考。 終止準則會對分解質量和分解速度產生極大的影響:如果SD值過大,雖然分解速度很快,但分解精度變粗,每個IMF會含有過寬的頻率成分,而不再是局部窄帶的單分量信號;反之,如果SD值過小,雖然能夠保證分解精度,但運算時間過長。對于不同工況下的工程信號,找到其合適的SD值是比較困難的,往往需要多次使用不同的SD值,人工地尋找合理的分解結果,但這種“試探”的過程,與小波分解中基函數的尋找相類似,都破壞了方法的自適應性。 而本文提出的FAEMD方法,使用順序統計濾波器繪制包絡線,窗口寬度w的產生是基于信號極值點。經過大量的實驗可知,在w尺度下,無須迭代或較少次迭代即可獲得一個IMF,直接避免了過度篩分導致的幅值同化。 3.2端點效應 在使用EMD處理信號時,端點效應是常遇到的一個問題,特別是在處理長度較短的信號時,端點效應會成為影響分析結果精度的原因之一。 EMD的一個關鍵步驟是根據極值點作三次樣條插值獲得包絡線。除非信號兩端均為極值點,否則不能簡單地認為端點就是極值點,這將導致插值產生的包絡線存在擬合誤差。由于EMD的篩分過程需要反復求取信號的包絡線,而端點處極值的不確定性導致每次求包絡時都存在擬合誤差。隨著誤差不斷累積,第一階IMF的端點處就會存在較大誤差。而第二階IMF的分解是建立在原信號減去第一階IMF的殘余項的基礎上進行的。所以第一階IMF中的誤差會傳遞給第二階IMF,導致第二階IMF有更大的誤差。以此類推,隨著分解的進行,誤差會不斷的積累并從端點向內部傳遞。嚴重時會使分解后的IMF失去意義。 構造一個簡單的信號函數如下: (11) 使用傳統的EMD方法處理仿真信號,分解結果如圖4所示。可以看到,分解結果,不能如實反映原信號的兩個單分量成分,而是出現許多不能解其物理含義的分量。這就是因為,邊界效應的影響不但在邊界處,還會隨著分解的進行,深入到信號內部,最終導致分解失敗。 圖3 仿真信號及其真實分量Fig.3 The simulated signal and its real components 圖4 仿真信號的傳統EMD分解結果Fig.4 The decomposition results by applying traditional EMD to the simulated signal in Fig. 3 因此,使用本文提出的FAEMD方法,處理圖3所示仿真信號,所得結果如圖5所示。從結果可以看出,FAEMD很好地從原信號中分離出兩個單分量成分。 圖5 仿真信號的FAEMD分解結果Fig.5 The decomposition results by applying FAEMD to the simulated signal in Fig. 3 3.3模式混疊 在進行信號的篩分過程中,有時會出現某個IMF分量中包含不同特征時間尺度或相近時間尺度的成分分步在不同的IMF中,這種情況稱為模式混疊。 圖6 含有間歇性分量的仿真信號Fig.6 The silumated signal with intermittency compoents 模式混疊歸因于信號的間歇現象,它的產生和極值點的選擇有關。圖6所示的仿真信號,由正弦波s1、同周期發生的高頻三角脈沖s2和趨勢項s3組成,這是一種簡單卻典型的可產生模式混疊的信號。由于信號通常含有間歇型分量、脈沖干擾或噪聲,使得基礎成分(頻率較低)的某個局部出現高頻振蕩。這種局部振蕩和基礎成分各自的包絡線會很不相同,因此在繪制整體包絡時,特別是在高低頻成分的交界處會造成混亂現象,也就是模式混疊。它不僅能扭曲信號的時頻分布,還能掩蓋各IMF分量蘊含的真實物理意義。 為了解決模式混疊問題,Wu等提出了一種噪聲輔助的數據分析方法——集合經驗模態分解(Ensemble Empirical Mode Decomposition,EEMD)。 EEMD的思路是通過對信號添加高斯白噪聲使得間歇的信號連續化,再進行EMD分解,最后對多次EMD分解結果做集合平均。由于白噪聲隨機產生,因此根據其統計特性,可以通過足夠多次的試驗抵消。 EEMD在使用時有兩個參數需要設置,即高斯白噪聲的幅值k和算法執行EMD的總次數M: (1)k越小,越有利于分解精度的提高,但當小到一定程度時,可能不足以引起信號局部極值點的變化,就不能解決模式混疊問題,k太大則會淹沒原始的信號特征,使得分解失去意義。建議k取原信號幅值標準差的0.01倍~0.5倍。 (2)M增大,顯然噪聲影響會減少,但M過大,會導致運算時間過長。當M取100~200次時,殘留噪聲引起的誤差一般不足1%。 圖7為EMD算法的分解結果,得到IMF分量C1-C4及趨勢項。只有C4和趨勢項具有真實的物理意義,分別代表了原始信號中的正弦分量和趨勢項。C1出現嚴重的失真,多出C2和C3這兩個無法解讀其含義的成分。C1和C2都含有高頻和低頻成分,可以認為它們都是s1和s2的模式混疊。由于傳統EMD無法將原信號中的間歇成分s2正確的分解出來,導致C1夾雜其他模態成分產生嚴重的失真,且剩余的間歇組分需要在后續的IMF中不斷的被部分分解出來,使得信號中多出了C2和C3這種無意義的成分。由于C4幅值比前3個IMF大得多,所以受到的干擾較小。 圖7 運用傳統EMD算法分解結果Fig.7 The decomposition results by applying traditional EMD to the simulated signal in Fig. 6. 圖8為FAEMD算法分解結果,得到IMF分量C1-C4及趨勢項。res依舊保持了原有的趨勢,C4較好的表現出S1的特征。前三項均反映出了S2的特征,有效的分解出了間歇信號。從能量的角度看,C1、C2和C3幅值依次減小,C2、C3可以分別認為是高頻信號在上次分解中產生的能量泄漏。它們都代表了C1中損失的能量,這解釋了C1中幅值小于3的原因。 圖8 運用FAEMD算法分解結果Fig.8 The decomposition results by applying FAEMD to the simulated signal in Fig. 6. FAEMD方法確實能高質量地分離出信號中的各個成分,有效抑制模態混疊的發生,客觀展現了信號的真實物理內涵。之所以可以有效地避免間歇信號與低頻信號混淆,是因為當高頻分量存在時,窗口寬度較小,此時濾波法可以較準確的刻畫低頻信號從而對其保留,這樣低頻的信號可以留在下次分解 使用EEMD方法,雖然也可以消除模式混疊現象,得到類似于FAEMD方法的結果,但相比比下,EEMD 不足之處明顯: (1)EEMD運算時間較長,由于EEMD中包含了多個EMD循環,運算時間加長,處理數據效率低下。 (2)每個循環中引入的噪聲不同,得到IMF數量可能不同,造成最終分析結果可能包含更多的IMF項,分散了真實分量應有的能量。 (3)由于EEMD所得的每個IMF是由M個EMD循環產生的IMF求平均所得,導致分析結果可能不能嚴格滿足原本IMF定義。FAEMD不會有噪聲引入的誤差問題,也不會因為噪聲幅值過大導致信號中微小的特征被淹沒。 3.4運算效率 為解決EMD方法中存在的諸多弊端,學者們先后提出了很多技術。雖然這些弊端均得以解決,但卻派生出另一個問題:運算時間大幅上漲。①EMD方法的篩分思想,需要不斷地進行迭代和插值,這需要耗費大量的計算資源;②為了避免端點效應而采用的邊界延拓技術增加了需要處理的數據量;③為了消除模態混疊效應而采用的集合經驗模態分解法,需要數百次的EMD循環,因而增加數百倍運算時間。 對于任何一個信號都可以將其分解成為若干個固有模態分解和一個殘余分量,其中每個固有模態分量包含了信號的各種不同頻段成分,且它們在頻域近似相互正交。為了評估EMD分解結果的優劣,文獻[1]引入了“正交因子(Orthogonality Index,OI)”。在工程實際中的分解,OI值為0是幾乎不可能的,只能通過改良EMD算法來獲取盡可能低的OI值,OI值越低,代表分解的效果越來。一般來說,低于0.1的OI值是可以接受的。 我們以實驗中測得的旋轉機械松動故障實驗信號作為測試數據,分別用EMD,EEMD和FAEMD三個MATLAB程序進行處理,得計算時間和正交因子。從表1中,可以清楚看到,FAEMD方法不僅具有最少的計算時間,而且得到了最低的正交因子,也就是最高質量的分解結果。 表1 傳統EMDs與FAEMD計算時間及正交因子 4結論 FAEMD算法不僅能夠顯著提高原始EMD的計算速度,同時能很好的避免EMD的四個主要問題:終止準則、邊界效應、模態混疊和運算效率。由于FAEMD自適應的IMF分解過程,可以有效避免原始EMD中IMF篩分過度所導致的幅值無意義問題,并無須考慮篩分的終止條件,實現了真正的自適應。由于FAEMD采用的全新包絡均值求解方法,信號在邊界處不會出現過沖或欠沖的擬合失真問題,避免了嚴重的端點效應。另外FAEMD在處理容易出現模態混疊現象的信號時,分解結果明顯優于現有最好的EEMD算法。 參 考 文 獻 [1] Huang N E, Shen Z, Long S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J]. Proceedings of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, 1998,454(1971): 903-995. [2] 陳林,宋海斌. 基于經驗模態分解的地震瞬時屬性提取 [J]. 地球物理學進展,2008,23(4): 1179-1185. CHEN Lin, SONG Hai-bin. Seismic instantaneous attribute extraction based on empirical mode decomposition [J]. Progress in Geophysics, 2008, 23(4): 1179-1185. [3] 李永勤,王清,鄧親愷. EMD及其在生物醫學信號處理中的應用研究[J]. 生物醫學工程學雜志, 2005, 22(5): 1058-1062. LI Yong-qin, WANG Qing, DENG Qin-kai. Research on EMD and its application in biomedical signal processing [J]. Journal of Biomedical Engineering, 2005, 22(5): 1058-1062. [4] 時培明, 丁雪娟, 李庚,等. 一種EMD改進方法及其在旋轉機械故障診斷中的應用[J]. 振動與沖擊, 2013, 32(4): 185-190. SHI Pei-ming, DING Xue-jun, LI Geng, et al. An improved method of EMD and its applications in rotating machinery fault diagnosis [J]. Journal of Vibration and Shock, 2013, 32(4): 185-190. [5] 蘇文勝, 王奉濤, 張志新, 等. EMD 降噪和譜峭度法在滾動軸承早期故障診斷中的應用 [J]. 振動與沖擊, 2010, 29(3): 18-21. SU Wen-sheng, WANG Feng-tao, ZHANG Zhi-xin, et al. Application of EMD denoising and spectral kurtosis in early fault diagnosis of rolling element bearings [J]. Journal of Vibration and Shock, 2010, 29(3): 18-21. [6] Rilling G, Flandrin P, One or two frequencies? the empirical mode decomposition answers[J]. IEEE Transactions on Signal Processing, 2008, 56(1): 85-95. [7] 郭明威, 倪世宏, 朱家海, 等. 振動信號中HHT/EMD端點延拓方法研究[J]. 振動與沖擊, 2012, 31(8): 62-66. GUO Ming-wei, NI Shi-hong, ZHU Jia-hai, et al. HHT/EMD end extension method in vibration signal analysis [J]. Journal of Vibration and Shock, 2012, 31(8): 62-66. [8] 蔡艷平, 李艾華, 石林鎖,等. EMD端點效應的改進型混沌延拓方法及其在機械故障診斷中的應用[J]. 振動與沖擊, 2011, 30(11): 46-52. CAI Yan-ping, LI Ai-hua, SHI Lin-suo, et al. Processing method for end effects of EMD based on improved chaos forecasting model and its application in machinery fault diagnosis [J]. Journal of Vibration and Shock, 2011, 30(11): 46-52. [9] Wu Z, Huang N E. Ensemble empirical mode decomposition: a noise-assisted data analysis method[J]. Advances in Adaptive Data Analysis, 2009,1(1): 1-41. [10] Rilling G, Flandrin P, Gonalves P, et al. Bivariate empirical mode decomposition[J]. Signal Processing Letters, IEEE, 2007, 14(12): 936-939. [11] Damerval C, Meignen S, Perrier V. A fast algorithm for bidimensional EMD. Signal Processing Letters[J]. IEEE, 2005,12(10): 701-704. [12] Sharif B, Adhami R R, Khan J F. Fast and adaptive bidimensional empirical mode decomposition using order-statistics filter based envelope estimation[J]. EURASIP Journal on Advances in Signal Processing,2008,2008:728356. Basic principle of a fast and adaptive empirical mode decomposition and its performance evaluation ZHOUYi,LIHong-guang (State Key Laboratory of Mechanical System and Vibraiton, Shanghai Jiaotong University, Shanghai 200240, China) Abstract:Empirical mode decomposition is an effective signal decomposition method, especially for non-stationary and non-linear signals. But it is found that there exist a lot of drawbacks along with the deepening of studies. Therefore, a one-dimensional signal processing method, called fast and adaptive empirical mode decomposition, was proposed according to the Bhuiyan’s study. It has been proved by numerous numerical simulations that the method can not only overcome the drawbacks of traditional methods to obtain high quality decomposition results but also enhance the computing efficiency. Key words:empirical mode decomposition; fast and adaptive empirical mode decomposition; numerical simulation 中圖分類號:Tp12;Tp13.3 文獻標志碼:A DOI:10.13465/j.cnki.jvs.2016.03.003 通信作者李鴻光 男,博士,教授,博士生導師,1972年生 收稿日期:2013-12-19修改稿收到日期:2014-01-28 基金項目:國家自然科學基金資助項目(11372176) 第一作者 周義 男,博士生,1985年生







