基于多特征融合的運動想象腦電信號分類研究

陸振宇，陸旭峰，楊瑞洪，常珊

（1.江蘇理工學(xué)院電氣信息工程學(xué)院生物信息與醫(yī)藥工程研究所，常州213001；2.揚州工業(yè)職業(yè)技術(shù)學(xué)院，揚州225127）

0 引言

腦機接口（Brain Computer Interface，BCI）是綜合神經(jīng)科學(xué)、信號處理、機器學(xué)習(xí)和信息技術(shù)等領(lǐng)域的跨學(xué)科技術(shù)。BCI 的目的是把大腦活動轉(zhuǎn)變成對設(shè)備的控制指令，具有恢復(fù)失去的感官和運動功能的潛力，如治療帕金森病的腦深度電刺激，截肢者和脊髓損傷病人的假肢、輪椅，還有中風和肌肉萎縮患者使用的鼠標和單詞拼寫設(shè)備[1]。BCI 已經(jīng)在游戲娛樂、機器人、生物識別和教育等應(yīng)用領(lǐng)域表現(xiàn)出強大的發(fā)展前景。

目前，運用較為廣泛的信號是當受試者自主想象某種特定動作所產(chǎn)生的神經(jīng)活動，稱之為基于運動想象（Motor Imagery，MI）的腦電信號。研究表明，當受試者進行手、腳、舌頭動作的運動想象，大腦皮層相對應(yīng)區(qū)域的μ頻帶（8-13Hz）表現(xiàn)顯著，而對側(cè)區(qū)域其節(jié)律能量下降，該現(xiàn)象被稱為事件相關(guān)同步/去同步（ERD/ERS）現(xiàn)象，這也為信號的特征提取提供了依據(jù)[2-4]。

BCI 的實現(xiàn)主要涉及信號處理和模式識別。信號處理包括數(shù)據(jù)的預(yù)處理，如濾波、偽跡去除，還有關(guān)鍵的特征提取技術(shù)。文獻中的特征提取方法一般有，頻域分析有自回歸（AR）模型、功率譜估計；時頻分析有小波變換、小波包變換、局部均值分解；空間域有主成分分析（PCA）、獨立分量分析（ICA）和共空間模式（CSP），用于多導(dǎo)聯(lián)的腦電信號[5-11]。模式識別是利用機器學(xué)習(xí)技術(shù)，對輸入特征進行分類。常用的二分類技術(shù)有線性判別分析（LDA）、神經(jīng)網(wǎng)絡(luò)（NN）、支持向量機（SVM）、基于隨機森林和Boosting 的集成技術(shù)。

通過對相關(guān)文獻的研究，信號變換和分解的方法多樣，一般在信號處理后選取能量、熵等作為特征，但變換方法及提取的特征單一化，無法多方面地表示腦電信號。本文選取提取算法易于理解、計算量小、效果良好的AR 模型、小波包變換、CSP 方法對信號進行處理，提取信號的多方面特征，采用SVM、LDA 分類器進行模式識別，并對分類結(jié)果評價。

1 實驗數(shù)據(jù)

本文實驗的數(shù)據(jù)是來源于2003 年BCI 競賽的Graz 數(shù)據(jù)集，由奧地利Graz 科技大學(xué)醫(yī)學(xué)系、生物醫(yī)學(xué)工程研究所提供[12]。實驗范式是受試者放松地坐在手扶椅上，根據(jù)屏幕顯示做出相應(yīng)的反饋。前2s 為準備時間，在第2s 時，蜂鳴器提示實驗開始，同時屏幕出現(xiàn)1s 的“+”符號，3-9s 是受試者按給出的箭頭方向進行對應(yīng)的運動想象。實驗采集了C3、CZ、C4 通道的腦電信號，采樣頻率為128Hz，并通過0.5-30Hz 的帶通濾波。整個數(shù)據(jù)集共有280 個樣本，其中左右手運動想象各140，訓(xùn)練和測試樣本各占一半，數(shù)據(jù)集分布均勻。

2 特征提取算法

2.1 AR模型

AR 模型可用如下公式表示：

式中，ε（n）是均值為零、方差為σ2的白噪聲，p 為模型階次，ap（i）為自回歸模型參數(shù)[4-5]。本文選取競賽數(shù)據(jù)C3、C4 通道，時間段在4-7s 的信號，基于最終預(yù)測誤差準則，確定階次為5，使用Burg 算法提取其AR 模型系數(shù)和系數(shù)差作為特征。

2.2 小波包變換

小波包變換是為了克服小波分解在高頻段的頻率分辨率較差，而在低頻段表現(xiàn)較好的基礎(chǔ)上提出的。它是一種更精細的信號分析的方法，提高了信號的時域分辨率[13]。

競賽數(shù)據(jù)采樣頻率為128Hz，依據(jù)Nyquist 采樣定理，腦電信號頻率范圍為1～64Hz，使用小波基函數(shù)對C3、C4 信號進行4 層小波包分解，分解頻寬為4Hz。分解后的第4 層小波包節(jié)點頻率范圍如表1 所示，具有代表性的μ節(jié)律信號為節(jié)點（4,3）。圖2 為采用bior3.3 小波基，對單個樣本的C3、C4 通道信號做小波包變換得到的時域和頻域圖。

表1 第4 層小波包分解的節(jié)點及頻帶范圍

圖2 單樣本節(jié)點（4,3）的bior3.3小波包變換

腦電信號經(jīng)小波包分解后，節(jié)點（i,j）信號的小波包能量由計算該節(jié)點的小波包系數(shù)平方和得到。

式中，i 為分解層數(shù)，j 為該層下的節(jié)點，di,j（n）為節(jié)點（i,j）的小波包系數(shù)。

根據(jù)能量守恒定律，分解后的各節(jié)點能量總和等于原信號的能量，則相對小波包能量為節(jié)點信號能量與該層節(jié)點總能量的比值：

式中，Eraw為原信號能量，Etotal為原信號分解i 層后各節(jié)點能量總和，pi,j則為節(jié)點（i,j）能量的歸一化，稱為相對能量。

信息熵可以表征信號的不確定度，結(jié)合小波包變換，可得出分解后的節(jié)點信號小波包熵為：

本文選取C3、C4 通道信號做小波包變換，以一次實驗所記錄的C3 信號為例，具體步驟如下：

（1）采用固定時間窗的方法，以128 點的時間窗，將長度為1152 點的信號分割成9 段；

（2）在一段128 點的信號上做小波包分解，小波基選擇bior3.3[4]，分解層數(shù)為4，選擇節(jié)點（4,3）的重構(gòu)信號，由式2 計算出小波包能量，式4 得出小波包相對能量，式5 得出小波包熵。

（3）C4 通道重復(fù)1、2 步驟，一個樣本的小波包變換結(jié)果為2×9 的矩陣，按樣本類別進行疊加平均，得到C3、C4 通道腦電信號的小波包能量、相對能量和熵的變化，圖3 為右手運動想象下的平均小波包能量和熵走勢圖。

圖3 右手運動想象的bior3.3平均小波包變換特征

由圖3 可知，4-7s 的C3、C4 通道小波包相對能量和小波包熵區(qū)分度明顯，走勢基本一致，數(shù)值大小不一樣。本文將其二維向量特征轉(zhuǎn)化成一維向量，一個樣本C3、C4 通道的小波包能量差、相對能量差、熵差△WPE 作為分類的特征，圖4 為小波包熵差的走勢圖。

圖4 左右手運動想象的bior3.3平均小波包熵差

2.3 共空間模式

CSP 是尋找空間濾波器，使得濾波處理后的數(shù)據(jù)與其中一類方差達到最大，而與另一類方差達到最小的算法?；谶\動想象的BCI，依賴特定的頻帶能量，對于兩類的左右手運動想象信號，經(jīng)過CSP 處理得到的特征向量增強了兩類的區(qū)別，可區(qū)分程度最大化[14]。

單次腦電信號樣本可由Ns×Nc的二維矩陣E 表示，其中，Ns是采樣點數(shù)，Nc是通道數(shù)，則歸一化協(xié)方差矩陣為：

式中ET表示E 的轉(zhuǎn)置，表示矩陣的秩。按樣本類別求疊加平均的協(xié)方差矩陣，然后總樣本的協(xié)方差矩陣為：

式中Cov1、Cov2 分別表示類別1、2 的協(xié)方差矩陣。對總樣本協(xié)方差矩陣作特征值分解，得到特征向量矩陣U、特征值對角陣λ：

將U、λ進行降序排列，求得白化矩陣：

對類別1 的協(xié)方差矩陣使用白化矩陣進行轉(zhuǎn)變：

對其作特征值分解，得到特征向量矩陣U1、特征值對角陣λ1：

同樣對U1、λ1進行降序排列，得到類別1 的投影矩陣W：

取W 的前m 行、后m 行，構(gòu)成空間濾波器W’，原信號E 經(jīng)過空間濾波器得到：

對E’按行求方差，并取對數(shù)，得到特征值f：

一般地還有，在方差的歸一化后取對數(shù)，得到特征值f：

本文選取4-7s 的運動想象腦電信號，對其做CSP變換，按式（15）、（16）分別計算，各自得到2 個特征值，用于區(qū)分左右手運動想象。

圖5 左右運動想象各70個樣本的特征值散點圖

3 模式判別

3.1 線性判別分析

線性判別分析（LDA）是線性二分類器，將輸入向量映射到一個超平面，該超平面將輸入空間劃分為兩個半空間：每個空間代表一種類別。決策邊界由超平面公式?jīng)Q定：

式中w 為法向量，w0為閾值，均由訓(xùn)練數(shù)據(jù)得出。對于新的輸入特征向量，代入式16，如果g（x）≥0，為類別1，反之為類別2。

3.2 支持向量機

支持向量機（SVM）是尋找兩類樣本間隔最大的分離超平面，以獲得最好的泛化能力的一種分類器。本文采用libsvm 工具，選擇C-SVM 類型，利用徑向基核函數(shù)實現(xiàn)腦電信號特征的非線性映射，使用網(wǎng)格搜索獲得最優(yōu)參數(shù)C、gama。

4 實驗結(jié)果分析

為了驗證本文提出方法的有效性，將單特征和多特征融合提取方法運用在整個數(shù)據(jù)集上，采用SVM、LDA 分類器對其分類，結(jié)果如表2、3，其中識別率是10折交叉驗證的平均值。

通過表2 可以直觀地發(fā)現(xiàn)：單特征選取方面，AR系數(shù)和系數(shù)差都能較好地被區(qū)分，識別率最高為87.5%；小波包熵差的效果高于能量差、相對能量差，原因在于小波包熵的數(shù)量級比后兩者大，從圖3 可知，小波包能量的范圍在0.25～0.75，而小波包熵是1.2～3.0；而共空間模式下的兩種特征在分類器上表現(xiàn)了較大差異，未歸一化的特征在LDA 上的表現(xiàn)最好，原因在于數(shù)據(jù)的數(shù)值較大，而對于SVM，無論特征是否歸一化，準確率都是73.93%，表現(xiàn)一般。

通過表2 的分析，選取分類效果好的AR 模型系數(shù)、小波包熵差、未歸一化的CSP 特征，進行排列組合，用于多特征融合，得到表3 的結(jié)果。其中，分類器方面，LDA 的識別率均比SVM 略高；特征融合方面，三者組合的準確率最高，達到了91.43%，而兩兩組合的準確率比各自的單特征高，均在90%左右。

表2 單特征的平均識別率（單位：%）

表3 多特征融合的分類準確率（單位：%）

5 結(jié)語

本文提出了多方面對運動想象腦電信號進行特征提取及融合，并使用分類器對其評價。AR 模型代表了腦電信號的功率譜估計，僅使用C3、C4 兩通道信號的系數(shù)，也可以得到較好的分類結(jié)果；小波包變換是定位到給定頻帶，計算熵值可以反映其能量的差異性，選擇C3、C4 通道信號的熵差作為特征，也得到了較好的效果；CSP 算法是提取腦電信號的空間特征，適用于較多導(dǎo)聯(lián)，且需要訓(xùn)練數(shù)據(jù)構(gòu)造濾波器，本文的競賽數(shù)據(jù)是3 導(dǎo)聯(lián)，使用CSP 僅得到兩個特征，分類效果不穩(wěn)定。而多特征融合的方式，則彌補了單特征的缺陷，在分類效果上也更好，這為運動想象腦電信號的特征工程提供了思路。在分類器的選擇上，SVM 雖然泛化能力高，但計算量大，對數(shù)據(jù)敏感，而LDA 實現(xiàn)簡單，運算速度快，效果也較好，適用于在線BCI。