徐雨，鄭威，程怡

(江蘇科技大學(xué)，鎮(zhèn)江 212000)

引言

腦電信號(hào)[1]是一種極微弱的生理信號(hào)，采集過(guò)程中易受肌電偽跡和眼電偽跡等生理偽跡影響。由于眼電偽跡幅值較大且與腦電信號(hào)頻率相近，對(duì)腦電信號(hào)產(chǎn)生很大干擾[2]，致使腦電信號(hào)產(chǎn)生畸變，影響分析結(jié)果。

為得到無(wú)干擾的腦電信號(hào)，國(guó)內(nèi)外團(tuán)隊(duì)研究了多種去除眼電偽跡的方法[3-4]。其中，有研究認(rèn)為在實(shí)驗(yàn)對(duì)象采集信號(hào)時(shí)不轉(zhuǎn)動(dòng)眼球，可避免眼電偽跡的產(chǎn)生，但實(shí)際情況中眼球活動(dòng)不可控制，所以在腦電信號(hào)采集過(guò)程中仍會(huì)有大量眼電偽跡產(chǎn)生[5]。目前，獨(dú)立成分分析(independent component analysis, ICA[6])、經(jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition, EMD[7])等算法可有效去除眼電偽跡。其中，ICA算法只能處理單一的非線性函數(shù)[8]，去偽跡時(shí)需要人工識(shí)別偽跡成分。EMD可分解非線性和非平穩(wěn)信號(hào)，但其在分解過(guò)程中存在模態(tài)混疊問(wèn)題[9]。

針對(duì)EMD和ICA在去除腦電信號(hào)中眼電偽跡時(shí)存在的弊端，本研究提出了自適應(yīng)噪聲完整經(jīng)驗(yàn)?zāi)B(tài)分解(complete ensemble empirical mode decomposition with adaptive noise, CEEMDAN)、(kernel independent compoment analysis，KICA)[10]和樣本熵相結(jié)合的方法用以去除腦電信號(hào)中眼電偽跡。該方法可以保留更多有用信息，更具有魯棒性。

1 基本原理

本研究提出的眼電偽跡去除方法流程見(jiàn)圖1。首先，對(duì)腦電信號(hào)進(jìn)行帶通濾波，對(duì)其進(jìn)行自適應(yīng)白噪聲完整經(jīng)驗(yàn)?zāi)B(tài)分解得到模態(tài)分量[11]；再對(duì)模態(tài)分量進(jìn)行KICA降維；隨后對(duì)得到的獨(dú)立分量進(jìn)行熵值計(jì)算，利用閾值判別并去除眼電偽跡，最后重構(gòu)信號(hào)。

圖1 去除腦電信號(hào)中眼電偽跡的流程圖

1.1 CEEMDAN算法

EMD分析時(shí)頻分布的能力較強(qiáng)，但在分解過(guò)程中會(huì)存在模態(tài)混疊問(wèn)題。為解決此問(wèn)題，王海龍等[12]提出在EMD中加入白噪聲的EEMD算法。但EEMD[13]不具有完備性，通過(guò)對(duì)EEMD分解后的經(jīng)驗(yàn)?zāi)B(tài)分量自適應(yīng)地加入白噪聲，解決了EEMD分解中無(wú)完備性的問(wèn)題。

CEEMDAN算法的步驟具體如下：

(1)構(gòu)造加入噪聲的序列：xi(t)=x(t)+ω0εi(t),i∈{1,2,…,I}，其中，x(t)為觀測(cè)信號(hào)，εi(t)為獨(dú)立的方差為1的高斯白噪聲，ω0為噪聲的標(biāo)準(zhǔn)差。最后得到的模態(tài)分量記為cq(I).

(4)設(shè)分解Q層，得出每個(gè)階段第q個(gè)剩余分量和第q+1個(gè)模態(tài)分量(q=1,…,Q-1)，再次進(jìn)行步驟(2)的計(jì)算：

1.2 核獨(dú)立成分分析

ICA算法是基于非線性函數(shù)期望值的一個(gè)固定形式[14]，由于對(duì)比函數(shù)的局限性，導(dǎo)致處理更復(fù)雜的非線性形式的盲源問(wèn)題時(shí)，精確度低[15-16]。KICA在ICA基礎(chǔ)上引入核函數(shù)[17]，其使用的函數(shù)空間適用于各種不同分布的源。與ICA相比，KICA具有更好的分離能力和更好的魯棒性[18-19]。

KICA算法的步驟如下[20]：

(1)輸入樣本集A=(a1,a2,…aI)，xi∈Rn，i=1,2,…，I，利用核函數(shù)代替向量的內(nèi)積，并對(duì)數(shù)據(jù)進(jìn)行白化處理。

(2)中心化數(shù)據(jù)。

(3)將RKHS中的非線性函數(shù)用作比較函數(shù)。

(4)最小化比較函數(shù),輸出分離矩陣W。

1.3 基于熵的偽跡判斷

樣本熵算法是一種判別時(shí)間序列復(fù)雜程度的算法[21-22]。與腦電信號(hào)相比，眼電信號(hào)比較規(guī)律，對(duì)應(yīng)熵值較低，通過(guò)熵值的判別可以去除熵值較低的眼電偽跡成分。與近似熵相比，樣本熵[23]的計(jì)算對(duì)數(shù)據(jù)長(zhǎng)度要求低，且具有更好的穩(wěn)定性[24]。樣本熵的計(jì)算如下：

(1)設(shè)長(zhǎng)度為P的觀測(cè)信號(hào)為{a(i),1≤i≤P}，即可得到一組e維的矢量:Ae(1),…,Ae(P-m+1)，其中Ae(i)={a(i),a(i+1),…,a(i+e-1)},1≤i≤P-e+1。其中，e為相空間維數(shù)。

(2)定義d[Ae(i),Ae(j)]為向量Ae(i)與Ae(j)之間的距離最大差值的絕對(duì)值。

(4)增加維數(shù)到e+1，循環(huán)步驟(1)到(3)，同理得到Ve+1(r)。樣本熵為：

(5)若N為有限值，樣本熵可為：

式中：o為相似容限。本實(shí)驗(yàn)運(yùn)用Gomez-Herrero提出的閾值判別式[25]識(shí)別偽跡。首先，定義變量c滿(mǎn)足1

φ(c+1)-φ(c)<φ(c)-φ(c-1)

(1)

式中，φ(c)表示KICA得到的獨(dú)立分量升序排序后的第c個(gè)分量的熵值，定義[P/2]為不大于P/2的最大整數(shù)。

若存在變量c滿(mǎn)足上述條件，c取滿(mǎn)足條件的最小整數(shù)，則φ(c)為獨(dú)立分量閾值，由于眼電偽跡的熵值比腦電信號(hào)的熵值更低，因此，判定排序前c個(gè)熵值對(duì)應(yīng)的分量為眼電偽跡，對(duì)其進(jìn)行置零。

2 實(shí)驗(yàn)

2.1 實(shí)驗(yàn)數(shù)據(jù)來(lái)源

本實(shí)驗(yàn)?zāi)X電信號(hào)采用BCI Competition IV的Data sets 2b 的數(shù)據(jù)集。實(shí)驗(yàn)對(duì)象為四名健康人員，采樣頻率250 Hz。實(shí)驗(yàn)對(duì)象都是右利手，有正常或者矯正正常視力。實(shí)驗(yàn)采集分為三個(gè)部分：(1)2 min睜眼和1 min閉眼(2)眼電信號(hào)采集(3)腦電信號(hào)采集。

由圖2可知，對(duì)象A的腦電信號(hào)中800和1 500采點(diǎn)處存在眼電偽跡信號(hào)；對(duì)象B的腦電信號(hào)中1 200采點(diǎn)處存在眼電偽跡信號(hào)；對(duì)象C的腦電信號(hào)中220采點(diǎn)處存在眼電偽跡；對(duì)象D的腦電信號(hào)中220、1 750采點(diǎn)處存在眼電偽跡。

圖2 含有偽跡的原始腦電信號(hào)

2.2 實(shí)驗(yàn)步驟

本實(shí)驗(yàn)具體步驟：

(1)對(duì)濾波后的腦電信號(hào)進(jìn)行自適應(yīng)白噪聲完整經(jīng)驗(yàn)?zāi)B(tài)分解。

(2)通過(guò)核獨(dú)立成分分析，對(duì)模態(tài)分解得到的模態(tài)分量進(jìn)行分析，得到核獨(dú)立分量。

(3)對(duì)獨(dú)立分量進(jìn)行樣本熵計(jì)算，利用閾值判別式判別眼電分量，若為偽跡，則將其置零去除。

(4)對(duì)去除偽跡后的信號(hào)進(jìn)行核獨(dú)立成分分析逆變換，隨后對(duì)逆變換后的信號(hào)進(jìn)行疊加，最終得到去除偽跡后的腦電信號(hào)。

圖3為自適應(yīng)白噪聲完整經(jīng)驗(yàn)?zāi)B(tài)分解后得到的分量圖。

圖3 CEEMDAN分解得到的模態(tài)分量Fig.3 The modal components obtained by CEEMDAN decomposition

核獨(dú)立成分分析獲得的分量見(jiàn)圖4。對(duì)圖4所得的各個(gè)分量進(jìn)行熵值計(jì)算，所得熵值見(jiàn)表1、表2、表3、表4。

表1 對(duì)象A的各分量的熵值

表2 對(duì)象B的各分量的熵值

表3 對(duì)象C的各分量的熵值

圖4 核獨(dú)立分量Fig.4 Kernel-independent components

表4 對(duì)象D的各分量的熵值

通過(guò)閾值判別式(1)，可以判定表1中，sig2、sig10、sig11、sig12屬于眼電偽跡；表2中，sig2、sig9、sig10、sig11屬于眼電偽跡；表3中，sig2、sig9、sig10、sig11；表4中，sig2、sig10、sig11、sig12屬于眼電偽跡，將眼電偽跡成分置零。腦電信號(hào)中眼電偽跡去除前和去除后的效果對(duì)比圖，見(jiàn)圖5—圖8。

圖5 對(duì)象A的眼電偽跡去除前后對(duì)比圖

圖6 對(duì)象B的眼電偽跡去除前后對(duì)比圖

圖7 對(duì)象C的眼電偽跡去除前后對(duì)比圖

圖8 對(duì)象D的眼電偽跡去除前后對(duì)比圖

2.3 實(shí)驗(yàn)結(jié)果及分析

本研究通過(guò)均方根誤差(RMSE)和相關(guān)系數(shù)(R)來(lái)判別去除偽跡的效果，其公式如下：

式中：x和y分別表示初始腦電信號(hào)和偽跡去除后的腦電信號(hào)。實(shí)驗(yàn)對(duì)象A、B、C、D的RMSE與相關(guān)系數(shù)見(jiàn)表5。本研究選取對(duì)象A的數(shù)據(jù)來(lái)作為參考分析數(shù)據(jù)，對(duì)其進(jìn)行15次實(shí)驗(yàn)，取實(shí)驗(yàn)結(jié)果的平均值。為證實(shí)本方法的有效性，與KICA、CEEMDAN-FastICA兩種方法進(jìn)行對(duì)比，結(jié)果見(jiàn)表6。

表5 各實(shí)驗(yàn)對(duì)象的RMSE和相關(guān)系數(shù)

表6 各算法RMSE和相關(guān)系數(shù)

由表6可知，CEEMDAN-FastICA方法的相關(guān)系數(shù)最小，說(shuō)明FastICA在分析信號(hào)時(shí)可能受外界影響，分離能力下降，使眼電偽跡未得到有效去除，最終重構(gòu)得到的腦電信號(hào)與原腦電信號(hào)相差較大。而CEEMDAN-KICA方法RMSE最小，相關(guān)系數(shù)最大，說(shuō)明KICA較好地解決了非線性混合信號(hào)的問(wèn)題，并且重構(gòu)后的信息誤差較小，且能夠較好地保留有用的數(shù)據(jù)信息。

3 結(jié)論

本研究提出CEEMDAN和KICA相結(jié)合的眼電偽跡去除方法，首先，對(duì)腦電信號(hào)進(jìn)行模態(tài)分解，再對(duì)分解后得到的模態(tài)分量進(jìn)行核獨(dú)立成分分析，得到獨(dú)立分量，隨后對(duì)獨(dú)立分量進(jìn)行樣本熵計(jì)算，利用閾值對(duì)熵值進(jìn)行判別，將眼電分量置零，最后重構(gòu)去除偽跡后的信號(hào)。通過(guò)實(shí)驗(yàn)表明，本研究提出的眼電偽跡去除方法更準(zhǔn)確，能有效去除眼電偽跡，更具有魯棒性。