基于VMD與CSP的腦電特征提取方法

劉 帥，烏日￣開西·艾依提

(新疆大學(xué)機(jī)械工程學(xué)院，新疆 烏魯木齊 830047)

1 引言

腦機(jī)接口(Brain Computer Interface， BCI)為人類與周圍環(huán)境進(jìn)行交互提供了一種不依賴于外周神經(jīng)和肌肉的新的輸出通道[1，2]。腦電圖(Electroencephalographic，EEG)具有無創(chuàng)、高時(shí)間分辨率、低成本等優(yōu)點(diǎn)，被廣泛應(yīng)用于腦機(jī)接口中。研究表明，當(dāng)人在想象左手或右手運(yùn)動(dòng)時(shí)，大腦的感覺運(yùn)動(dòng)皮層Mu節(jié)律和Beta節(jié)律能量會(huì)產(chǎn)生變化，這一現(xiàn)象被稱為事件相關(guān)去同步(Event-related Desynchronization， ERD)和事件相關(guān)同步(Event-related Synchronization， ERS)[3]。近年來，運(yùn)動(dòng)想象腦機(jī)接口技術(shù)被用來訓(xùn)練神經(jīng)肌肉損傷的患者，幫助其恢復(fù)肢體的運(yùn)動(dòng)感知功能[4，5]。由于運(yùn)動(dòng)想象腦電信號(hào)是由大腦活動(dòng)產(chǎn)生的，信噪比低，容易受到外界噪聲的影響，難以獲得較高的識(shí)別精度。如何對(duì)腦電信號(hào)進(jìn)行有效的特征提取和分類成為一個(gè)難點(diǎn)問題。

目前，國(guó)內(nèi)外研究者提出了若干運(yùn)動(dòng)想象腦電特征提取和分類算法，常用的特征提取算法有自回歸模型(Auto Regressive，AR)、小波包變換(Wavelet Package Transform，WPD)、功率譜密度(Power Spectral Density，PSD)、共空間模式(Common Spatial Pattern，CSP)等[6-8]。常用的分類算法有支持向量機(jī)(Support Vector Machine，SVM)、線性判別分析(Linear Discriminant Analysis，LDA)等[9]。楊默涵等[10]提出了一種基于總體經(jīng)驗(yàn)?zāi)B(tài)分解(Ensemble Empirical Mode Decomposition，EEMD)的運(yùn)動(dòng)想象腦電處理方法，將分解后的分量利用希爾伯特變換提取邊際譜及瞬時(shí)能譜，選取較優(yōu)頻帶波幅特征作為特征，取得了較好的分類效果。Shalu Chaudhary等[11]采用兩種時(shí)頻分析方法將腦電信號(hào)轉(zhuǎn)化為圖像，輸入到基于AlexNet模型的深度卷積神經(jīng)網(wǎng)絡(luò)(Deep Convolutional Neural Networks，DCNN)進(jìn)行特征提取與分類，發(fā)現(xiàn)連續(xù)小波變換(Continuous Wavelet Transform， CWT)轉(zhuǎn)化的圖像具有更高的識(shí)別率。Yu Z等[12]提出了一種時(shí)間約束稀疏組空間模式(Temporally Constrained Sparse Group Spatial Patterns，TSGSP)特征提取算法，通過同時(shí)優(yōu)化CSP中的濾波器頻帶和時(shí)間窗口提高運(yùn)動(dòng)想象腦電信號(hào)的分類精度。

為了進(jìn)一步提高在少數(shù)通道下的運(yùn)動(dòng)想象腦電信號(hào)的分類準(zhǔn)確率，本文提出一種變分模態(tài)分解(Variational Mode Decomposition，VMD)結(jié)合CSP的特征提取方法。將腦電信號(hào)進(jìn)行VMD分解，得到多個(gè)具有物理意義的固有模態(tài)函數(shù)(Intrinsic Mode Function，IMF)，通過巴氏距離計(jì)算各階固有模態(tài)函數(shù)與原始信號(hào)之間的相關(guān)程度，選取合適的固有模態(tài)函數(shù)組合成為新的信號(hào)矩陣作為輸入，再通過CSP構(gòu)造空間濾波器，得到各IMF的空間分布特征，最后使用支持向量機(jī)進(jìn)行分類。

2 數(shù)據(jù)集描述

本文采用兩個(gè)數(shù)據(jù)集，第一個(gè)是BCI Competition II的Data set III的數(shù)據(jù)集[13]，數(shù)據(jù)集記錄了一位25歲的健康女性，實(shí)驗(yàn)要求被試坐在屏幕前，根據(jù)屏幕上隨機(jī)出現(xiàn)的左右箭頭提示進(jìn)行相應(yīng)的左右手運(yùn)動(dòng)想象。數(shù)據(jù)采用雙導(dǎo)聯(lián)方式記錄，采集了C3，Cz和C4三個(gè)通道的EEG信號(hào)。采樣頻率為128Hz，對(duì)信號(hào)進(jìn)行0.5-30Hz濾波。實(shí)驗(yàn)共記錄7組，每組包含40次實(shí)驗(yàn)，每次實(shí)驗(yàn)持續(xù)9s，共280次實(shí)驗(yàn)，選其中的140次實(shí)驗(yàn)作為訓(xùn)練樣本，剩余的140次實(shí)驗(yàn)作為測(cè)試樣本。單次實(shí)驗(yàn)的采集范式如圖1所示。

圖1 Data set III單次實(shí)驗(yàn)采集范式

第二個(gè)數(shù)據(jù)集來自三位身體健康的志愿者，均為右利手，且均為第一次參加腦電實(shí)驗(yàn)。使用Neuracle公司的NeuSenW16通道無線腦電采集系統(tǒng)進(jìn)行實(shí)驗(yàn)數(shù)據(jù)的采集，腦電帽的電極按照國(guó)際標(biāo)準(zhǔn)10-20系統(tǒng)放置，僅采集了C3、C4兩個(gè)通道的數(shù)據(jù)。參考電極位于頭頂部，接地電極位于前額。采樣頻率設(shè)置為1000Hz，初始帶通濾波為0.5-100Hz，陷波濾波50Hz。

在實(shí)驗(yàn)過程中，受試者佩戴腦電帽被要求坐在距離電腦屏幕大約1米的椅子上，雙手自然放置在桌子上，保持放松狀態(tài)。在想象過程中，要求被試盡量避免眨眼、眼動(dòng)、身體以及頭部移動(dòng)。單次實(shí)驗(yàn)數(shù)據(jù)采集范式如圖2所示。試驗(yàn)開始前2s為準(zhǔn)備階段，屏幕顯示空白；2s后有一段短暫的提示音并且屏幕顯示“+”字形，提示被試實(shí)驗(yàn)馬上開始；4-7s屏幕會(huì)隨機(jī)出現(xiàn)向左或向右的箭頭，被試需要根據(jù)提示進(jìn)行左右手的運(yùn)動(dòng)想象；7s后有2-4s的隨機(jī)休息時(shí)間；為了防止被試疲勞，每隔20次實(shí)驗(yàn)休息2min。實(shí)驗(yàn)共分為3組，每組100次實(shí)驗(yàn)，共300次實(shí)驗(yàn)，其中200次實(shí)驗(yàn)作為訓(xùn)練樣本，100次實(shí)驗(yàn)作為測(cè)試樣本。

圖2 實(shí)驗(yàn)室單次實(shí)驗(yàn)數(shù)據(jù)采集范式

3 改進(jìn)的特征提取方法

因個(gè)體差異較大，每位受試者在運(yùn)動(dòng)想象時(shí)，ERD/ERS發(fā)生的頻率段不相同。傳統(tǒng)的共空間模式是根據(jù)固定頻率段的腦電信號(hào)進(jìn)行特征提取，且需要大量通道的數(shù)據(jù)。針對(duì)這一問題，本文通過結(jié)合VMD和CSP用于改善特征提取算法。改進(jìn)后的特征提取流程如圖3所示。首先，將濾波后的信號(hào)進(jìn)行VMD分解成K個(gè)IMF分量，然后通過計(jì)算各IMF分量與原始信號(hào)的巴氏距離，選擇巴氏距離最小的3階IMF分量進(jìn)行重構(gòu)，采用CSP濾波提取特征，得到特征向量F。

3.1 變分模態(tài)分解

變分模態(tài)分解通過迭代搜索變分模型的最優(yōu)解來確定每個(gè)模態(tài)分量的中心頻率和有限帶寬，可以有效抑制經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)方法中存在的端點(diǎn)效應(yīng)和模態(tài)混疊現(xiàn)象，具有堅(jiān)實(shí)的理論基礎(chǔ)，且對(duì)采樣和噪聲的魯棒性更強(qiáng)[14]。

圖3 VMD與CSP相結(jié)合的特征提取流程

在VMD算法中，定義固有模態(tài)函數(shù)(IMF)為一個(gè)調(diào)幅-調(diào)頻信號(hào)，通過混合一個(gè)預(yù)估中心頻率，將每個(gè)模態(tài)的頻譜調(diào)制到相應(yīng)的基頻帶

(1)

計(jì)算梯度平方L2范數(shù)對(duì)信號(hào)解調(diào)，得到各個(gè)模態(tài)函數(shù)的帶寬，則相約束變分問題如下：

(2)

式中，{uk}={u1，u2，…uk}為分解后的k個(gè)模態(tài)分量，{ωk}={ω1，ω2，…ωk}為每個(gè)模態(tài)對(duì)應(yīng)的中心頻率，f為輸入的腦電信號(hào)。

使用二次懲罰因子α和Lagrange乘法算子λ，將約束變分問題轉(zhuǎn)變成非約束變分問題，增廣型的拉格朗日表達(dá)式如下

L({uk}，{ωk}，λ)

(3)

利用交替方向乘子法(Alternate Direction Method of Multipliers，ADMM)，并結(jié)合Parseval/Plancherel 傅里葉等距變換，優(yōu)化得到各模態(tài)分量和中心頻率，求解增廣拉格朗日函數(shù)的最優(yōu)解，交替尋優(yōu)迭代un+1，ωn+1和λn+1的表達(dá)式如下：

(6)

迭代終止條件

(7)

最終得到K個(gè)有限帶寬的IMF分量。

3.2 模態(tài)函數(shù)選擇

利用巴氏距離用于測(cè)量IMF分量與原始信號(hào)的相似性。對(duì)于相同定義域X上的離散概率分布p和q，巴氏距離為

DB(p，q)=-ln(BC(p，q))

(8)

(9)

其中，p表示第K個(gè)IMF分量，q表示原始信號(hào)，BC(p，q)為離散概率分布 Bhattacharyya系數(shù)。

3.3 傳統(tǒng)共空間模式算法

共空間模式是一種對(duì)兩分類任務(wù)下的空域?yàn)V波特征提取算法，能夠從多通道的腦電信號(hào)中提取每一類的空間分布成分。共空間模式的原理是利用矩陣對(duì)角化，找到一組最優(yōu)空間濾波器對(duì)信號(hào)進(jìn)行投影，使得一類方差最大化，另一類方差最小化，從而得到具有較高區(qū)分度的特征向量[15]。

單次運(yùn)動(dòng)想象腦電實(shí)驗(yàn)數(shù)據(jù)由一個(gè)N×M的矩陣X組成，其中N表示通道數(shù)，M表示采樣點(diǎn)數(shù)，其歸一化后的協(xié)方差矩陣計(jì)算如下

(10)

(11)

對(duì)混合協(xié)方差矩陣Vc進(jìn)行特征分解

(12)

式中，Uc為特征向量矩陣，λc為特征值矩陣。構(gòu)造白化矩陣P

(13)

(14)

式中，Sl和Sr具有共同的特征矢量E，且對(duì)應(yīng)的特征值之和為1，則Sl最大特征值對(duì)應(yīng)的特征向量是Sr最小特征值對(duì)應(yīng)的特征向量。對(duì)Sl和Sr進(jìn)行特征分解可得

Sl=EλlETSr=EλrET

(15)

式中，對(duì)角矩陣λl和λr之和為單位矩陣I。由此構(gòu)造投影矩陣W

W=ETP

(16)

則對(duì)輸入矩陣X進(jìn)行濾波可得Z=WX，取Z的前m列和后m列構(gòu)成特征矩陣Zp，p=1，2，…，2m。

(17)

3.4 支持向量機(jī)

支持向量機(jī)是一種有效的二分類模型，具有較強(qiáng)的泛化能力。其基本原理是將輸入的向量映射到高維的特征空間，并在此空間尋找一個(gè)最優(yōu)分類超平面，使兩類樣本的分類間隔最大[16]。本文將改進(jìn)的CSP特征提取出來的特征向量作為輸入，使用LIBSVM工具包進(jìn)行分類。訓(xùn)練集用來訓(xùn)練模型，測(cè)試集用來評(píng)估最優(yōu)超參數(shù)訓(xùn)練出來的模型，同時(shí)對(duì)模型的參數(shù)進(jìn)行尋優(yōu)。采用網(wǎng)格搜索算法尋找最佳懲罰系數(shù)c和核參數(shù)g組合。

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

4.1 數(shù)據(jù)預(yù)處理

根據(jù)ERD/ERS現(xiàn)象，在做肢體運(yùn)動(dòng)想象時(shí)，感覺運(yùn)動(dòng)區(qū)域的節(jié)律變化主要發(fā)生在Mu節(jié)律(8-13Hz)和Beta節(jié)律(14-30Hz)，本文使用6階Butterworth濾波器對(duì)原始EEG信號(hào)進(jìn)行8-30Hz的濾波。

利用連續(xù)小波變換對(duì)信號(hào)進(jìn)行時(shí)頻分析，通過時(shí)頻圖觀察運(yùn)動(dòng)想象期間ERD發(fā)生的時(shí)間段。想象左右手運(yùn)動(dòng)時(shí)的時(shí)頻圖與腦地形圖如圖4和圖5所示。想象左手運(yùn)動(dòng)時(shí)，在8-12Hz頻段對(duì)側(cè)的C4通道3.5-7s時(shí)間段的能量顯著降低，發(fā)生ERD現(xiàn)象；而同側(cè)的C3通道3.5-7s時(shí)間段的能量顯著增加，發(fā)生ERS現(xiàn)象。與之相反，在想象右手運(yùn)動(dòng)時(shí)，對(duì)側(cè)的C3通道3.5-7s發(fā)生ERD現(xiàn)象，同側(cè)的C4通道3.5-7s發(fā)生ERS現(xiàn)象。對(duì)于BCI Competition II的Data set III的數(shù)據(jù)集，本文選擇單次實(shí)驗(yàn)的3.5-7s時(shí)間段的腦電信號(hào)進(jìn)行后續(xù)的特征提取。數(shù)據(jù)集中每次實(shí)驗(yàn)共采集9s數(shù)據(jù)，采樣頻率128Hz，每次實(shí)驗(yàn)每個(gè)通道共有1152個(gè)采樣點(diǎn)。選取其中3.5-7s的數(shù)據(jù)，此時(shí)數(shù)據(jù)為448×3×140。

圖4 想象左手運(yùn)動(dòng)時(shí)頻圖與腦地形圖

圖5 想象右手運(yùn)動(dòng)時(shí)頻圖與腦地形圖

4.2 特征提取與分類

將預(yù)處理后的腦電信號(hào)進(jìn)行VMD自適應(yīng)分解，得到K個(gè)IMF分量，每階IMF分量都有各自的中心頻率和有限帶寬。然后計(jì)算各IMF分量與原始信號(hào)的巴氏距離來選擇合適的IMF分量。本文K取值為7。以想象左手為例，C3通道的腦電信號(hào)經(jīng)VMD分解后的IMF波形及各分量對(duì)應(yīng)的頻譜圖如圖6所示。從頻譜圖中可以看出，自適應(yīng)分解后的IMF分量主要集中在8-30Hz范圍，不同于EMD分解時(shí)出現(xiàn)的相同成分的信息出現(xiàn)在不同IMF分量中，VMD分解后的各階分量均為有限帶寬，有效克服了EMD中的模態(tài)混迭現(xiàn)象。

圖6 VMD分解后的各階IMF分量及各分量對(duì)應(yīng)的頻譜圖

根據(jù)式(8)計(jì)算各IMF分量與原始信號(hào)的巴氏距離。以單次實(shí)驗(yàn)數(shù)據(jù)為例，表1為C3和C4通道經(jīng)VMD分解后各IMF分量與原始信號(hào)的巴氏距離。巴氏距離越小，對(duì)應(yīng)的模態(tài)分量與原始信號(hào)越相似。觀察表1可知，C3通道巴氏距離最小的三階IMF分量分別為IMF2、IMF3和IMF4；Cz通道巴氏距離最小的三階IMF分量分別為IMF2、IMF1和IMF3；C4通道巴氏距離最小的三階IMF分量分別為IMF2、IMF4和IMF7。

表1 IMF分量與原始信號(hào)的巴氏距離

本文選取每個(gè)通道數(shù)據(jù)分解后與原始信號(hào)巴氏距離最小的三階IMF分量重新構(gòu)造新的輸入矩陣作為特征提取的輸入，原輸入矩陣重新構(gòu)造成一個(gè)448×9的矩陣，再通過CSP濾波提取特征得到特征向量，最終將特征向量輸入到支持向量機(jī)中進(jìn)行分類。

4.3 實(shí)驗(yàn)對(duì)比

為了證明本文方法的有效性，在相同腦電數(shù)據(jù)集情況下，與其它文獻(xiàn)的方法進(jìn)行了比較，對(duì)比結(jié)果如表2所示?？梢钥闯?，在使用同一數(shù)據(jù)集的情況下，相比于其它三種算法，本文所提方法的分類結(jié)果取得了91.43%的分類準(zhǔn)確率。說明本文所提出的特征提取算法能夠有效提高運(yùn)動(dòng)想象腦電的分類準(zhǔn)確率。

表2 本文方法與其它方法分類準(zhǔn)確率對(duì)比

為了進(jìn)一步驗(yàn)證該方法的有效性，在實(shí)測(cè)數(shù)據(jù)中進(jìn)行了驗(yàn)證，與傳統(tǒng)的特征提取算法CSP進(jìn)行了比較。如表3所示，對(duì)實(shí)驗(yàn)室實(shí)測(cè)的三名被試的C3和C4通道的腦電信號(hào)分別進(jìn)行CSP特征提取和VMD融合CSP的特征提取，并用支持向量機(jī)進(jìn)行分類。從表3可以看出，相比于傳統(tǒng)的CSP特征提取算法，經(jīng)VMD分解后的信號(hào)再進(jìn)行CSP濾波提取特征能夠獲得更高的分類識(shí)別率，三名被試的分類準(zhǔn)確率分別提高了6%、8%、6%。

表3 實(shí)測(cè)數(shù)據(jù)的分類識(shí)別率比較

5 結(jié)論

針對(duì)運(yùn)動(dòng)想象腦電信號(hào)個(gè)體差異較大，傳統(tǒng)的CSP特征提取算法需要大量通道的數(shù)據(jù)的問題，本文提出了一種VMD結(jié)合CSP的特征提取算法。通過巴氏距離可以選擇與腦電信號(hào)最相關(guān)的固有模態(tài)函數(shù)構(gòu)造信號(hào)矩陣進(jìn)行特征提取，最后使用SVM分類。在公共數(shù)據(jù)集上取得了91.43%的分類準(zhǔn)確率，在實(shí)測(cè)數(shù)據(jù)上相比于傳統(tǒng)的共空間模式特征提取方法，取得了更好的分類效果。

實(shí)驗(yàn)結(jié)果表明該方法在少數(shù)通道的情況下能夠有效提取運(yùn)動(dòng)想象腦電信號(hào)的特征，取得了較好的分類準(zhǔn)確率，證明了該方法的有效性，為運(yùn)動(dòng)想象腦機(jī)接口應(yīng)用奠定了基礎(chǔ)。在下一階段，本文提出的方法將用于實(shí)現(xiàn)在線腦機(jī)接口。