在基于頭皮腦電(EEG)信號的腦-機接口(BCI)研究中, 用戶個體差異性和背景噪聲的復雜性是影響BCI系統穩定性的兩個主要因素。因此需要針對不同個體進行BCI系統參數優化, 其中包括對時域、空域濾波器參數的優化設計和分類器參數的學習。本文以提高BCI系統的準確性為目標, 提出了一種結合獨立分量分析空域濾波器(ICA-SF)優化設計和EEG多子帶特征的BCI信息處理新方法。基于所提方法, 對4位受試者在不同時間采集的三類運動想象EEG(MI-EEG)進行分析。實驗結果表明, 在同一受試者的自交叉測試和不同受試者數據集之間的互交叉驗證中, 多子帶特征結合方法所得到的平均識別率比僅使用單頻帶所得的平均識別率普遍提高, 識別率最大提升可達6.08%和5.15%。
引用本文： 康莎莎, 周蚌艷, 吳小培. 基于獨立分量優化子帶特征的三類運動想象分類. 生物醫學工程學雜志, 2016, 33(2): 208-215. doi: 10.7507/1001-5515.20160037 復制
引言
腦-機接口(brain-computer interface，BCI)是一種不依賴正常的腦神經和肌肉輸出通路的腦-機通信技術，它能為思維功能正常的肢體殘疾群體提供一種特殊的人-機交互模式，相關研究具有重要的科學意義和應用價值[1-3]。作為神經活動的信息載體，頭皮腦電(electroencephalogram，EEG)信號能實時反映思維狀態的變化，因此在非植入式BCI系統中得到了廣泛的應用。然而多通道EEG信號的空間分辨率和信噪比都非常低，并且EEG信號特征的個體差異性也非常明顯，這些不利因素給EEG信號的有效分析和精確解碼帶來了很大困難。EEG信號處理和特征獲取技術是BCI系統實現的關鍵，相關研究已經成為該領域普遍關注的焦點。本文以運動想象BCI(motor imagery BCI，MI-BCI)的系統實現為基本目標，重點研究基于獨立分量分析(independent component anaylsis，ICA)的運動想象EEG(motor imagery EEG，MI-EEG)特征提取新方法。當受試者進行肢體運動想象(motor imagery，MI)時，BCI系統通過對頭皮EEG信號的實時處理和特征獲取，實現MI類型的識別和控制命令生成。由于EEG信號中運動相關節律成分的信噪比非常低，并且其生成過程對環境因素和受試者自身狀態均比較敏感，因此EEG信號動態分析和處理方法非常關鍵[4]。在現有的EEG信號處理方法中，綜合利用多道EEG信號時-頻-空三域信息的空域濾波技術被認為是最有效的EEG信號處理方法[5]。目前最常用的空域濾波方法有ICA和共同空間模式(common spatial pattern，CSP)[6-14]。其中ICA在多通道EEG信號分析中非常有效，因為ICA空域濾波器(ICA spatial filter，ICA-SF)不受信號頻譜混迭的限制，濾波的同時能有效保護有用的細節成分，很大程度上彌補了傳統時域、頻域方法的不足[15-16]。
本文以提高BCI系統的識別率為目標，提出了一種結合ICA-SF優化設計和EEG多子帶特征的BCI信息處理新方法，具體內容包括：①基于ICA-SF和支持向量機(support vector machine，SVM)分類器分析單頻帶和多子帶下三類MI-EEG的分類效果，在此基礎上，實現了針對個體的運動相關子帶優化組合；②基于4位受試者的MI-EEG數據集，對ICA-SF的有效性和多子帶特征組合方法進行驗證，得到了較高的分類識別率。
1 運動想象試驗
根據神經電生理知識，肢體MI能誘發運動神經皮層mu和beta節律的功率譜變化，這一現象被稱之為事件相關去同步(event-related desynchronization，ERD)和事件相關同步化(event-related synchronization，ERS)[17-19]。ERD/ERS現象是BCI系統特征分類的一個非常重要的神經電生理特征依據，其中不同類型MI誘發的ERD/ERS現象在運動皮層具有不同的神經活動分布。當想象左手或右手運動時，肢體對側運動皮層(C4/C3電極位置)的mu/beta節律出現抑制(ERD)，而同側運動皮層(C3/C4電極位置)的mu/beta節律則相對增強(ERS)。當進行想象腳運動時，中央運動皮層(Cz電極位置)的mu/beta節律出現抑制(ERD)。
為了緩解受試者疲勞度，提高采集的MI-EEG數據質量，本文的MI試驗設計在BCI Competition 2005 dataⅢa(BCI2005)試驗范式的基礎上進行了改進，具體包括以下三個方面：
(1)?減少了導聯數。其目的是減少各電極之間的干擾和增加受試者舒服度，相比BCI2005的60導聯MI-EEG數據，本文的試驗設計在保證實驗效果的同時僅采集了16導聯數據，其中主要包含運動皮層相關的導聯，如圖 1所示。
(2)?延長了休息時間。將BCI2005試驗中2 s的休息時間延長到4 s，用以緩解受試者疲勞度，提高數據質量。否則受試者會因不適而做出一些較大幅度的動作，致使電極松動或產生其他干擾信號(眼電、肌電等)。
(3)?減少了每組數據中的單次試驗(single trial)數量。本文中的每組試驗只采集了75個單次試驗。單次試驗次數的適當減少可保證MI-EEG數據的整體可靠性。
圖1
自主實驗的電極分布(上)和試驗設計(下)
Figure1.
Electrode distribution (up) and experimental design of independent experiment (down)
EEG信號采集設備為NeuroScan 40導放大器，依照國際標準10-20系統，采用16導電極，分別為VEOU，VEOL，Fp1，Fp2，FC3，FCZ，FC4，C3，CZ，C4，CP3，CPZ，CP4，O1，OZ，O2，參考電極為A1，A2。采樣頻率為250 Hz，數據精度是16 bit，采集時對原始EEG進行了50 Hz陷波濾波和0.5~100 Hz帶通濾波。受試者為實驗室4位研究生S1(女，25歲)、S2(男，29歲)、S3(女，23歲)和S4(女，22歲)，身體健康，均為右利手。實驗環境無電磁屏蔽。單次試驗持續10 s，大致分為三個時段：0~1 s的提示，1~6 s的MI和6~10 s的休息。具體過程是：每單次試驗開始，計算機發出一個短促的“滴”提示聲，提醒受試者要保持注意力集中，1 s后計算機屏幕顯示方向箭頭：向右、向下或向左的方向箭頭，提示受試者要進行右手、腳或左手的MI，每類MI持續5 s，第6秒箭頭消失，此后4 s的持續時間計算機保持黑屏，受試者則處于休息放松狀態，安靜地等待下一個單次試驗開始。其中箭頭方向出現的順序是隨機的，每類MI各25次，每組MI-EEG數據中含75個單次試驗，每次試驗采集A、B兩組MI-EEG數據。需要說明的是，綜合考慮運算效率和識別率，本文只采用了16導聯中位于運動皮層附近的9導聯MI-EEG數據(FC3、FCZ、FC4、C3、C4、CZ、CP3、CPZ和CP4)[20]。
2 方法
2.1 ICA算法
ICA是一種基于高階統計信息的盲源分離方法。該算法假定源是非高斯分布，分離后的源信號具有最大的獨立性，ICA基本原理如圖 2所示。由于EEG信號具有非高斯性，且各個局部神經皮層產生的神經活動可近似認為相互獨立，所以可利用ICA算法從頭皮EEG信號分離出具有明確神經電生理意義的神經活動[21]。
圖2
ICA算法原理圖
Figure2.
Principle chart of ICA algorithm
假定n路混合信號x(t)=[x1(t), …, xn(t)]T是由多個隱含EEG源信號s(t)=[s1(t), …, sn(t)]T的線性瞬時混合而成，如式(1)所示：
| $ \boldsymbol{x}\left(t \right)=\boldsymbol{A}\boldsymbol{s}\left(t \right) $ |
式(1)中混合矩陣A和源信號s(t)未知。與式(1)對應的分離模型如式(2)所示：
| $ \boldsymbol{u}\left(t \right)=\boldsymbol{W}\boldsymbol{x}\left(t \right)=\boldsymbol{W}\boldsymbol{A}\boldsymbol{s}\left(t \right) $ |
其中，W是分離矩陣，u(t)=[u1(t), …, un(t)]T是對源信號s(t)的近似估計，且u(t)中的各個分量間要盡可能相互獨立。本文用于多通道MI-EEG數據分析的ICA算法是實驗室自編Matlab代碼[22]，算法中采用的是信息極大獨立性判決和自然梯度優化算法[8-9]。
2.2 ICA-SF設計
在MI-BCI系統實現中，一個好的空域濾波器對提升MI-EEG信噪比和改善BCI分類效果起著非常關鍵的作用。本文采用實驗室近期提出的基于ICA空域濾波的三類運動相關獨立分量(motor-related independent components，MRICs)檢測新方法，具體實現流程如圖 3所示。
圖3
MRICs檢測濾波器流程圖
Figure3.
Flow chart of MRICs detection filter
運行ICA算法之前，對原始MI-EEG數據進行8~30 Hz數字帶通濾波，該頻段包含了運動相關的mu和beta節律成分[18]。此外，沒有進行其他偽跡處理。多通道單次MI-EEG數據經ICA方法后分解成多個獨立分量(independent component，IC)，其代表MI-EEG各個隱含源。同時得到分離矩陣W，其列向量則為ICA-SF，而混合矩陣A的列向量則表示IC在全部頭皮電極的空間模式分布(投影強度)。文獻[21]給出的基于混合矩陣A的三類ICA-SF有效性判定，可初步剔除不符合要求的單次訓練樣本，保證最終用于BCI系統實現的ICA-SF檢測性能。該判據準則的基本思路是：由訓練集的單次MI-EEG數據所得的混合矩陣A中，必須包含類似如圖 4所示三類MRICs的空間分布模式的三個列向量。如果滿足這一條件，則根據這3個列向量在A中的位置，選擇分離矩陣W中相應的列作為有效的三類MRICs檢測濾波器組：{wl(i) wr(i) wf(i), i=1, …, M≤N}，l，r和f分別代表左手、右手和腳，M為有效濾波器個數，N為訓練集中單次試驗數。在此基礎上，將所得有效濾波器組分別用于全部訓練集的單次MI-EEG進行MI標簽估計，然后根據識別率高低進一步判斷各濾波器組的性能。
圖4
三類MRICs的空間模式分布
Figure4.
Spatial distribution pattern of three-class of MIRCs
通過下面的例子可進一步說明本文采用的ICA-SF設計思想。文中的數據采集試驗每組包含75個單次試驗。根據如圖 3所示算法流程，理想情況下可得到75個基于單次MI-EEG的ICA-SF和相應的識別率。如果單次MI-EEG數據存在采集質量不高的“壞”數據，那么所設計的ICA-SF會失效或導致識別率較低，如圖 5所示。圖中存在識別率為零的情況，表明此次ICA-SF設計失敗；另外，還有一些識別率很低的情況，同樣表明對應的ICA-SF檢測性能較差，此兩種情況對應的單次MI-EEG數據不能用于設計ICA-SF。
根據已報道的MI-BCI研究論文，未見有對訓練樣本質量進行評估和選擇的討論，通常是采用全部單次訓練樣本的拼接用于ICA-SF的設計。根據本文的前期研究工作，如果不加選擇地將訓練集中一些“壞”的單次數據用于設計ICA-SF，會導致其性能下降。另外，從統計理論上看，雖然大數據樣本有利于保證ICA-SF設計性能的穩定性，但是也會使計算過程過于緩慢，不利于在線BCI的實現。因此，本文依據單次識別率大小，從合格的單次試驗數據中，選擇前L個(L的經驗值一般為10)單次MI-EEG數據拼接成用于ICA-SF設計的訓練樣本，最終得到用于MI-BCI系統實現的三類MI濾波器組：[wl wr wf]。
2.3 多子帶特征結合方法設計
多道MI-EEG經空域濾波后，通常需要選擇最優運動相關頻帶對運動相關節律進一步增強，但不同受試者的最優運動節律頻帶會因個體差異性而存在較大差別。如圖 6所示是S1和S2某組MI-EEG數據在C3、CZ和C4導聯的功率譜，其中S1在子帶10~13 Hz有明顯的譜峰，而S2在子帶10~14 Hz和19~23 Hz的譜峰都比較明顯。如果手動選擇S2的最優濾波子帶可能要經過多次嘗試才能找到，這顯然不利于實用MI-BCI系統的設計。
圖6
S1和S2與運動相關導聯對應的功率譜
Figure6.
Corresponding power spectrum of motor-related of S1 and S2
為了進一步提高MI-BCI的實用性和識別率，本文提出基于ICA-SF的多子帶特征結合方法，如圖 7所示。該方法的三類MI-EEG數據的識別步驟如下：
圖7
多子帶特征結合方法流程圖
Figure7.
Flow chart of multiple sub-band feature combination
(1)?基于如圖 3所示算法對ICA-SF設計樣本進行優化選擇；
(2)?選取識別率最高的前10個單次MI-EEG數據進行拼接，然后用來設計最終的左手、右手和腳的ICA-SF；
(3)?最優運動相關頻帶選擇：把8~30 Hz濾波頻段分為5個子帶，利用單頻帶ICA-SF方法對每個子帶進行節律增強，每個子帶下都能得到相應的分類結果[23-24]；
(4)?選取識別率高的若干子帶的空域特征進行組合，以組合后的特征向量訓練SVM分類器；
(5)?測試集輸入訓練的SVM分類器完成最后的分類識別。
3 特征提取和分類
根據如圖 7所示方法，ICA-SF從原始MI-EEG中檢測出與MI相關的MRICs成分，然后經過5個子頻帶濾波器，輸出相應的MRICs子帶成分。針對每個子帶輸出，建立相應的特征向量用作SVM分類器的訓練集和分類集。為了避免過擬合現象，采用了5×5交叉驗證。特征向量建立過程如下：
假設x為原始多通道MI-EEG數據，首先對x進行8~30 Hz帶通預處理，消除部分干擾成分；然后利用ICA算法優化設計的濾波器組[wl wr wf]提取三類MRICs成分，即：zi=wix，其中i=l, r, f；特征向量取之于zi的MI時間段，通常選擇MI提示后0.5~5 s的MI-EEG數據段[22]。在已報道的文獻中，特征向量通常為式(3)所示的標準化方差[25-27]。但是，單一方差特征通常不能完整地描述MI-EEG數據的多域特征，會導致三類MI識別率不高。因此本文采用式(4)所示的特征向量組，即直接將三個MRICs的MI時段(0.5~5 s)的樣本組合成特征向量，為了降低特征維數，可對特征樣本進行4~5倍下采樣。
| $ {F_i}=\frac{{{\mathop{\rm var}} \left({{z_i}} \right)}}{{\sum\limits_{i=l, r, f} {{\mathop{\rm var}} \left({{z_i}} \right)} }}, {F_{feature}}=[{F_l}\; {F_r}\; {F_f}] $ |
| $ {F_{feature}}=[{z_l}\; {z_r}\; {z_f}] $ |
4 實驗結果
本文數據集包含S1的4次8組數據和S2、S3、S4的各2次4組數據，共有20組MI-EEG數據集。每組有效數據點為187 500(250×10×75)。因此，S1的數據樣本總數為8×187 500，其他三位受試者S2、S3和S4的數據樣本總數均為4×187 500。為了驗證所提算法的有效性，論文基于20組MI-EEG數據集進行了兩種形式的交叉測試，即同一受試者數據集之間的自交叉測試和不同受試者數據集之間的互交叉驗證，為了表述方便，文中分別稱之為“SS_TEST”和“DS_TEST”。在用于學習分類的特征選擇方面，采用識別率最高和次高的兩個子帶特征進行組合，以驗證ICA-SF優化設計和多子帶特征結合方法的有效性。
4.1 同一受試者的自交叉測試(SS_TEST)
分別對4位受試者的20組數據進行SS_TEST，比較采用MRICs單頻帶和多子帶特征對BCI識別率的影響。以受試者S1的數據為例，給出一些代表性的識別率測試結果，最后給出全部受試者的SS_TEST結果。
受試者S1在4次MI試驗中共采集的8組MI-EEG數據集，鑒于同一受試者在同一次實驗采集中導聯位置誤差較小，因此同一次實驗采集的A、B兩組MI-EEG數據集間交叉測試識別率可能較高。為了驗證是否存在此類現象，本小節把S1的4次試驗的8組MI-EEG數據進行了分組，即分成A、B兩大組，A代表每次試驗的A組數據，B為每次試驗的B組數據。
S1-A和S1-B兩大組的SS_TEST結果如圖 8所示，橫坐標指明了各個對應的訓練集。在S1-A的SS_TEST中，訓練集S1-1-A對應的識別率描述為：在單頻帶和兩個子帶下，第1、2列是訓練集S1-1-A測試同一次采集的B組(S1-1-B)數據而得到的兩個識別率(分別對應采用單頻帶和兩子帶特征)。接下來的第3~8列是訓練集S1-1-A對隔天采集的A組其他MI-EEG數據(S1-2-A、S1-3-A、S1-4-A)的測試分類結果，類似的SS_TEST在S1_B組數據集中同樣進行。通過比較可以看出，基于兩子帶聯合特征的三類MI識別率總體上要高于基于單頻帶特征的識別率。但是，對于同一個受試者而言，無論兩組數據集是否來自同一次數據采集試驗，SS_TEST識別率并沒有體現出明顯的規律性差別。為了更直觀地比較，圖 9給出了S1的8組MI-EEG數據分別作為訓練集的交叉測試結果平均值，從中明顯可以看到，基于多子帶特征的MI識別率明顯高于基于單頻帶特征的識別率。
圖8
S1作為訓練集的SS_TEST識別率
Figure8.
Recognition rates of training set S1 by SS_TEST
圖9
訓練集S1的SS_TEST平均識別率
Figure9.
Average recognition rates of training set S1 by SS_TEST
如表 1所示，給出了對其他3位受試者數據集測的全部SS_TEST結果。由于S2~S4只含有兩次實驗的4組數據。因此，在交叉測試時，將各受試者數據分為｛S2-1，S2-2｝、｛S3-1，S3-2｝、｛S4-1，S4-2｝。然后在單頻帶和聯合子帶特征兩種情況下進行SS_TEST。實驗結果顯示，對于不同受試者的SS_TEST試驗，兩個子帶下MI識別率均比單頻帶下識別率要高。其中識別率最高提升為6.08% (S4-1)，最低為1.45%(S4-2)，平均提升了3.64%。
表1
S2、S3和S4作為訓練集的SS_TEST的平均識別率(%)
Table1.
Average recognition rates of training sets S2, S3 and S4 by SS_TEST (%)
4.2 不同受試者間的互交叉測試(DS_TEST)
為了更好地說明本文方法的有效性，此小節在不同受試者之間進行了組間交叉測試(DS_TEST)。為了表述方便，重新統一了MI-EEG數據集：S1、S2、S3和S4各2次4組試驗數據，共16組。由于受試者個體存在一定的差異性，不同受試者的最佳MI頻段也有所不同。所以在DS_TEST中，針對不同受試者使用的最優濾波頻帶也需適當調整。
如表 2所示是S1作為訓練集，S2、S3和S4作為測試集的DS_TEST分類結果，其結果表明基于聯合子帶特征的MI識別率比單頻帶下的識別率要高。如，在訓練集S1-3-A下，單頻帶下的組間測試平均識別率為77.2%，兩個子帶下的平均識別率為80.8%，相比單頻帶下的平均識別率高出3.65%；而在測試集S1-3-B下，兩個子帶下的平均識別率更是比單頻帶下的平均識別率高出5.15%。
表2
S1作為訓練集的DS_TEST識別率(%)
Table2.
Recognition rates of training set S1 by DS_TEST (%)
如圖 10所示是S2、S3和S4作為訓練集，S1作為測試集的左手、右手和腳的MI分類結果。圖中受試者的平均識別率能清晰地分析出兩種頻帶方式下分類結果的變化規律，本小節得到了與4.1小節的相同結論：三類MI在單頻帶下的平均識別率沒有兩個子帶下的平均識別率高，驗證了本文方法的可行性和有效性，為進一步的BCI系統提供了一種頻帶優化方法。
圖10
S1作為測試集的DS_TEST平均識別率
Figure10.
Average recognition rates of testing set S1 by SS_TEST
5 總結與討論
在基于空域濾波的MI-BCI研究中，由于ERD/ERS現象的空間分布特性，以及EEG信號在不同受試者身上所體現的高度復雜性與變異性，針對不同受試者的特征優化選擇對于MI-BCI系統的實現尤為重要。雖然MI-EEG數據自動特征選擇是BCI系統實現的基本任務之一，但是一些傳統方法的穩定性和適用性依然有需要改進之處。在與本文相關的研究中發現，多通道MI-EEG數據的空域濾波器技術對BCI系統整體性能起到非常關鍵的作用，雖然ICA空域濾波方法理論上能夠提取出真實的運動相關神經活動，但是ICA方法的穩定性和對隨機干擾的敏感性很大程度上限制了該方法的性能。
本文從ICA-SF的優化設計著手，盡可能保證MRICs檢測濾波器的性能，在此基礎上實現運動節律子帶的優化組合。獲得了較好的三類MI識別結果。本文僅選擇了最優和次優的兩個子帶進行驗證，一方面是出于對算法計算復雜度的考慮，另外，更多子帶的優化組合效果相比兩個子帶而言，對識別率改善并不十分明顯。原因可能是與子頻帶劃分策略有一定的關系，這也是本文后續研究的重點內容之一。
引言
腦-機接口(brain-computer interface，BCI)是一種不依賴正常的腦神經和肌肉輸出通路的腦-機通信技術，它能為思維功能正常的肢體殘疾群體提供一種特殊的人-機交互模式，相關研究具有重要的科學意義和應用價值[1-3]。作為神經活動的信息載體，頭皮腦電(electroencephalogram，EEG)信號能實時反映思維狀態的變化，因此在非植入式BCI系統中得到了廣泛的應用。然而多通道EEG信號的空間分辨率和信噪比都非常低，并且EEG信號特征的個體差異性也非常明顯，這些不利因素給EEG信號的有效分析和精確解碼帶來了很大困難。EEG信號處理和特征獲取技術是BCI系統實現的關鍵，相關研究已經成為該領域普遍關注的焦點。本文以運動想象BCI(motor imagery BCI，MI-BCI)的系統實現為基本目標，重點研究基于獨立分量分析(independent component anaylsis，ICA)的運動想象EEG(motor imagery EEG，MI-EEG)特征提取新方法。當受試者進行肢體運動想象(motor imagery，MI)時，BCI系統通過對頭皮EEG信號的實時處理和特征獲取，實現MI類型的識別和控制命令生成。由于EEG信號中運動相關節律成分的信噪比非常低，并且其生成過程對環境因素和受試者自身狀態均比較敏感，因此EEG信號動態分析和處理方法非常關鍵[4]。在現有的EEG信號處理方法中，綜合利用多道EEG信號時-頻-空三域信息的空域濾波技術被認為是最有效的EEG信號處理方法[5]。目前最常用的空域濾波方法有ICA和共同空間模式(common spatial pattern，CSP)[6-14]。其中ICA在多通道EEG信號分析中非常有效，因為ICA空域濾波器(ICA spatial filter，ICA-SF)不受信號頻譜混迭的限制，濾波的同時能有效保護有用的細節成分，很大程度上彌補了傳統時域、頻域方法的不足[15-16]。
本文以提高BCI系統的識別率為目標，提出了一種結合ICA-SF優化設計和EEG多子帶特征的BCI信息處理新方法，具體內容包括：①基于ICA-SF和支持向量機(support vector machine，SVM)分類器分析單頻帶和多子帶下三類MI-EEG的分類效果，在此基礎上，實現了針對個體的運動相關子帶優化組合；②基于4位受試者的MI-EEG數據集，對ICA-SF的有效性和多子帶特征組合方法進行驗證，得到了較高的分類識別率。
1 運動想象試驗
根據神經電生理知識，肢體MI能誘發運動神經皮層mu和beta節律的功率譜變化，這一現象被稱之為事件相關去同步(event-related desynchronization，ERD)和事件相關同步化(event-related synchronization，ERS)[17-19]。ERD/ERS現象是BCI系統特征分類的一個非常重要的神經電生理特征依據，其中不同類型MI誘發的ERD/ERS現象在運動皮層具有不同的神經活動分布。當想象左手或右手運動時，肢體對側運動皮層(C4/C3電極位置)的mu/beta節律出現抑制(ERD)，而同側運動皮層(C3/C4電極位置)的mu/beta節律則相對增強(ERS)。當進行想象腳運動時，中央運動皮層(Cz電極位置)的mu/beta節律出現抑制(ERD)。
為了緩解受試者疲勞度，提高采集的MI-EEG數據質量，本文的MI試驗設計在BCI Competition 2005 dataⅢa(BCI2005)試驗范式的基礎上進行了改進，具體包括以下三個方面：
(1)?減少了導聯數。其目的是減少各電極之間的干擾和增加受試者舒服度，相比BCI2005的60導聯MI-EEG數據，本文的試驗設計在保證實驗效果的同時僅采集了16導聯數據，其中主要包含運動皮層相關的導聯，如圖 1所示。
(2)?延長了休息時間。將BCI2005試驗中2 s的休息時間延長到4 s，用以緩解受試者疲勞度，提高數據質量。否則受試者會因不適而做出一些較大幅度的動作，致使電極松動或產生其他干擾信號(眼電、肌電等)。
(3)?減少了每組數據中的單次試驗(single trial)數量。本文中的每組試驗只采集了75個單次試驗。單次試驗次數的適當減少可保證MI-EEG數據的整體可靠性。
圖1
自主實驗的電極分布(上)和試驗設計(下)
Figure1.
Electrode distribution (up) and experimental design of independent experiment (down)
EEG信號采集設備為NeuroScan 40導放大器，依照國際標準10-20系統，采用16導電極，分別為VEOU，VEOL，Fp1，Fp2，FC3，FCZ，FC4，C3，CZ，C4，CP3，CPZ，CP4，O1，OZ，O2，參考電極為A1，A2。采樣頻率為250 Hz，數據精度是16 bit，采集時對原始EEG進行了50 Hz陷波濾波和0.5~100 Hz帶通濾波。受試者為實驗室4位研究生S1(女，25歲)、S2(男，29歲)、S3(女，23歲)和S4(女，22歲)，身體健康，均為右利手。實驗環境無電磁屏蔽。單次試驗持續10 s，大致分為三個時段：0~1 s的提示，1~6 s的MI和6~10 s的休息。具體過程是：每單次試驗開始，計算機發出一個短促的“滴”提示聲，提醒受試者要保持注意力集中，1 s后計算機屏幕顯示方向箭頭：向右、向下或向左的方向箭頭，提示受試者要進行右手、腳或左手的MI，每類MI持續5 s，第6秒箭頭消失，此后4 s的持續時間計算機保持黑屏，受試者則處于休息放松狀態，安靜地等待下一個單次試驗開始。其中箭頭方向出現的順序是隨機的，每類MI各25次，每組MI-EEG數據中含75個單次試驗，每次試驗采集A、B兩組MI-EEG數據。需要說明的是，綜合考慮運算效率和識別率，本文只采用了16導聯中位于運動皮層附近的9導聯MI-EEG數據(FC3、FCZ、FC4、C3、C4、CZ、CP3、CPZ和CP4)[20]。
2 方法
2.1 ICA算法
ICA是一種基于高階統計信息的盲源分離方法。該算法假定源是非高斯分布，分離后的源信號具有最大的獨立性，ICA基本原理如圖 2所示。由于EEG信號具有非高斯性，且各個局部神經皮層產生的神經活動可近似認為相互獨立，所以可利用ICA算法從頭皮EEG信號分離出具有明確神經電生理意義的神經活動[21]。
圖2
ICA算法原理圖
Figure2.
Principle chart of ICA algorithm
假定n路混合信號x(t)=[x1(t), …, xn(t)]T是由多個隱含EEG源信號s(t)=[s1(t), …, sn(t)]T的線性瞬時混合而成，如式(1)所示：
| $ \boldsymbol{x}\left(t \right)=\boldsymbol{A}\boldsymbol{s}\left(t \right) $ |
式(1)中混合矩陣A和源信號s(t)未知。與式(1)對應的分離模型如式(2)所示：
| $ \boldsymbol{u}\left(t \right)=\boldsymbol{W}\boldsymbol{x}\left(t \right)=\boldsymbol{W}\boldsymbol{A}\boldsymbol{s}\left(t \right) $ |
其中，W是分離矩陣，u(t)=[u1(t), …, un(t)]T是對源信號s(t)的近似估計，且u(t)中的各個分量間要盡可能相互獨立。本文用于多通道MI-EEG數據分析的ICA算法是實驗室自編Matlab代碼[22]，算法中采用的是信息極大獨立性判決和自然梯度優化算法[8-9]。
2.2 ICA-SF設計
在MI-BCI系統實現中，一個好的空域濾波器對提升MI-EEG信噪比和改善BCI分類效果起著非常關鍵的作用。本文采用實驗室近期提出的基于ICA空域濾波的三類運動相關獨立分量(motor-related independent components，MRICs)檢測新方法，具體實現流程如圖 3所示。
圖3
MRICs檢測濾波器流程圖
Figure3.
Flow chart of MRICs detection filter
運行ICA算法之前，對原始MI-EEG數據進行8~30 Hz數字帶通濾波，該頻段包含了運動相關的mu和beta節律成分[18]。此外，沒有進行其他偽跡處理。多通道單次MI-EEG數據經ICA方法后分解成多個獨立分量(independent component，IC)，其代表MI-EEG各個隱含源。同時得到分離矩陣W，其列向量則為ICA-SF，而混合矩陣A的列向量則表示IC在全部頭皮電極的空間模式分布(投影強度)。文獻[21]給出的基于混合矩陣A的三類ICA-SF有效性判定，可初步剔除不符合要求的單次訓練樣本，保證最終用于BCI系統實現的ICA-SF檢測性能。該判據準則的基本思路是：由訓練集的單次MI-EEG數據所得的混合矩陣A中，必須包含類似如圖 4所示三類MRICs的空間分布模式的三個列向量。如果滿足這一條件，則根據這3個列向量在A中的位置，選擇分離矩陣W中相應的列作為有效的三類MRICs檢測濾波器組：{wl(i) wr(i) wf(i), i=1, …, M≤N}，l，r和f分別代表左手、右手和腳，M為有效濾波器個數，N為訓練集中單次試驗數。在此基礎上，將所得有效濾波器組分別用于全部訓練集的單次MI-EEG進行MI標簽估計，然后根據識別率高低進一步判斷各濾波器組的性能。
圖4
三類MRICs的空間模式分布
Figure4.
Spatial distribution pattern of three-class of MIRCs
通過下面的例子可進一步說明本文采用的ICA-SF設計思想。文中的數據采集試驗每組包含75個單次試驗。根據如圖 3所示算法流程，理想情況下可得到75個基于單次MI-EEG的ICA-SF和相應的識別率。如果單次MI-EEG數據存在采集質量不高的“壞”數據，那么所設計的ICA-SF會失效或導致識別率較低，如圖 5所示。圖中存在識別率為零的情況，表明此次ICA-SF設計失敗；另外，還有一些識別率很低的情況，同樣表明對應的ICA-SF檢測性能較差，此兩種情況對應的單次MI-EEG數據不能用于設計ICA-SF。
根據已報道的MI-BCI研究論文，未見有對訓練樣本質量進行評估和選擇的討論，通常是采用全部單次訓練樣本的拼接用于ICA-SF的設計。根據本文的前期研究工作，如果不加選擇地將訓練集中一些“壞”的單次數據用于設計ICA-SF，會導致其性能下降。另外，從統計理論上看，雖然大數據樣本有利于保證ICA-SF設計性能的穩定性，但是也會使計算過程過于緩慢，不利于在線BCI的實現。因此，本文依據單次識別率大小，從合格的單次試驗數據中，選擇前L個(L的經驗值一般為10)單次MI-EEG數據拼接成用于ICA-SF設計的訓練樣本，最終得到用于MI-BCI系統實現的三類MI濾波器組：[wl wr wf]。
2.3 多子帶特征結合方法設計
多道MI-EEG經空域濾波后，通常需要選擇最優運動相關頻帶對運動相關節律進一步增強，但不同受試者的最優運動節律頻帶會因個體差異性而存在較大差別。如圖 6所示是S1和S2某組MI-EEG數據在C3、CZ和C4導聯的功率譜，其中S1在子帶10~13 Hz有明顯的譜峰，而S2在子帶10~14 Hz和19~23 Hz的譜峰都比較明顯。如果手動選擇S2的最優濾波子帶可能要經過多次嘗試才能找到，這顯然不利于實用MI-BCI系統的設計。
圖6
S1和S2與運動相關導聯對應的功率譜
Figure6.
Corresponding power spectrum of motor-related of S1 and S2
為了進一步提高MI-BCI的實用性和識別率，本文提出基于ICA-SF的多子帶特征結合方法，如圖 7所示。該方法的三類MI-EEG數據的識別步驟如下：
圖7
多子帶特征結合方法流程圖
Figure7.
Flow chart of multiple sub-band feature combination
(1)?基于如圖 3所示算法對ICA-SF設計樣本進行優化選擇；
(2)?選取識別率最高的前10個單次MI-EEG數據進行拼接，然后用來設計最終的左手、右手和腳的ICA-SF；
(3)?最優運動相關頻帶選擇：把8~30 Hz濾波頻段分為5個子帶，利用單頻帶ICA-SF方法對每個子帶進行節律增強，每個子帶下都能得到相應的分類結果[23-24]；
(4)?選取識別率高的若干子帶的空域特征進行組合，以組合后的特征向量訓練SVM分類器；
(5)?測試集輸入訓練的SVM分類器完成最后的分類識別。
3 特征提取和分類
根據如圖 7所示方法，ICA-SF從原始MI-EEG中檢測出與MI相關的MRICs成分，然后經過5個子頻帶濾波器，輸出相應的MRICs子帶成分。針對每個子帶輸出，建立相應的特征向量用作SVM分類器的訓練集和分類集。為了避免過擬合現象，采用了5×5交叉驗證。特征向量建立過程如下：
假設x為原始多通道MI-EEG數據，首先對x進行8~30 Hz帶通預處理，消除部分干擾成分；然后利用ICA算法優化設計的濾波器組[wl wr wf]提取三類MRICs成分，即：zi=wix，其中i=l, r, f；特征向量取之于zi的MI時間段，通常選擇MI提示后0.5~5 s的MI-EEG數據段[22]。在已報道的文獻中，特征向量通常為式(3)所示的標準化方差[25-27]。但是，單一方差特征通常不能完整地描述MI-EEG數據的多域特征，會導致三類MI識別率不高。因此本文采用式(4)所示的特征向量組，即直接將三個MRICs的MI時段(0.5~5 s)的樣本組合成特征向量，為了降低特征維數，可對特征樣本進行4~5倍下采樣。
| $ {F_i}=\frac{{{\mathop{\rm var}} \left({{z_i}} \right)}}{{\sum\limits_{i=l, r, f} {{\mathop{\rm var}} \left({{z_i}} \right)} }}, {F_{feature}}=[{F_l}\; {F_r}\; {F_f}] $ |
| $ {F_{feature}}=[{z_l}\; {z_r}\; {z_f}] $ |
4 實驗結果
本文數據集包含S1的4次8組數據和S2、S3、S4的各2次4組數據，共有20組MI-EEG數據集。每組有效數據點為187 500(250×10×75)。因此，S1的數據樣本總數為8×187 500，其他三位受試者S2、S3和S4的數據樣本總數均為4×187 500。為了驗證所提算法的有效性，論文基于20組MI-EEG數據集進行了兩種形式的交叉測試，即同一受試者數據集之間的自交叉測試和不同受試者數據集之間的互交叉驗證，為了表述方便，文中分別稱之為“SS_TEST”和“DS_TEST”。在用于學習分類的特征選擇方面，采用識別率最高和次高的兩個子帶特征進行組合，以驗證ICA-SF優化設計和多子帶特征結合方法的有效性。
4.1 同一受試者的自交叉測試(SS_TEST)
分別對4位受試者的20組數據進行SS_TEST，比較采用MRICs單頻帶和多子帶特征對BCI識別率的影響。以受試者S1的數據為例，給出一些代表性的識別率測試結果，最后給出全部受試者的SS_TEST結果。
受試者S1在4次MI試驗中共采集的8組MI-EEG數據集，鑒于同一受試者在同一次實驗采集中導聯位置誤差較小，因此同一次實驗采集的A、B兩組MI-EEG數據集間交叉測試識別率可能較高。為了驗證是否存在此類現象，本小節把S1的4次試驗的8組MI-EEG數據進行了分組，即分成A、B兩大組，A代表每次試驗的A組數據，B為每次試驗的B組數據。
S1-A和S1-B兩大組的SS_TEST結果如圖 8所示，橫坐標指明了各個對應的訓練集。在S1-A的SS_TEST中，訓練集S1-1-A對應的識別率描述為：在單頻帶和兩個子帶下，第1、2列是訓練集S1-1-A測試同一次采集的B組(S1-1-B)數據而得到的兩個識別率(分別對應采用單頻帶和兩子帶特征)。接下來的第3~8列是訓練集S1-1-A對隔天采集的A組其他MI-EEG數據(S1-2-A、S1-3-A、S1-4-A)的測試分類結果，類似的SS_TEST在S1_B組數據集中同樣進行。通過比較可以看出，基于兩子帶聯合特征的三類MI識別率總體上要高于基于單頻帶特征的識別率。但是，對于同一個受試者而言，無論兩組數據集是否來自同一次數據采集試驗，SS_TEST識別率并沒有體現出明顯的規律性差別。為了更直觀地比較，圖 9給出了S1的8組MI-EEG數據分別作為訓練集的交叉測試結果平均值，從中明顯可以看到，基于多子帶特征的MI識別率明顯高于基于單頻帶特征的識別率。
圖8
S1作為訓練集的SS_TEST識別率
Figure8.
Recognition rates of training set S1 by SS_TEST
圖9
訓練集S1的SS_TEST平均識別率
Figure9.
Average recognition rates of training set S1 by SS_TEST
如表 1所示，給出了對其他3位受試者數據集測的全部SS_TEST結果。由于S2~S4只含有兩次實驗的4組數據。因此，在交叉測試時，將各受試者數據分為｛S2-1，S2-2｝、｛S3-1，S3-2｝、｛S4-1，S4-2｝。然后在單頻帶和聯合子帶特征兩種情況下進行SS_TEST。實驗結果顯示，對于不同受試者的SS_TEST試驗，兩個子帶下MI識別率均比單頻帶下識別率要高。其中識別率最高提升為6.08% (S4-1)，最低為1.45%(S4-2)，平均提升了3.64%。
表1
S2、S3和S4作為訓練集的SS_TEST的平均識別率(%)
Table1.
Average recognition rates of training sets S2, S3 and S4 by SS_TEST (%)
4.2 不同受試者間的互交叉測試(DS_TEST)
為了更好地說明本文方法的有效性，此小節在不同受試者之間進行了組間交叉測試(DS_TEST)。為了表述方便，重新統一了MI-EEG數據集：S1、S2、S3和S4各2次4組試驗數據，共16組。由于受試者個體存在一定的差異性，不同受試者的最佳MI頻段也有所不同。所以在DS_TEST中，針對不同受試者使用的最優濾波頻帶也需適當調整。
如表 2所示是S1作為訓練集，S2、S3和S4作為測試集的DS_TEST分類結果，其結果表明基于聯合子帶特征的MI識別率比單頻帶下的識別率要高。如，在訓練集S1-3-A下，單頻帶下的組間測試平均識別率為77.2%，兩個子帶下的平均識別率為80.8%，相比單頻帶下的平均識別率高出3.65%；而在測試集S1-3-B下，兩個子帶下的平均識別率更是比單頻帶下的平均識別率高出5.15%。
表2
S1作為訓練集的DS_TEST識別率(%)
Table2.
Recognition rates of training set S1 by DS_TEST (%)
如圖 10所示是S2、S3和S4作為訓練集，S1作為測試集的左手、右手和腳的MI分類結果。圖中受試者的平均識別率能清晰地分析出兩種頻帶方式下分類結果的變化規律，本小節得到了與4.1小節的相同結論：三類MI在單頻帶下的平均識別率沒有兩個子帶下的平均識別率高，驗證了本文方法的可行性和有效性，為進一步的BCI系統提供了一種頻帶優化方法。
圖10
S1作為測試集的DS_TEST平均識別率
Figure10.
Average recognition rates of testing set S1 by SS_TEST
5 總結與討論
在基于空域濾波的MI-BCI研究中，由于ERD/ERS現象的空間分布特性，以及EEG信號在不同受試者身上所體現的高度復雜性與變異性，針對不同受試者的特征優化選擇對于MI-BCI系統的實現尤為重要。雖然MI-EEG數據自動特征選擇是BCI系統實現的基本任務之一，但是一些傳統方法的穩定性和適用性依然有需要改進之處。在與本文相關的研究中發現，多通道MI-EEG數據的空域濾波器技術對BCI系統整體性能起到非常關鍵的作用，雖然ICA空域濾波方法理論上能夠提取出真實的運動相關神經活動，但是ICA方法的穩定性和對隨機干擾的敏感性很大程度上限制了該方法的性能。
本文從ICA-SF的優化設計著手，盡可能保證MRICs檢測濾波器的性能，在此基礎上實現運動節律子帶的優化組合。獲得了較好的三類MI識別結果。本文僅選擇了最優和次優的兩個子帶進行驗證，一方面是出于對算法計算復雜度的考慮，另外，更多子帶的優化組合效果相比兩個子帶而言，對識別率改善并不十分明顯。原因可能是與子頻帶劃分策略有一定的關系，這也是本文后續研究的重點內容之一。
