基于改進EMD方法的系統誤差分離方法研究

姚 剛，趙建軍，姚躍亭，任 喜

(海軍航空工程學院, 山東 煙臺 264001)

誤差為測量值與真值之差。測量系統中所有測量結果均含誤差，其大小直接影響測量系統精度。測量誤差分為隨機誤差與系統誤差。隨機誤差具有不確定性，其發生、變化規律雖無法確定，但可用統計方法研究。而系統誤差按確定性規律變化，其正負、大小可預知。為研究、提高測量系統精度，需對系統誤差深入分析。系統誤差分離為研究重要內容之一。對復雜系統，其誤差具有模型難識別、對處理結果影響難估計的復雜特性。

系統誤差一般可分為常值誤差、線性漂移誤差、周期性誤差及復雜規律變化誤差等。系統誤差分離可視為對誤差信號中各頻帶分離過程。測量誤差分離方法現使用較多的有最小二乘回歸法、中值濾波法、傅里葉級數逼近法、小波及小波包分解法等[1]。各種方法的使用均存在限制，即經典最小二乘回歸法的模型階數較難確定；中值濾波方法精度難以提高；傅里葉級數逼近法為對誤差進行諧波分解；小波基函數選取會直接影響小波及小波包分解法分離質量，增加使用難度。基于上述原因，本文提出用改進EMD方法對系統誤差分析、處理，并進行仿真及實例驗證。

1 EMD基本理論

1.1 EMD方法概述

Huang等[2-3]在瞬時頻率概念基礎上提出的信號處理方法-基于EMD(Empirical Mode Decomposition)的時頻分析方法可處理非線性、非平穩信號，分解時基函數不確定，只依賴數據本身，具有良好的自適應特點。EMD方法將數據視為由一個或多個IMF(Intrinsic Mode Functions，固有模態函數)組成。由于經驗模態分解方法的優點，在生物、海洋、地球科學、天文學及工程技術等領域中廣泛應用。實際的系統誤差往往非平穩、非線性，EMD為處理此類復雜多變的誤差提供了新的方法。

經EMD方法分解所得每個IMF，具有的性質為：① 在整個信號上，極值點個數與過零點數目相等或最多相差一個；② 局部均值為零，即由局部極大值構成的上包絡線與由局部極小值構成的下包絡線的平均值為零。實際中，上下包絡均值無法為零，通常認為包絡均值滿足IMF均值為零的條件為

(1)

1.2 EMD分解步驟

EMD算法假設對于何信號均由若干有限的IMF組成，每個IMF通過方法獲得[4-5]。設原信號為x(t)，找出整個信號x(t)的極值點，并通過三次樣條擬合出信號的上下包絡線，分別為emax(t),emin(t)，計算上下包絡線均值為

(2)

將原始信號減去包絡線均值，獲得新信號為

(3)

重復上述步驟，依次獲得imf2,…,imfn，直至所得余項rn為單調信號或小于預先給定閾值，分解結束。

原函數經EMD分解得：

(4)

2 基于波形相似度的端點延拓方法

2.1 端點效應產生原因

對誤差信號進行EMD分解時，需解決端點效應問題。端點效應產生原因即在對信號求上下包絡線時需對極值點進行擬合。若直接將信號端點作為信號的極值點，可使擬合所得包絡覆蓋整個信號。雖滿足經驗模態分解要求條件，但會導致端點附近包絡出現較大誤差，且隨篩分過程的不斷反復，邊界效應產生的誤差會向信號內部傳播，影響整個分解過程。對較長數據，可采取在分解過程中通過拋棄兩個端點附近數據方法削弱端點效應。而對短數據而言，只能采取端點延拓方式抑制。

2.2 改進的端點延拓方法

目前已有抑制端點效應的方法，包括鏡像延拓法、基于神經網絡的延拓法、混沌延拓、基于多項式擬合的延拓方法、基于極值時間尺度的LS-SVM延拓方法等[6-10]。雖對端點效應抑制有一定效果，但亦存在各自問題。鏡像延拓[6]方法處理短數據時效果較差，因其會截去部分數據。神經網絡[7]及混沌方法[8]的不足在于速度太慢。而多項式擬合方法[9]適應性較差，對某些數據結果較好，對另一些數據效果會不理想。基于LS-SVM的端點延拓方法[10]易受內積函數、核函數、精度參數及懲罰參數經驗選擇影響。為此，本文提出基于波形相似度的端點延拓方法。

端點效應抑制原理為尋找信號之外離端點最近的極值點，與信號本身的極值點一起進行樣條插值，構成信號的上下包絡線，再截取所需分解數據。本質上，延拓即為通過各種數據處理方法，對端點處信號變化趨勢進行預測。基于波形相似度的端點延拓方法基本思想即能在信號內部找到一段波形與端點處的波形相似，據相似程度建立數據關系對原始信號進行延拓。可最大限度維護端點處波形的變化趨勢，從而能有效抑制端點效應的產生。

2.3 波形相似度定義

波形相似度為表明兩波形相似程度的度量。對信號x(t)，y(t),(t=1,2…n)，定義兩波形相似度為

(5)

式中：cov(·)為x(t)與y(t)協方差；σ(x)，σ(y)分別為x(t)與y(t)的方差。

定理1：當且僅當隨機變量Y與X之間存在線性關系Y=a+bX時，相關系數|ρ(X,Y)|=1。

據定理1知：兩波形相似度較高時(ρ≥0.8)，可近似將兩者視為線性關系，用線性回歸方法尋找兩者間線性關系。

2.4 延拓步驟

圖1 待匹配波示意圖

對原始信號延拓時，左右兩端方法相同，以左端點為例說明延拓步驟：

(1) 找出待匹配的波，并確定其長度。設待匹配的波為ω，起點為左端點值，向右包括極大值與極小值。當第一個極大值距左端點距離Lmax1較第一個極小值距左端點距離Lmin1小時，則待匹配波ω的終點為第一個極小值點min1，反之則為第一個極大值點max1。待匹配波ω終點為min1情況見圖1。

(2) 計算整個數據段中各波段與ω的相似度ρ。為提高運算速度，僅計算極值點附近波的相似度。

(3) 找出相似度最高的子波段ω_match，即相似度|ρ|最大的波。

(4) 用線性回歸方法計算兩波段數據關系，即

Y=a+bX

(6)

式中：X為ω幅值向量，Y為ω_match幅值向量。

(5) 分別向左延拓一個極大值(tmax0,vmax0)與一個極小值(tmin0,vmin0)。設在ω_match左端最近的極大值、極小值點為(tmaxp,vmaxp)與(tminp,vminp)。延拓的極大值max0與極小值min0的幅值據式(6)計算：

vmax0=(maxp-a)/b

(7)

vmin0=(minp-a)/b

(8)

延拓極值對應時刻由極值距端點距離確定。設(tmaxp,maxp)，(tminp,minp)距離ω_match左端點距離分別為Lmaxp，Lminp，則延拓的極大值、極小值對應時刻分別為

tmax0=-LmaxP

(9)

tmin0=-LminP

(10)

3 改進的EMD分解停止準則

分解停止準則為決定EMD分解停止的標準。終止條件太苛刻時會產生過多IMF分量，會干擾對信號的正確判斷；終止條件寬松時則會使分解不夠細致，無法識別信號的變化規律。對整個EMD分解過程停止需有適當的判斷依據。通常采用兩個判斷標準，即當最后一個IMF或殘余分量幅值小于預設值時，整個分解過程停止；或當殘余分量變成單調函數或常數時，不能再篩分出IMF時停止。用主元分析方法(PCA)的累計方差貢獻率對分解停止時機確定[11]時需對IMF及殘余分量進行主元分析。即將原始變量投影到新坐標系，使新主元互不相關，從而可用較少維度代表原始變量特征。而用EMD方法所得IMF已具有近似正交性及不相關性，因此可直接進行累計方差貢獻率計算。

3.1 IMF近似正交性

由EMD分解過程知，正交性在實際意義上是滿足。但理論上尚未給出嚴格的數學證明。通過EMD分解所得各IMF在局部均應相互正交，原因為每個IMF均由原信號與其極大、極小值包絡的局部均值之差獲得，故有

(11)

式(11)嚴格意義上講不準確，因均值由信號極大值、極小值點由曲線擬合求得，與真實均值存在一定誤差，雖泄露不可避免，但較小。IMF正交性可通過后驗方法給出。原信號重構信號可表示為

(12)

將殘余分量rn(t)視為第n+1個IMF分量，對上式兩邊平方得：

(13)

因此對原信號x(t)的正交性指標(Index of Orthogonality)可定義[2]為

(14)

式中：T為信號總長度。

若分解是正交的，則式(14)中交叉項為0。通過大量試驗驗證[2-3]得出，對一般信號其正交性指標通常不超過1%，對短信號在極端情況下可能達5%，即EMD分解所得各IMF為近似正交。

3.2 基于累計方差貢獻率的分解停止準則

方差貢獻率概念源于主元分析方法。方差貢獻率及累計方差貢獻率為衡量主成分重要性指標。在改進的EMD方法中，定義第i個IMF方差貢獻率為

(15)

該值越大，說明第i個IMF所含信息量越多。其中σ(·)為信號方差。

累計方差貢獻率(Cumulative Percent Variance，CPV)為前k個IMF方差貢獻率之和。實際應用時需大于一個期望閾值，從而可更好反映變量波動，即

(16)

式中：CL為確定閾值。

據各IMF分量的近似正交性，每個IMF分量攜帶信息量的大小可用方差貢獻率衡量。由累計方差貢獻率確定分解停止時機。計算每次分解IMF分量imf1,imf2,…imfn及殘余分量rn的方差貢獻率，求得每個分量及殘余分量的貢獻率。當殘余分量貢獻率小于給定閾值時則停止分解。分量方差貢獻率小于4×10-3時不必再分解[11]。因此本文采用兩個分解停止判斷標準，即累積方差貢獻率大于設定的閾值時，整個分解過程停止；或殘余分量變成單調函數或常數時，不能再篩分出IMF時停止。

4 仿真實例

4.1 多種類型系統誤差疊加

系統誤差表現形式可為單一常值誤差、線性漂移誤差或周期性誤差，也可為其組合。本文以兩種形式誤差疊加說明改進EMD方法的有效性。設系統誤差有三個誤差分量構成，其中y1(t)為線性誤差(常值誤差可視為線性誤差的特殊形式)，y2(t),y3(t)為周期性誤差。各誤差數學表達式為

y(t)=y1(t)+y2(t)+y3(t)

(17)

式中：y1(t)=ax+b,a=0.5,b=1;y2(t)=2sin(2πt);y3(t)=sin(10πt),(t=0:0.001:3)

為驗證改進EMD方法的有效性，分別對上述信號直接進行EMD分解、鏡像延拓、本文方法分解。采用基于累積方差貢獻率的停止準則，所得結果見圖2、圖3、圖4(其中實線為IMF分量，虛線為原始誤差分量)。

由三圖及表1看出，基于累積方差貢獻的停止準則下，三種方式均能獲得兩個IMF分量與一個殘余分量。未對端點進行延拓時的EMD分解所得分量與原始分量誤差較小，三個分量中線性誤差最大，此因端點效應向里蔓延之結果；基于鏡像延拓的分解方法，由于鏡像點選擇不合理原因產生一虛假分量，誤差較大，尤其線性誤差淹沒其中；基于波形相似度的端點延拓方法所得三個分量在端點處分解效果較好，與各原始信號誤差最小。

圖2 EMD分解

表1 各分解方法均方差

4.2 隨機誤差與系統誤差共存

實際所測系統中，誤差可由測量值減去真值獲得，誤差成分必含隨機誤差與系統誤差。進行系統誤差分離時，需對隨機誤差進行處理。本節仿真信號在式(17)基礎上加入N(0,0.62) 的正態白噪聲獲得。仿真時分別對信號進行直接分解及降噪后分解，實驗結果見圖5、圖6。圖5為未對隨機誤差處理的結果，分解時因隨機誤差頻率較高，產生較多的IMF高頻分量，使分解過程受到干擾，從而未很好分離出系統誤差。圖6先對原始信號進行隨機誤差處理，再進行基于波形相似度的延拓分解，不但較完整的分離出各項系統誤差，且端點效應亦得到較好抑制；但由于隨機誤差干擾，端點效應未能完全消除。

圖5 直接進行EMD分解

圖6 降噪后基于波形相似度延拓方法分解

5 實例分析

本文用動態測角儀誤差數據[12]驗證改進EMD方法的分離效果。考慮采樣點數較少，對誤差數據先進行三次樣條插值，所得誤差曲線見圖7。利用改進EMD方法分解所得結果見圖8。由二圖中看出，原始系統誤差包括線性漂移誤差與周期性誤差。經改進EMD方法分離得出4個誤差成分，其中imf1,imf2,imf3為具有周期性誤差諧波，所得殘余分量res為線性漂移誤差。系統誤差得到較好分離，可為測量儀器的誤差分析提供數據支持。

圖7 動態測角儀誤差圖

圖8 改進EMD方法分解

6 結 論

由于EMD方法為新的數據處理方法，不但可處理線性、平穩數據，且可自適應處理非線性及非平穩數據。本文將改進的EMD方法用于系統誤差分離過程。仿真實驗、實例分析表明，基于波形相似度的端點延拓方法可有效抑制端點效應，累積方差貢獻率可適當確定分解停止時機，從而較好完成系統誤差的分離，具有一定應用前景。

[1] 李世平，付宇，張進.一種基于EMD的系統誤差分離方法[J].中國測試，2011, 37(3): 9-13.

LI Shi-ping，FU Yu，ZHANG Jin.Systematic error separation method based on EMD[J].China Measurement & Test, 2011，37(3): 9-13.

[2] Huang N E, Shen Z, Long S R, et al.The empirical mode decomposition method and the Hilbert spectrum for non-stationary time series analysis[C].Proc.R.Soc.London, U.K., 1998: 903-995.

[3] Huang N E.Review of empirical mode decomposition[C].Proceedings of SPIE.Orlando, USA, 2001: 71-80.

[4] Rilling G, Flandrin P,Goncalves P.On empirical mode decomposition and its algorithms[C].IEEE-EURASIP Workshop on Nonlinear Signal and Image Processing, Grado, Italy, 2003:811-815.

[5] Zhang J,Yan R Q, Gao R X, et al.Performance enhancement of ensemble empirical mode decomposition[J].Mechanical Systems and Signal Processing, 2010(24):2104-2123.

[6] 王婷.EMD算法研究及其在信號去噪中的應用[D].哈爾濱:哈爾濱工程大學,2010.

[7] 邵晨曦，王劍，范金鋒，等.一種自適應的EMD端點延拓方法[J].電子學報,2007,35(10):1944-1948.

SHAO Chen-xi, WANG Jian, FAN Jin-feng, et al.A self- adaptive method dealing with the end issue of EMD[J].Acta Electronica Sinica, 2007, 35(10):1944-1948.

[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] 許寶杰,張建民,徐小力,等.抑制EMD端點效應方法的研究[J].北京理工大學學報，2006, 26(3):196-200.

XU Bao-jie, ZHANG Jian-min, XU Xiao-li, et al.A study on the method of restraining the ending effect of empirical mode decomposition[J].Transactions of Beijing Institute of Technology, 2006, 26(3): 196-200.

[10] 郭明威,倪世宏,朱家海，等.振動信號中的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.

[11] 李雪耀，鄒曉杰，張汝波，等.譜熵和主成分分析用于EMD分解研究[J].哈爾濱工程大學學報,2009,30(7):797-803.

LI Xue-yao, ZOU Xiao-jie, ZHANG Ru-bo, et al.Research on empirical mode decomposition based on spectrum entropy methods and principal component analysis[J].Journal of Harbin Engineering University, 2009, 30(7):797-803.

[12] 費業泰，盧榮勝.動態測量誤差修正原理與技術[M].北京：中國計量出版社,2001.