基于小波變換共空間模式的腦電信號解碼

曲思霖，王從慶，李建亮，展文豪，張 民

(1.南京航空航天大學(xué) 自動化學(xué)院，江蘇 南京 210016；2.中國航天員科研訓(xùn)練中心 人因工程國防科技重點實驗室, 北京 100094)

0 引言

在空間站中，由于失重，航天員不能像在地面一樣自如地控制機器完成指定操作，而腦機接口技術(shù)可以解碼腦電，利用腦電信號(EEG signal)控制空間機械臂運動。腦機接口(brain computer interface,BCI)[1]是一種不依賴于外圍神經(jīng)和肌肉組織的大腦正常輸出通路來實現(xiàn)大腦與外部環(huán)境之間的直接通信的技術(shù)。通過腦機接口可以將操作者的意圖傳遞給機器，達到人腦控制機器運動的效果。

腦電信號的特征提取方法主要有時域分析法、頻域分析法、時頻分析法、空域分析法。張紹榮等[2]提出基于稀疏貝葉斯邏輯回歸運動想象腦電信號分類模型，提取方差、峰度、峭度、6階自回歸模型(auto regression,AR)系數(shù)、通帶頻率以及小波能量作為特征。Mahmood等[3]采用共空間模式方法進行特征提取，用支持向量機(support vector machine,SVM)分類。Park等[4]采用特征向量中心性特征選擇方法，通過小波包分解、共空間模式提取特征，并采用內(nèi)核極限學(xué)習(xí)機分類。Sun等[5]采用共空間模式提取特征，將概率神經(jīng)網(wǎng)絡(luò)(probabilistic neural network,PNN)分成6類腦電信號。

由于運動想象信號主要表現(xiàn)為特定頻率的波段信號幅值變化，只采用時域特征無法體現(xiàn)，因此，常用時頻分析法和空域分析法進行分析?；诳焖俑道锶~變換的功率譜密度(power spectral density,PSD)表示信號功率在各頻率點的分布情況，該方法簡單易行，但事件和頻率不能同時滿足最佳效果[6]。小波包分解(wavelet packet decomposition,WPD)是進行信號時頻分析和處理的常用工具，通過伸縮平移對信號進行多尺度細化，達到高頻處時間細分、低頻處頻率細分的效果，適用于非平穩(wěn)信號處理[7]。共空間模式(common space pattern，CSP)是一種有效提取運動想象腦電信號特征的方法，于1999年首次用于運動想象信號分類[8]，目前在腦機接口系統(tǒng)中已得到廣泛應(yīng)用[9]。通過共空間模式計算空間濾波器，使不同類別的腦電信號特征差異最大化。

以上所提及的相關(guān)文獻中分別提取了運動想象信號的時域特征、頻域特征、空域特征，但沒有將3類特征結(jié)合起來，且未考慮由于受試者個人條件不同選擇不同特征的情況。本文針對小波包分解腦電信號，采用遞歸特性消除方法(recursive feature elimination,RFE)選擇特征，將經(jīng)過共空間模式濾波得到的特征向量輸入長短期記憶網(wǎng)絡(luò)進行分類，識別4類腦電信號，平均識別準(zhǔn)確率達93.69%。將識別的腦電信號作為控制信號，控制虛擬環(huán)境中的空間機械臂的中臂和大臂順時針或逆時針運動。

1 腦電信號特征提取與分類

1.1 信號預(yù)處理

采用基于負熵最大的FastICA方法對腦電信號進行盲源分離。由于相鄰?fù)ǖ啦杉男盘栔g相互干擾，且這些信號相互獨立，因此，常用獨立成分分析(independent component analysis,ICA)的方法對腦電信號進行分離[10]。

基于神經(jīng)生理經(jīng)驗知識，選擇大腦皮層區(qū)域感知區(qū)22通道采集的腦電信號。運動想象信號采集通道如圖1所示。

圖1 運動想象信號采集通道分布

以C3、C4、Cz通道信號為例，未經(jīng)獨立成分分析方法分離的運動想象信號如圖2所示。經(jīng)過FastICA盲源分離后得到C3、C4、Cz獨立成分如圖3所示。

圖2 C3、C4、Cz通道腦電信號

圖3 經(jīng)過FastICA盲源分離后C3、C4、Cz通道腦電信號

1.2 小波包分解

小波包中提供多種類型的小波基。對跳變較多的信號，Haar-Wsalsh基比較適用；局部三角函數(shù)基適用于含正弦分量或高頻振蕩的信號；對于運動想象腦電信號這類由分段多項式結(jié)構(gòu)組成的離散型的信號，常采用Daubechies小波基[11]。本文采用Db3小波(濾波器長度為3)對4類運動想象信號進行3層分解。

小波包分解原理如下[12]：

(1)

式中：s為尺度因子；b為平移因子；f(t)為源信號；ψ(·)為小波基函數(shù)。離散表達式為

(2)

式中：cl為逼近系數(shù)；gj,l為細節(jié)系數(shù)；φ(·)為尺度函數(shù)。

(3)

(4)

(5)

(6)

式中：h0為低通濾波器；h1為高通濾波器。

小波包分解如圖4所示。本文中數(shù)據(jù)集采集運動想象信號經(jīng)過Butterworth濾波后的頻率為0.5～100 Hz,因此第3層小波包節(jié)點對應(yīng)的頻率如表1所示。

圖4 小波包3層分解圖

表1 第3層小波包節(jié)點對應(yīng)的頻率

1.3 特征選擇

當(dāng)進行運動想象時，不同的受試者子頻帶變化不同。針對不同的受試者選擇不同的通道和頻帶信號可以去除每個受試者不重要的子頻帶，使剩余子頻帶有更好的可分性，可提高識別準(zhǔn)確率[13]。特征選擇具體步驟如下。

步驟1計算第3層中每個節(jié)點的功率值：

(7)

式中：l表示第l個通道；j表示第j層；n表示第n個節(jié)點；r表示采樣點個數(shù)。

步驟4將被選入V′中子頻帶的系數(shù)重新構(gòu)成一個新的矩陣E165×q，q為每個節(jié)點系數(shù)個數(shù)。

1.4 共空間模式提取特征

二分類共空間模式的原理為求取空間濾波器使得一類信號的投影矩陣方差最大，另一類信號的投影矩陣方差最小[14]，算法如下。

第k類腦電信號Xk可表示為

(8)

式中：N為腦電儀的通道數(shù)；T為總采樣點數(shù)。

求取腦電信號的規(guī)范化協(xié)方差矩陣Rk：

(9)

復(fù)合的規(guī)范化協(xié)方差矩陣RC：

RC=R1+R2。

(10)

對規(guī)范化協(xié)方差矩陣RC進行特征值分解：

(11)

式中：λC為特征值矩陣；UC為特征向量法矩陣。

構(gòu)造白化變換矩陣Q：

(12)

利用白化變換矩陣Q對協(xié)方差矩陣R1和R2進行處理：

S1=QR1QTS2=QR2QT。

(13)

式中：S1、S2具有相同的特征向量。構(gòu)造對角矩陣λ1、λ2和特征向量矩陣B，對S1、S2進行主分量分解：

S1=Bλ1BTS2=Bλ2BT，λ1+λ2=I。

(14)

若λ1中特征值按降序排列，則λ2中特征值按升序排列，即當(dāng)S1特征值最大時，S2特征值最小，保證了2類信號差異最大化。

求取投影矩陣，即空間濾波器M：

M=BTQ。

(15)

通過空間濾波器M得到第k類腦電信號特征矩陣Zm×T：

Zm×T=MT×XN×T。

(16)

式中：N為腦電儀的通道數(shù)；T為總采樣點數(shù)；m為生成空間濾波器時選取特征個數(shù)。本文中m=30,計算平均頻帶功率：

(17)

將4類運動想象信號兩兩組合，作為CSP濾波器的輸入，得到6個空間濾波器。將2類信號組合，通過空間濾波器的方法增大信號與另3類信號的差異性，使4類信號具有可分性，為后續(xù)分類提供良好的基礎(chǔ)。CSP提取4類腦電特征如圖5所示。

圖5 CSP提取腦電特征示意圖

1.5 腦電信號分類

本文選擇長短期記憶(long short term memory,LSTM)網(wǎng)絡(luò)進行分類。長短期記憶網(wǎng)絡(luò)是一種改進的循環(huán)神經(jīng)網(wǎng)絡(luò)算法,它使用一種被稱為 LSTM的記憶單元來判別哪些信息應(yīng)該被保留,控制信息從前一時刻到下一時刻進行傳輸,是目前應(yīng)用最為廣泛的具有記憶功能的網(wǎng)絡(luò)[15]。

每個空間濾波器選取30維特征，信號與共空間模式濾波器相乘，計算平均頻帶功率，6個空間濾波器可組成30×6的特征矩陣。選擇長短期記憶網(wǎng)絡(luò)搭建分類模型。LSTM的記憶單元共100個神經(jīng)元，輸出層選擇Softmax作為激活函數(shù)，得到4類分類結(jié)果，共迭代30次。訓(xùn)練模型損失函數(shù)選擇交叉熵函數(shù)，優(yōu)化器選擇Adam優(yōu)化器。

2 實驗部分

本文介紹了用于驗證的公開數(shù)據(jù)集，并將本文方法在數(shù)據(jù)集上測試得到的識別準(zhǔn)確率與其他特征提取和分類方法比較。

本文選用的數(shù)據(jù)集為奧地利大學(xué)提供的第4屆腦-機接口大賽的data set A。該數(shù)據(jù)集中包含9位受試者4類運動想象的數(shù)據(jù)。每個受試者進行288次單次實驗，其中每類運動想象各72次。每個實驗持續(xù)時間約為8 s。在第2 s時，屏幕出現(xiàn)上、下、左、右箭頭，受試者分別進行舌頭運動想象、雙腳運動想象、左手運動想象、右手運動想象。3～6 s時，受試者保持該動作想象。信號采樣頻率為250 Hz，并且通過0.5～100 Hz的帶通濾波器和50 Hz的陷波濾波。實驗共采集25個通道信號，包括22個腦電通道和3個眼電(electro-oculogram,EOG)通道[16]。

選擇前3位受試者的3～6 s的22個腦電采集通道的訓(xùn)練數(shù)據(jù)隨機劃分，其中75%作為訓(xùn)練集，25%作為測試集。與其他方法對比結(jié)果如表2所示。

表2 不同的特征提取方法和分類方法的識別準(zhǔn)確率對比

表2中LDA為線性判別分析。由表2可知，在采用WPD和CSP進行特征提取并用LSTM進行分類時，采用RFE選擇特征的識別準(zhǔn)確率略高于未選擇特征的識別準(zhǔn)確率。在采用WPD和CSP進行特征提取并用RFE進行特征選擇時，LSTM分類效果好于LDA。若未采用CSP提取特征，則幾乎無法正確分類。由于采用RFE選擇特征時，需要遍歷所有特征，且對準(zhǔn)確率影響較小，因此當(dāng)特征數(shù)量較大時，不建議采用。本文將不同類別信號兩兩組合進行CSP空域濾波器計算，即需要計算k(k-1)/2個空域濾波器(其中k為類別數(shù))，若類別過多，則需要計算的空域濾波器會大量增加。

3 空間機械臂的腦電解碼控制

將上述分類結(jié)果作為控制信號，通過串口在Unity 3D中控制空間機械臂運動，仿真實驗流程如圖6所示。

圖6 仿真實驗流程

空間機械臂是由底座(基座)、大臂、中臂、小臂、末端手爪組成，將其安裝在空間站上，主要完成空間在軌抓取目標(biāo)物體任務(wù)。航天員利用腦電解碼指令控制空間站外的空間機械臂。當(dāng)腦電信號為左手運動想象時，控制空間機械臂大臂逆時針轉(zhuǎn)動。當(dāng)腦電信號為右手運動想象時，控制械臂大臂順時針轉(zhuǎn)動。當(dāng)腦電信號為舌頭運動想象時，控制機械臂中臂順時針轉(zhuǎn)動。當(dāng)腦電信號為雙腳運動想象時，控制機械臂中臂逆時針轉(zhuǎn)動。各類腦電信號對應(yīng)控制指令及作用如表3所示，空間機械臂初始位置如圖7所示，空間機械臂按指令轉(zhuǎn)動效果如圖8所示。

表3 4類腦電信號對應(yīng)控制指令及機械臂動作

圖7 空間機械臂初始位置

圖8 空間機械臂按指令轉(zhuǎn)動后位置

4 結(jié)論

本文提出一種基于小波包分解的共空間模式方法提取腦電信號特征，采用遞歸特性消除方法選擇特征，通過長短期記憶網(wǎng)絡(luò)進行分類。首先，用快速獨立成分分析(FastICA)方法對22通道的腦電信號進行盲源分離。然后，用Db3小波分解腦電信號，計算每個節(jié)點的功率值，采用RFE選擇特征，將選擇的子頻帶和通道組成特征矩陣,通過共空間模式對腦電信號進行空間濾波，并采用長短期記憶神經(jīng)網(wǎng)絡(luò)對提取的腦電特征進行分類。將本文方法與其他特征提取方法在data set A公開數(shù)據(jù)集上進行測試。測試結(jié)果表明，本文方法識別率均優(yōu)于常用的幾種方法。最后，采用解碼的腦電信號作為控制信號，控制空間機械臂中臂和大臂順時針、逆時針轉(zhuǎn)動，驗證了本文方法可以用于解碼控制。

該方法可用于幫助航天員控制空間站外的空間機械臂在軌操作。腦電信號易受周圍環(huán)境干擾，而且可實現(xiàn)任務(wù)較少，人與人之間具有較大差異性，如何將訓(xùn)練好的模型通過遷移學(xué)習(xí)的方法解決腦機接口系統(tǒng)的局限性[17]是下一步研究的重點。