周杰，楊國雨，徐濤

杭州電子科技大學(xué)智能控制與機器人研究所，浙江杭州310018

前言

對于因嚴(yán)重運動障礙而喪失正常交流能力的患者來說，急需重新建立起一種新的與外界進(jìn)行交流的方法，腦機接口（Brain Computer Interface,BCI）正是這樣的一種方法。它通過信號處理、模式識別以及其他過程將與患者意圖相關(guān)的腦電信號（Electroencephalograph,EEG）轉(zhuǎn)換為控制信號，從而控制機器人手臂或輪椅等機器［1］。

由于EEG具有較高的時間分辨率，相對較低的成本，并且更方便患者使用的特點，BCI系統(tǒng)廣泛采用EEG技術(shù)檢測大腦活動［2］。其中基于運動想象腦電的BCI對中風(fēng)患者有重要意義，它可以幫助恢復(fù)受損神經(jīng)，在康復(fù)訓(xùn)練中發(fā)揮重要的作用［3］。在運動想象任務(wù)中，要求受試者通過想象諸如手或腳移動的特定動作來生成EEG信號［4］。然后，采用機器學(xué)習(xí)算法對EEG信號進(jìn)行分類以將其轉(zhuǎn)換為命令。由于數(shù)十億個神經(jīng)元之間的復(fù)雜互連，記錄的EEG信號本質(zhì)上是非線性、非平穩(wěn)的。因此，如何從復(fù)雜的運動想象腦電中提取可以有效識別運動任務(wù)的特征是非常重要的。目前，共同空間模式（Common Spatial Pattern,CSP）方法是提取運動想象腦電相關(guān)特征最常用和最有效的方法之一，并已廣泛用于運動想象BCI系統(tǒng)［5-6］。

CSP是一種針對兩類運動想象EEG的處理算法，其通過構(gòu)建空間濾波器使得兩類方差差異最大化，在兩類運動想象任務(wù)問題上得到很高的分類準(zhǔn)確率。對于多類運動想象分類任務(wù)，已經(jīng)提出一對多 CSP（One Versus Rest CSP,OVR-CSP）［7］和一對一CSP（One Versus One CSP,OVO-CSP）［8］，然而分類結(jié)果并不理想。這在一定程度上是因為當(dāng)運動想象EEG未選取運動想象相關(guān)頻率范圍時，對CSP特征執(zhí)行的分類通常會導(dǎo)致準(zhǔn)確性較差［9］。因此出現(xiàn)了一些改進(jìn)的CSP，例如子帶共同空間模式［9］、濾波器帶通用空間模式［10］。這些方法的思想是將CSP與頻率信息結(jié)合，使原始EEG分解成多個頻帶，并在每個頻帶上使用CSP來選擇最具有區(qū)別性的特征。

但是結(jié)合頻域信息的CSP仍存在不足，EEG本質(zhì)上是一個高度非平穩(wěn)的時間序列信號，時間序列信息對于提高CSP特征提取和后續(xù)分類的性能也至關(guān)重要［11］。遞歸神經(jīng)網(wǎng)絡(luò)（Recurrent Neural Network,RNN）是一類人工神經(jīng)網(wǎng)絡(luò)，它通過單元之間的連接形成一個有向循環(huán)創(chuàng)建了網(wǎng)絡(luò)的內(nèi)部狀態(tài)，從而可以有效利用輸入信號的時間信息進(jìn)行時間序列分析。近年來，RNN被廣泛應(yīng)用于時間序列的信號處理中，并已成為語音識別［12］、機器人運動識別［13］、圖像文本識別［14］、自然語言處理［15］等領(lǐng)域的有效模型。然而，傳統(tǒng)的RNN在訓(xùn)練過程中存在梯度爆炸或梯度消失的問題，導(dǎo)致RNN無法有效利用早期的歷史信息。針對這個問題，Hochreiter等［16］將記憶單元引入到傳統(tǒng)的RNN中以存儲長時間信息，提出長短期記憶網(wǎng)絡(luò)（Long Short-Term Memory,LSTM）。

為了有效利用EEG的空間頻率與時間序列信息，獲得更好的分類效果，本文結(jié)合離散小波變換（Discrete Wavelet Transform,DWT）、CSP和LSTM神經(jīng)網(wǎng)絡(luò)的特點，提出一種基于空間頻率與時間序列信息的多類運動想象EEG特征提取方法，并采用BCI競賽III中Dataset III-a的K3受試者數(shù)據(jù)進(jìn)行分類驗證。與傳統(tǒng)腦電特征提取算法相比，分類準(zhǔn)確率有較大提高。

1 方法

1.1 DWT

小波變換是一種時頻分析方法，由于其多分辨率特性，適合對EEG等非平穩(wěn)信號進(jìn)行特征提取。多分辨率分析的步驟如圖1所示。初始DWT由連續(xù)應(yīng)用低通濾波器和高通濾波器組成，離散信號經(jīng)DWT分解為多個小波分量。其中，h［n］表示高通濾波器，g［n］表示低通濾波器，x［n］表示同時通過高低通濾波器的EEG信號注釋，↓2表示兩倍下采樣。重復(fù)該過程，直到每個通道的EEG信號被分解成n層小波細(xì)節(jié)，即D1,D2,…,Dn，以及近似系數(shù)An。

圖1 離散小波多分辨率分析結(jié)構(gòu)圖Fig.1 Block diagram of multi-resolution analysis of discrete wavelet transform(DWT)

1.2 OVR-CSP

CSP方法的原理：首先，兩類實對稱矩陣（協(xié)方差矩陣）同時對角化，并應(yīng)用主成分分析方法排除兩種任務(wù)的共同部分，提取不同部分。進(jìn)而，再通過空間因子和相應(yīng)的空間過濾器將特定信號成分提取出來。為了將CSP應(yīng)用到多分類問題常采用一對多方法，該方法是把其中一類模式當(dāng)成一類，而余下的所有模式當(dāng)成另外的一類，從而形成一個兩類的CSP，這樣依次對每一類信號計算對應(yīng)的CSP空間濾波器。以一個帶有不同類別標(biāo)簽的4類運動想象腦電數(shù)據(jù)為例。該算法的具體過程表示如下：

設(shè)Xi∈RN×T,i=1,2,3,4分別表示4類運動想象原始EEG，其中N表示腦電數(shù)據(jù)集的通道數(shù)，而T表示每個通道的采樣點數(shù)，則可以得到EEG的歸一化協(xié)方差矩陣：

其中，是矩陣Xi的轉(zhuǎn)置矩陣，trace是矩陣對角線元素的總和。對每一類想象動作多次實驗得到的空間協(xié)方差矩陣取平均得到平均協(xié)方差矩陣：

其中，Ci表示第i類的試驗次數(shù)，表示所有第i類實驗協(xié)方差矩陣之和，則4類的平均協(xié)方差之和為對進(jìn) 行 特 征 值 分 解，構(gòu)建白化矩陣W1：

其中，Λ代表分解后的特征值矩陣，U0是與特征值矩陣相對應(yīng)的特征向量矩陣。把類別1分為一類，剩下的類別2、3、4合并為一類（計算其他類時，同理），則可將平均協(xié)方差矩陣轉(zhuǎn)換為：

可以證明S1和S2具有相同的特征向量，即Σ1=Y1Λ1YΤ1,Σ2=Y1Λ2YT1,Λ1+Λ2=Ι，其中I是單位矩陣。由于兩個矩陣相應(yīng)特征值之和總是1，因此對于S1具有最大特征值的特征向量對于S2具有最小特征值，反之亦然。因此該算法按照特征值的大小，將矩陣U1中的特征向量依次由大到小排列，然后分別選擇對應(yīng)于Λ1和Λ2中m個最大特征值的2m個特征向量，形成一個新的矩陣U′1，它們對于區(qū)分兩個類是最優(yōu)的。依據(jù)矩陣U′1和白化矩陣W1構(gòu)造類別1的空間濾波器：同理可得其余類別的空間濾波器W4，將每個運動想象腦電信號通過4個空間濾波器濾波后得到的4個信號組成新的信號：j∈1,…,n,其中n表示EEG總數(shù)。最后，將能量信號的方差規(guī)范化和取對數(shù)后為最終的特征：

1.3 LSTM

RNN模型的輸入由兩個輸入組成：一個是當(dāng)前時刻的輸入，另一個是前一時刻隱藏層的輸出。因此，RNN具有記憶功能，可以有效利用輸入信號的時間信息進(jìn)行時間序列分析。但在傳統(tǒng)的RNN中，在梯度反向傳播階段期間，梯度信號可能多次（與時間步數(shù)一樣多）乘以隱藏層神經(jīng)元之間連接的權(quán)重矩陣。如果矩陣中的權(quán)重很小，則可能導(dǎo)致消失梯度，在數(shù)據(jù)中學(xué)習(xí)長期依賴關(guān)系的信息也會變得更加困難。相反，如果這個矩陣中的權(quán)重很大，則可能導(dǎo)致梯度爆炸。

這些問題是提出LSTM模型的主要動機，它引入了記憶單元，如圖2所示。記憶單元由4個主要元件組成：具有自回歸連接的神經(jīng)元（與自身的連接）、輸入門、忘記門和輸出門。自回歸連接的權(quán)重為1.0，并且確保除外部干擾的情況下，記憶單元的狀態(tài)可以從一個時間步到另一個時間步保持恒定。輸入門可以允許輸入信號改變或保持記憶單元的狀態(tài)。另一方面，輸出門允許記憶單元的狀態(tài)對其他神經(jīng)元產(chǎn)生影響或阻止它。最后，忘記門可以控制記憶單元的自回復(fù)連接，允許單元根據(jù)需要記住或忘記其以前的狀態(tài)［17-18］。

設(shè)Xt為t時刻的輸入矢量，對于本文中使用的LSTM單元結(jié)構(gòu)，正向傳播的矢量公式如下：

圖2 LSTM記憶單元Fig.2 Memory units of long short-term memory(LSTM)

其中，邏輯Sigmoid函數(shù)σ用作門的激活函數(shù)，雙曲正切函數(shù)tanh通常用作記憶單元輸入輸出激活函數(shù)h，符號?表示按元素乘。

與輸入的時間序列x1,x2,…,xn類似，LSTM神經(jīng)網(wǎng)絡(luò)的輸出也是一個時間序列y1,y2,…,yn。針對運動想象EEG的分類問題，本文將輸出序列在所有時間步上進(jìn)行平均，得到最終的輸出y，并將y輸入一層全連接隱含層，其中隱含層的輸出激活函數(shù)是tanh函數(shù)，如圖3所示。

1.4 識別方法步驟

圖3 LSTM分類模型Fig.3 LSTM classification model

本文基于運動想象EEG的空間頻率與時間序列信息，提出一種多類運動想象EEG的特征提取方法，識別方法步驟如下：（1）首先將訓(xùn)練集原始EEG數(shù)據(jù)進(jìn)行離散小波變換提取運動想象任務(wù)相關(guān)子帶小波系數(shù)，接著利用OVR-CSP進(jìn)行特征提取。對n類運動想象任務(wù)構(gòu)建n個兩類空間濾波器，并將其組成一個多類空間濾波器。（2）將訓(xùn)練集原始EEG數(shù)據(jù)通過滑動矩形窗（長度為1 s，重疊50%）提取出多段EEG，并將這些EEG數(shù)據(jù)段按照時間序列組成時間序列EEG。然后，對時間序列EEG中的每段信號進(jìn)行離散小波分解，選取運動想象EEG相關(guān)頻帶的小波系數(shù)，將訓(xùn)練數(shù)據(jù)經(jīng)過多類空間濾波器濾波并計算特征，從而得到降維的DWT-CSP時間序列特征。（3）將訓(xùn)練集的DWT-CSP時間序列特征作為LSTM網(wǎng)絡(luò)的輸入，在LSTM之后再添加一層全連接隱含層，將隱含層的輸出輸入Softmax分類器，得到訓(xùn)練集的分類結(jié)果。采用沿時間反向傳播算法訓(xùn)練網(wǎng)絡(luò)，反向傳播包括兩個方向：一個是沿時間的反向傳播，即從當(dāng)前時刻開始，計算每個時刻的誤差項；另一個則是將誤差項向上一層傳播。（4）將測試集的EEG數(shù)據(jù)依次經(jīng)過DWT、多類空間濾波器濾波和訓(xùn)練好的LSTM網(wǎng)絡(luò)，從而得到時間步上平均的DWT-CSP-LSTM特征，特征經(jīng)Softmax分類器得到測試集的分類結(jié)果。

2 實驗及結(jié)果分析

2.1 實驗數(shù)據(jù)描述

本文采用2005年BCI競賽III-a數(shù)據(jù)集K3受試者腦電數(shù)據(jù)，該受試者執(zhí)行4類運動想象任務(wù)，分別為左手、右手、腳以及舌頭。實驗使用美國Neuro公司的Neuroscan腦電采集設(shè)備，該設(shè)備包含64導(dǎo)聯(lián)電極的腦電帽，以及采樣頻率為250 Hz的放大器。此外，該腦電采集設(shè)備將采集到的腦電數(shù)據(jù)以GDF格式存儲。

該數(shù)據(jù)集的實驗范式如圖4所示。受試者處于相對安靜的實驗環(huán)境，并以較放松的狀態(tài)坐在實驗椅上。為了防止由于信息反饋對系統(tǒng)產(chǎn)生影響，所有任務(wù)的顯示是隨機的。當(dāng)實驗開始后，受試者需要保持安靜，當(dāng)t=2 s時，BCI系統(tǒng)發(fā)出響聲提示受試者試驗正式開始，此時系統(tǒng)顯示器上出現(xiàn)一個“+”，提示受試者需要開始集中注意力。當(dāng)t=3 s時，顯示器上出現(xiàn)具有上下左右4個方向的箭頭，其中上下左右方向分別與左手、右手、腳、舌頭想象運動對應(yīng)，受試者按照不同方向箭頭做出相應(yīng)的運動想象。該數(shù)據(jù)集箭頭只保持顯示1 s時間，因而在t=4 s，該箭頭從顯示器上消失，但受試者仍需保持執(zhí)行運動想象任務(wù)直至t=7 s。

2.2 參數(shù)選取

圖4 實驗范式Fig.4 Timing of the paradigm

在基于空間頻率和時間序列信息的特征提取中，需要確定DWT的分解層數(shù)，CSP特征提取時每類空間濾波器所選取的特征向量數(shù)目，以及LSTM神經(jīng)網(wǎng)絡(luò)隱含層數(shù)與每層單元數(shù)等參數(shù)，這些參數(shù)的選取將直接影響本文特征提取的分類正確率。本文采用訓(xùn)練樣本數(shù)據(jù)5折交叉驗證的方法來優(yōu)化選取這些參數(shù)。

由于運動想象EEG的相關(guān)節(jié)律為α節(jié)律（8～12 Hz）和β節(jié)律（14～30 Hz），且該數(shù)據(jù)集的采樣頻率是250 Hz，本文選取Daubechies類db4小波，對每段EEG進(jìn)行4層小波分解，各層子帶頻帶范圍如表1所示。從表1可以看出D3和D4子帶的頻率范圍與運動想象EEG相關(guān)節(jié)律的頻率范圍相近，因而選取D3和D4子帶的小波系數(shù)組成新的腦電小波信號[ ]D3,D4。

OVR-CSP算法則是構(gòu)造多個空間濾波器，每個空間濾波器選擇其中一類作為第一類，其余類作為第二類。因而每個空間濾波器可選擇2m特征向量組成空間濾波器，或者選擇第一類最大的m特征值對應(yīng)的特征向量組成空間濾波器。K3受試者數(shù)據(jù)的分類正確率隨m或2m變化的實驗結(jié)果如圖5所示。從圖5可以看出，無論是由m或2m特征向量組成的空間濾波器，隨著m值從1增大到16，各受試者數(shù)據(jù)的分類正確率都是按先逐漸增大，達(dá)到最大值后開始降低的趨勢變化。但是由2m特征向量組成的空間濾波器在m=2時達(dá)到最大值92.23%，提取出的特征維數(shù)為16，而m特征向量組成的空間濾波器在m=7時達(dá)到最大值91.86%，提取出的特征維數(shù)為28，考慮特征維數(shù)和準(zhǔn)確率兩個因素，本文選取最優(yōu)特征向量數(shù)m=2的2m特征向量組成空間濾波器。

表1 腦電信號4層小波分解Tab.1 Frequency band of electroencephalogram(EEG)signals using 4-layer wavelet decomposition

圖5 K3受試者取不同m值時的分類準(zhǔn)確率Fig.5 Classification accuracy of subject K3 in different m

LSTM神經(jīng)網(wǎng)絡(luò)選擇合適數(shù)量的隱藏層沒有一般規(guī)則，具有少量神經(jīng)元的神經(jīng)網(wǎng)絡(luò)可能不足以對復(fù)雜函數(shù)建模。另一方面，具有太多神經(jīng)元的神經(jīng)網(wǎng)絡(luò)可能導(dǎo)致網(wǎng)絡(luò)過擬合，從而失去對測試集的泛化能力。找到最佳隱藏層數(shù)最常見的方法是實驗和試錯，不同隱含層層數(shù)和單元數(shù)的參數(shù)實驗結(jié)果如圖6所示。依據(jù)圖6的對比實驗結(jié)果，本文選擇單隱含層25個LSTM單元的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)進(jìn)行特征提取。

2.3 結(jié)果分析

圖6 LSTM單元數(shù)對分類精度的影響Fig.6 Effect of the number of LSTM units on classification accuracy

本文應(yīng)用基于空間頻率和時間序列信息的特征提取方法在BCI競賽III中Dataset III-a的K3受試者數(shù)據(jù)集上進(jìn)行4類運動想象任務(wù)分類測試實驗。圖7為CSP、DWT-CSP算法提取受試者K3的特征分布圖，為了對特征進(jìn)行可視化，本文只選擇由左手想象運動與其余想象運動形成的空間濾波器提取的一組二維特征，從圖7可知，單獨采用CSP算法從3名受試者腦電數(shù)據(jù)得到的特征類別間區(qū)分度較小。而經(jīng)過DWT-CSP方法得到的特征可以看出具有較為明顯的差異性，這說明頻率信息能夠增大類別間的區(qū)分度。當(dāng)采用更加完備的特征向量時，不同類別間的差異性會更加明顯。

為了進(jìn)一步探索空間、頻率以及時間序列信息對于運動想象EEG分類性能的影響，本文還對CSP特征、DWT-CSP特征以及CSP-LSTM特征與本文算法基于相同Softmax分類器的分類結(jié)果進(jìn)行對比，如圖8所示。其中，單獨采用CSP算法準(zhǔn)確率為78.33%，而結(jié)合空間與頻率信息的DWT-CSP特征的分類準(zhǔn)確率達(dá)到85.55%，結(jié)合空間與時間序列信息的CSP-LSTM特征的分類準(zhǔn)確率達(dá)到84.96%，這說明頻率或時間序列信息對于運動想象EEG分類的有效性，并且本文算法的分類準(zhǔn)確率達(dá)到了92.23%，進(jìn)一步說明空間、頻率、時間序列信息之間的互補性。

圖7 CSP、DWT-CSP算法提取受試者K3的特征分布圖Fig.7 Distribution of CSP and DWT-CSP features of subject K3

圖8 二維特征分布Fig.8 Two-dimensional feature distribution

最后，將本文所提出的方法與其他文獻(xiàn)中方法所得的分類準(zhǔn)確率進(jìn)行對比。文獻(xiàn)［19］結(jié)合CSP與CNN神經(jīng)網(wǎng)絡(luò)進(jìn)行特征提取，并采用Softmax分類器進(jìn)行分類，正確率達(dá)91.46%；文獻(xiàn)［20］采用CSP方法進(jìn)行特征提取，并應(yīng)用設(shè)計的BP神經(jīng)網(wǎng)絡(luò)對提取的特征數(shù)據(jù)進(jìn)行分類，正確率達(dá)76.67%；文獻(xiàn)［21］是基于CSP、Hilbert變換和SVM的特征提取與分類算法，準(zhǔn)確率達(dá)91.11%?？梢钥闯?，本文所提出的方法分類正確率比其他文獻(xiàn)有所提高，準(zhǔn)確率達(dá)92.23%。

3 結(jié)論

從非線性、非平穩(wěn)的運動想象EEG中提取特征是BCI系統(tǒng)的關(guān)鍵環(huán)節(jié)，針對在多類運動想象EEG分類時，CSP特征提取方法缺少頻率和時間序列信息的問題，本文將CSP方法與DWT、LSTM神經(jīng)網(wǎng)絡(luò)算法結(jié)合起來對多類運動想象EEG進(jìn)行特征提取，并采用Softmax分類器進(jìn)行任務(wù)分類。在2005年BCI競賽III-a的K3的4類運動想象腦電數(shù)據(jù)集上，相比單獨采用CSP、DWT-CSP、CSP-LSTM方法，融合空間、頻率、時間序列三者信息的DWT-CSP-LSTM方法達(dá)到了更好的分類結(jié)果，表明空間、頻率、時間序列信息的互補性和有效性，為BCI系統(tǒng)中EEG的特征提取提供一條新的思路。