鄭佳佳，郭濱

（長(zhǎng)春理工大學(xué) 電子信息工程學(xué)院，長(zhǎng)春 130021）

把大腦皮層或頭皮表面的微弱的生物電通過腦電圖描記儀畫成曲線圖，這個(gè)曲線圖反映了腦電的活動(dòng)情況，也叫做腦電圖［1］（electroencephalogram，EEG）。腦電檢測(cè)是遵照國際10-20電極放置系統(tǒng)，把腦電測(cè)試儀放于頭部表皮上來察看大腦活動(dòng)情況［2］的檢測(cè)方法。

由于各種干擾信號(hào)會(huì)影響腦電信號(hào)的采集，包括眼電、心電、工頻干擾等，為了更好的通過腦電信號(hào)來分析病情，在臨床醫(yī)學(xué)診斷中，給采集到的腦電信號(hào)進(jìn)行去噪是必不可少的。腦電信號(hào)去噪處理方法有時(shí)域方法中的W-igner分布、小波包分析、小波分析、HHT分析，以及獨(dú)立分量分析等。生物醫(yī)學(xué)領(lǐng)域常用具有多分辨的小波變換，并且時(shí)頻局部化特性較好，但其存在頻率分辨率較差的高頻段和時(shí)間分辨率較差的低頻段等不足，在小波分解過程中，為了解決以上的問題，本文應(yīng)用了小波包分析。

黃鍔提出一種自適應(yīng)時(shí)頻分析的方法［3-5］。該方法又稱希爾伯特-黃變換，它的重要環(huán)節(jié)是經(jīng)驗(yàn)?zāi)B(tài)分解，自適應(yīng)地分解了信號(hào)中的局部變化特征，然而，如果極值點(diǎn)不均勻的出現(xiàn)在信號(hào)中，分解結(jié)果圖中會(huì)出現(xiàn)一個(gè)IMF分量上的多尺度信號(hào)，這就是模態(tài)混疊現(xiàn)象。這樣分解得到的信號(hào)就不是所需要的，不利于接下來的研究。故使用可以有效抑制模態(tài)混疊現(xiàn)象的聚合經(jīng)驗(yàn)?zāi)B(tài)分解（EEMD），之后對(duì)含噪聲的IMF進(jìn)行Hilbert變換，剖析Hilbert譜。將含噪的IMF分量進(jìn)行小波包分析，最后疊加各分量，得到去噪后腦電信號(hào)。

1 基于EEMD算法的HHT方法

白噪聲具有零均值，且頻率譜均勻分布的特征，如果把腦電加入到不同的白噪聲背景上，最后對(duì)結(jié)果求取總體均值，就可以消除噪聲，這就是EEMD的改進(jìn)部分。

EEMD算法步驟如下：

（1）給目標(biāo)信號(hào)S(t)加入一組高斯白噪聲n1(t)，獲得一個(gè)總體S1(t)；

（2）對(duì)S1(t)進(jìn)行EMD分解，得到n個(gè)imf分量，即：

（3）重復(fù)前二個(gè)步驟，每一次都加入新的白噪聲；

（4）得到各自的imf分量組；

（5）將相應(yīng)的imf總體均值作為結(jié)果：

其中，imfi(t)為第i個(gè)imf分量的總體均值；T為總的平均數(shù)；imfij(t)是白噪聲加入第j次后得到的imf；

（6）將相應(yīng)的剩余R的總體均值作為結(jié)果：

EEMD算法流程如圖1所示。

圖1EEMD算法

2 小波包分析

小波分析是基于正交小波基函數(shù)，將信號(hào)先處理成頻率較高和頻率較低的兩部分，之后僅對(duì)頻率低的進(jìn)一步分解，不再分析頻率高的，所以小波變換不能分解和表示高頻段含有有用信息的信號(hào)，如生物醫(yī)學(xué)信號(hào)等。而小波包分析不僅對(duì)低頻段分解，而且還分析高頻段，因此可以對(duì)中高頻段含有的用信息的信號(hào)進(jìn)行時(shí)頻局部分析。小波包是小波分析基礎(chǔ)結(jié)合最優(yōu)基選擇的結(jié)果。根據(jù)各個(gè)被分解的頻帶段信號(hào)特征，自適應(yīng)的選擇最優(yōu)基函數(shù)去匹配信號(hào)，這樣信號(hào)的分析就能提高，因此，小波包分析被廣泛的應(yīng)用到該領(lǐng)域中［6，7］。

2.1 小波包分解

在小波多分辨分析理論的基礎(chǔ)上，將尺度函數(shù)φ(t)記為u0(t)，小波函數(shù)φ(t)記為u1(t)，則二尺度方程：

式中，u0(t)=φ(t)所對(duì)應(yīng)的小波包為{un(t)}n∈Z，而H(k)，G(k)是共軛濾波器系數(shù)［8］。

假設(shè)S(t)是時(shí)間信號(hào)，pij(t)表示第j層上的第i個(gè)小波包，則二進(jìn)小波包分解為：

2.2 小波包重構(gòu)

式中，j=J-1,J-2,...,1,0;i=2j,2j-1,...,2,1;J=logN2，H是與尺度函數(shù)相關(guān)的一個(gè)小波包重構(gòu)濾波器，G為與小波函數(shù)相關(guān)的另一個(gè)小波包重構(gòu)濾波器。

小波包過程如圖2所示。

圖2 小波包分解過程

經(jīng)過小波包分解后，S代表每一層被分解的高頻和低頻部分，其良好的時(shí)頻分辨率正適合非平穩(wěn)信號(hào)EEG的頻域分析［9］。

3 消噪過程

帶有噪聲的腦電信號(hào)模型：

其中，S(t)是帶噪信號(hào)；f(t)為純凈腦電信號(hào)；e(t)為噪聲；ε為噪聲強(qiáng)度。

基于改進(jìn)的EMD結(jié)合小波包去噪的主要內(nèi)容為：含噪的信號(hào)經(jīng)過EEMD分解得到從高頻到低頻率的IMF分量，經(jīng)觀察可得噪聲主要分布在高頻部分，同時(shí)經(jīng)過Hilbert譜分析，可以確定高頻段含有噪聲（前三個(gè)IMF分量中），不處理第4個(gè)到第9個(gè)的IMF分量，只需要通過小波包對(duì)含噪聲的高頻段IMF進(jìn)行去噪處理，最后把新的前三個(gè)IMF分量和后面的6個(gè)原始分量相加，即為處理過噪聲的腦電信號(hào)。

具體步驟如下所示：

（1）對(duì)帶噪腦電信號(hào)S(t)進(jìn)行EEMD分解，得到n個(gè)IMF分量；

（2）計(jì)算每個(gè)IMF分量的Hilbert譜：

（3）通過小波包分析對(duì)前三個(gè)IMF分量進(jìn)行處理；

（4）把留下來的低頻部分IMF與經(jīng)過小波包處理的高頻部分IMF相加起來，可以得到去過噪聲的腦電信號(hào)；

（5）用均方根誤差（RMSE）和信噪比（SNR）來說明消噪效果，表達(dá)式如下：

其中，xi為i時(shí)刻在原始信號(hào)上的采樣值；xi是去除噪聲后的值；N是信號(hào)序列長(zhǎng)度。信噪比SNR反映了去除噪聲的信號(hào)中噪聲所占的比值，SNR和均方根誤差RMSE是呈反相關(guān)系的，研究需要信噪比大的，信號(hào)的去噪效果就更好。

4 實(shí)驗(yàn)結(jié)果與分析

實(shí)驗(yàn)數(shù)據(jù)來源于身體健康的在校研究生，在放松閉眼狀態(tài)下，1000赫茲為采樣的頻率，采樣的時(shí)間是3秒。實(shí)驗(yàn)中有噪聲的腦電信號(hào)與加入的白噪聲標(biāo)準(zhǔn)差之比為5，總體的平均次數(shù)是100。腦電信號(hào)在F3通道里獲取，然后進(jìn)行EEMD分解，最終得到9個(gè)IMF分量和一個(gè)剩余量R，所得IFM分量的時(shí)域圖如圖3、4所示。由圖可知，不同的時(shí)間特征尺度對(duì)應(yīng)各階IMF分量，并且隨著階次的增長(zhǎng)，頻率從高頻到低頻出現(xiàn)。

圖3EEMD分解后IMF分量時(shí)域圖

圖4EEMD分解后IMF分量時(shí)域圖

由于腦電信號(hào)主要集中在0-30HZ，而高頻部分主要分布在前3層的IMF分量中，所以噪聲主要存在于這三層，然后進(jìn)行前3個(gè)IMF分量的Hilbert譜分析可以得出這3層含有噪聲（如圖5），因?yàn)辄c(diǎn)分布比較分散，有少量的能量。因此，只對(duì)這三層IMF分量進(jìn)行改善，利用小波包去噪處理。消噪效果如圖6所示。

圖5 Hilbert譜

圖6 小波包去噪后的IMF1，IMF2，IMF3

最后，把通過小波包去噪后的IMF1～I(xiàn)MF3和其他保留的IMF分量進(jìn)行疊加，可以得到改進(jìn)EMD結(jié)合小波包去噪后的腦電圖，和小波包消噪、改進(jìn)EMD去噪相對(duì)比，三種方法的消噪結(jié)果如圖8、9、10所示。

圖7 原始信號(hào)

圖8 小波包去噪

圖9 改進(jìn)EMD去噪

圖10 改進(jìn)EMD結(jié)合小波包去噪

三種去噪方法的2個(gè)評(píng)價(jià)指標(biāo)如表1所示。

表1 3種消噪方法評(píng)價(jià)指標(biāo)比較

由表1可知，改進(jìn)EMD結(jié)合小波包方法比單獨(dú)使用小波包方法具有更高的信噪比和更低的均方根誤差，故消噪效果更佳。這是因?yàn)楸疚氖褂玫姆椒ú粌H可以抑制模態(tài)混疊現(xiàn)象，而且可以得到更多高頻部分的有用信息。

5 結(jié)語

本文把改進(jìn)EMD的HHT方法結(jié)合小波包來進(jìn)行腦電信號(hào)的去噪。運(yùn)用EEMD在EEG信號(hào)中加入白噪聲，接著對(duì)高頻部分對(duì)應(yīng)的IMF分量進(jìn)行HHT譜分析，可以得到含噪聲的高頻分量，最后對(duì)其進(jìn)行小波包去噪處理，可以更大限度的得到有用信息，有效的消除噪聲。經(jīng)驗(yàn)證，改進(jìn)EMD結(jié)合小波包去噪效果更好。有利于未來推動(dòng)臨床醫(yī)學(xué)方法的研究及一些疾病的診斷與預(yù)防。

參考文獻(xiàn)

［1］ 王兆源，周龍旗.腦電信號(hào)的分析方法［J］.第一軍醫(yī)大學(xué)學(xué)報(bào)，2000，20（2）：189-190.

［2］ 季忠，秦樹人，彭麗玲.腦電信號(hào)的現(xiàn)代分析方法［J］.重慶大學(xué)學(xué)報(bào)，2002，25（9）：108-112.

［3］ 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］.J Proc Rsoc Londna，1998，454（1971）：903-995.

［4］ 韓立喜.心音信號(hào)去噪及特征值提取的研究［D］.長(zhǎng)春：長(zhǎng)春理工大學(xué)，2013.

［5］ 黃友朋，趙山，許凡，方彥軍.EEMD排列熵與PCA-GK的滾動(dòng)軸承聚類故障診斷［J］.河南科技大學(xué)學(xué)報(bào)：自然科學(xué)版，2017，38（02）：17-24，30.

［6］ Yang BH，Yan GZ，Yan RG，et al.Adaptive subject-based feature extraction in brain-computer interfaces using wavelet packet best basis decomposition［J］.Medical Engineering&Physics，2007，29（1）：48-53.

［7］ Mcfarland DJ，Anderson CW，Muller KR，et al.BCI meeting 2005-workshop on BCI signal processing：feature extraction and translation［J］.IEEE Trans.Neural Syst.Rehabil.Eng.，2006，14（2）：135-138.

［8］ 張國華，袁中凡，李彬彬.心音信號(hào)特征提取小波包算法研究［J］.振動(dòng)與沖擊，2008，27（7）：47-49.

［9］ 楊幫華，陸文宇，何美燕，等.腦機(jī)接口中基于WPD和CSP 的特征提取［J］.儀器儀表學(xué)報(bào)，2012，33（11）：2560-2565.