吳林彥，魯昊，高諾，王濤

（山東建筑大學(xué)信息與電氣工程學(xué)院，山東濟(jì)南 250100）

1 引 言

腦－機(jī)接口（BCI）的概念始于1973年，它是在人或動(dòng)物腦與外部設(shè)備間建立的直接連接通路的一種新穎的人機(jī)接口方式。作為一個(gè)新興的交叉學(xué)科，BCI技術(shù)已經(jīng)引起了國(guó)內(nèi)、國(guó)際上許多科技工作者的重視，成為生物醫(yī)學(xué)、通信、計(jì)算機(jī)技術(shù)等領(lǐng)域交叉中的研究熱點(diǎn)［1］。腦電信號(hào)的特征提取對(duì)于腦電信號(hào)的處理至關(guān)重要，直接影響著腦電信號(hào)識(shí)別的準(zhǔn)確率。

現(xiàn)有的腦電信號(hào)特征提取方法主要有以下幾種：功率譜分析法，自適應(yīng)自回歸（AAR）模型法，獨(dú)立分量分析（ICA），共同空間模式（CSP），小波包分析法［2］。功率譜分析法可以反映腦電信號(hào)的能量變化，但是會(huì)損失信號(hào)的時(shí)間信息［3］。AAR模型法所需計(jì)算量較小且不需要相關(guān)頻帶的先驗(yàn)知識(shí)，但是對(duì)高度非平穩(wěn)信號(hào)效果不佳［4］。ICA算法對(duì)采集的腦電信號(hào)數(shù)據(jù)進(jìn)行獨(dú)立分量提取，一般很難得到有價(jià)值的結(jié)果［5］。CSP算法的基本原理是基于兩個(gè)協(xié)方差矩陣的同時(shí)對(duì)角化來設(shè)計(jì)最優(yōu)的空間濾波器，以分辨左手想象運(yùn)動(dòng)和右手想象運(yùn)動(dòng)腦電特征，因此不需要對(duì)被試者選擇特定的頻帶進(jìn)行特征提?。?］。小波包分析作為一種時(shí)頻分析方法，適用于非平穩(wěn)信號(hào)的處理。每級(jí)分解對(duì)低、高頻信息同時(shí)進(jìn)行分解，有利于提取更加豐富的特征。由于小波變換的變焦距特性，很容易將類別間差距最大的部分突出表示，從而將不同類之間的差異放大，這有助于提高識(shí)別正確率［7］。

由上述方法的優(yōu)劣對(duì)比，本研究選取目前較為常用且效果較好的CSP算法和小波包分析方法對(duì)競(jìng)賽數(shù)據(jù)的腦電信號(hào)進(jìn)行分析，并采用支持向量機(jī)（SVM）進(jìn)行分類，對(duì)兩種方法的正確率進(jìn)行了比較。隨后對(duì)由Emotivepoc＋系統(tǒng)采集的數(shù)據(jù)利用兩種算法進(jìn)行驗(yàn)證，并分析了本實(shí)驗(yàn)中兩種算法的特點(diǎn)。

2 算法基本理論

2.1 共同空間模式算法

CSP算法作為一種在多任務(wù)條件下，從多通道數(shù)據(jù)集中提取某一特定任務(wù)信號(hào)成分的方法［6］，步驟如下：

若將實(shí)驗(yàn)測(cè)得的腦電數(shù)據(jù)表示為N×T的矩陣E，其中N為腦電測(cè)量的通道數(shù)，T為腦電測(cè)量時(shí)每個(gè)通道的采樣點(diǎn)數(shù)，那么歸一化的腦電數(shù)據(jù)協(xié)方差矩陣可以表示為：

若用CL和CR表示想左和想右兩種情況下的空間協(xié)方差矩陣。那我們可以應(yīng)用主成分分析法，對(duì)協(xié)方差矩陣的和C按式2進(jìn)行特征分解：

其中，φ∈R3×3為特征矩陣對(duì)角陣，F(xiàn)C∈R3×3為對(duì)應(yīng)的特征向量矩陣。那么，白化矩陣為：

則協(xié)方差矩陣CL和CR可以變換為：

由此可求得廣義特征值構(gòu)成的對(duì)角陣D和特征向量矩陣B使得：

由此可以構(gòu)造如下空間濾波器：

將訓(xùn)練集的運(yùn)動(dòng)想象矩陣XL、XR經(jīng)過構(gòu)造的相應(yīng)濾波器W濾波可得特征ZL、ZR為：

根據(jù)CSP算法在多電極采集腦電信號(hào)特征提取的定義，本研究選取fL和fR為想象左和想象右的特征向量，定義如下：

對(duì)于測(cè)試數(shù)據(jù)Xi來說，其特征向量fi提取方式如式12，并與fL和fR進(jìn)行比較以確定第i次想象為想左或想右。

2.2 小波包原理

小波包分析是一種窗口面積固定而其形狀可改變，即時(shí)間窗和頻率窗均可變的時(shí)頻局部化分析方法。小波函數(shù)是小波包函數(shù)族中的一個(gè)，因此可以說小波包是小波函數(shù)的推廣［8］，其定義如下：

令正交小波基的濾波器系數(shù)分別為hn和gn，并將尺度函數(shù) φ（t）改記為 u2n（t），將小波函數(shù) ψ（t）改記為 u2n＋1（t），那么原來小波分析中關(guān)于 φ（t）和ψ（t）的雙尺度方程可以變換為：

由上式13、14所定義的函數(shù)集合｛un（t）｝，n∈Z為由u0≠φ所定義的小波包。

我們對(duì)于小波包的分析步驟可以分為小波包分解和特征提取兩部分。設(shè)原始信號(hào)為f（n），則小波包分解算法可表示為［9］：

式 15中，j＝0，1，2，…，J，j為分解層數(shù)，J為最高分解層數(shù)，Aj表示第j層上小波分解近似系數(shù)，Dj表示第 j層上小波分解逼近系數(shù)，h和g分別是正交低通分解濾波器和高通分解濾波器。表示濾波器h，g的共軛反轉(zhuǎn)；Aj*表示 Aj和的卷積，C（Aj*）表示對(duì)卷積的二元下抽樣，C為二元下抽樣算子。

本研究采用Shannon熵作為標(biāo)準(zhǔn)采用小波包熵對(duì)給定的腦電信號(hào)進(jìn)行特征提?。?0］。Shannon熵定義如下：

其中A代表信號(hào)，Aj代表信號(hào)A在一個(gè)正交小波包基上的投影系數(shù)。

3 實(shí)驗(yàn)數(shù)據(jù)分析

3.1 實(shí)驗(yàn)數(shù)據(jù)采集

本研究的數(shù)據(jù)來源分為兩部分：一部分為奧地利GRAZ大學(xué)腦機(jī)接口實(shí)驗(yàn)中心提供的運(yùn)動(dòng)想象腦電數(shù)據(jù)，該數(shù)據(jù)為2003年 BCI競(jìng)賽數(shù)據(jù)。另一部分?jǐn)?shù)據(jù)由本實(shí)驗(yàn)室采用Emotiv公司的epoc＋設(shè)備測(cè)量得到。

2003年 BCI競(jìng)賽數(shù)據(jù)實(shí)驗(yàn)的基本構(gòu)成見圖1（左）：第一階段為準(zhǔn)備期，要求受試者安靜放松坐在椅子上，無動(dòng)作。第二階段為提示階段，計(jì)算機(jī)屏幕會(huì)顯示該時(shí)間段的起始時(shí)刻，要求受試者按照提示進(jìn)行左、右手的想象動(dòng)作。第三階段是記錄階段，在此過程會(huì)完整記錄腦電變化過程。實(shí)驗(yàn)數(shù)據(jù)分為兩部分，即訓(xùn)練集和測(cè)試集，各為140次實(shí)驗(yàn)，其中左、右手想象實(shí)驗(yàn)各為70次。

本實(shí)驗(yàn)室的運(yùn)動(dòng)想象實(shí)驗(yàn)設(shè)置基本構(gòu)成見圖1（右），與競(jìng)賽數(shù)據(jù)相比，提示階段延長(zhǎng)到了2 s，記錄階段為第4 s到第9 s。每組實(shí)驗(yàn)數(shù)據(jù)共200次實(shí)驗(yàn)，其中左、右手想象實(shí)驗(yàn)各為100次。

圖1 競(jìng)賽數(shù)據(jù)實(shí)驗(yàn)時(shí)序示意圖（左）；Emotiv epoc＋系統(tǒng)測(cè)試數(shù)據(jù)實(shí)驗(yàn)時(shí)序示意圖（右）Fig 1 The scheduling of com petition data experiment（left）；the scheduling of Emotiv epoc＋system’s data experiment（right）

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

由于運(yùn)動(dòng)想象時(shí)的事件相關(guān)去同步（ERD）和事件相關(guān)同步（ERS）現(xiàn)象主要出現(xiàn)在μ（8～12 Hz）頻段，因此本實(shí)驗(yàn)采用2～24 Hz的數(shù)字濾波器對(duì)數(shù)據(jù)進(jìn)行濾波［9］。

由于小波包分析是一種時(shí)間窗和頻率窗均可變的時(shí)頻局部化分析方法。因此本研究只對(duì)要采用CSP方法分析數(shù)據(jù)進(jìn)行加窗處理。

由運(yùn)動(dòng)想象數(shù)據(jù)實(shí)驗(yàn)的時(shí)序示意圖可知，被測(cè)試人員的想象時(shí)間是從3 s開始，9 s結(jié)束。同時(shí)，由于想左想右的平均瞬時(shí)功率會(huì)隨著想象時(shí)間的延長(zhǎng)趨于一致，因此想象時(shí)間選取的時(shí)間段越短，測(cè)試的效果就越好?；谶@兩點(diǎn)，本研究對(duì)要進(jìn)行CSP特征提取的數(shù)據(jù)進(jìn)行了加窗處理，由圖1（左）可知，可選取競(jìng)賽數(shù)據(jù)每次實(shí)驗(yàn)中3.5～4.5 s的數(shù)據(jù)進(jìn)行特征提??；對(duì)于Emotiv epoc＋系統(tǒng)測(cè)得的數(shù)據(jù)，由圖1（右）可知，可選取測(cè)試數(shù)據(jù)每次實(shí)驗(yàn)中4.5～5.5 s的數(shù)據(jù)進(jìn)行特征提取。

3.3 數(shù)據(jù)特征提取及分類

按照第1節(jié)所述原理，本研究對(duì)時(shí)間加窗后的腦電信號(hào)分別進(jìn)行標(biāo)準(zhǔn)CSP特征提取和小波包變換特征提取。對(duì)于小波包變換特征提取來說，由于本研究采用的競(jìng)賽數(shù)據(jù)采樣率為1000 Hz，Emotiv epoc＋系統(tǒng)測(cè)得數(shù)據(jù)采樣率為128 Hz。因此對(duì)于競(jìng)賽數(shù)據(jù)，首先進(jìn)行降頻處理，先將其數(shù)據(jù)降為256 Hz，再進(jìn)行5層小波包變換；而對(duì)于Emotiv epoc＋系統(tǒng)測(cè)得數(shù)據(jù)，本實(shí)驗(yàn)進(jìn)行4層小波包變換。

以5層小波包變換為例，5層分解樹示意圖見圖2?？捎?jì)算出數(shù)據(jù)采樣率為256 Hz的數(shù)據(jù)變換后μ節(jié)律對(duì)應(yīng)的頻段為（5，2）節(jié)點(diǎn)。同理，數(shù)據(jù)采樣率為128 Hz的數(shù)據(jù)4層小波變化時(shí)μ節(jié)律對(duì)應(yīng)的頻段為（4，2）節(jié)點(diǎn)。所以可以根據(jù)式16計(jì)算圖2中節(jié)點(diǎn)的小波包分解系數(shù)Sij的熵：

其中，Xij表示頻帶中小波包系數(shù)Sij的時(shí)域重構(gòu)離散信號(hào)。

支持向量機(jī)（SVM）是基于統(tǒng)計(jì)學(xué)習(xí)理論的一種機(jī)器學(xué)習(xí)方法，針對(duì)二分類問題，SVM有更強(qiáng)的適應(yīng)性 、更好的分類能力以及更高的計(jì)算效率。SVM方法常用的核函數(shù)有線性核、多項(xiàng)式核以及徑向基核等，選用不同的核函數(shù)可以構(gòu)造不同的SVM分類器，本研究選用徑向基核函數(shù)：

圖2 小波包變換的5層小波樹分解示意圖Fig 25 layer wavelet tree decomposition of wavelet packet transform

其中，xi表示第i個(gè)輸入樣本，γ為核參數(shù)。對(duì)于支持向量機(jī)來說，核參數(shù)γ和誤差懲罰因子C是影響性能的主要參數(shù)。參數(shù)γ影響空間變換后的數(shù)據(jù)分布，而參數(shù)C決定了SVM的收斂速度和推廣能力。本研究采用網(wǎng)格搜索法確定參數(shù)γ和C，并采用10段交叉驗(yàn)證的方法進(jìn)行驗(yàn)證，即200次數(shù)據(jù)中，其中180次已知該次實(shí)驗(yàn)是想左還是想右，這部分?jǐn)?shù)據(jù)作為訓(xùn)練集，另外20次沒有給出運(yùn)動(dòng)方式，作為測(cè)試集［11］。

4 結(jié)果分析

根據(jù)上述的分析，本實(shí)驗(yàn)對(duì)競(jìng)賽腦電信號(hào)數(shù)據(jù)的處理過程如下：對(duì)原始腦電信號(hào)進(jìn)行2～24 Hz帶通濾波。隨后對(duì)數(shù)據(jù)分別進(jìn)行5層小波包變換特征提取和加窗后的CSP特征提取，將提取的特征利用SVM算法進(jìn)行分類，對(duì)競(jìng)賽數(shù)據(jù)樣本a和數(shù)據(jù)樣本b分別求得分類結(jié)果并計(jì)算準(zhǔn)確率，結(jié)果見表1。

由表1可知，CSP算法提取相對(duì)較好的數(shù)據(jù)時(shí)正確率能達(dá)到80%以上，而且偶然性較大。但是相對(duì)于CSP算法而言，小波包分析處理所得特征值分類效果更好也更穩(wěn)定，驗(yàn)證的最低正確率為97.5%。

表1 BCI2003競(jìng)賽數(shù)據(jù)CSP＋SVM與小波包變換＋SVM的比較Table 1 The comparison between CSP＋SVM and wavelet packet transform＋SVM in analyzing com petition data

由于競(jìng)賽數(shù)據(jù)噪聲較小，考慮到腦機(jī)接口信號(hào)采集的普適性，本研究采用同樣的方法對(duì)本實(shí)驗(yàn)室利用Emotiv epoc＋系統(tǒng)所測(cè)數(shù)據(jù)進(jìn)行分析處理。其過程為：將測(cè)量數(shù)據(jù)進(jìn)行2～24 Hz帶通濾波，隨后分別進(jìn)行CSP特征提取和以4層小波包分析為基礎(chǔ)的特征提取，并用SVM算法進(jìn)行分類，采用10段交叉驗(yàn)證的方法進(jìn)行驗(yàn)證。對(duì)腦電信號(hào)數(shù)據(jù)進(jìn)行4層小波包分析后可以得到圖3的特征向量以區(qū)分兩類運(yùn)動(dòng)想象的特征。圖3（左）為想象左右手運(yùn)動(dòng)時(shí)的小波包熵隨時(shí)間變化過程，圖3（右）為CSP特征提取時(shí)左右手的腦電信號(hào)平均功率。

圖3 小波包熵隨時(shí)間變化過程（左）；CSP特征提取時(shí)腦電信號(hào)平均功率（右）Fig 3 WPE vs.time（left）；Average power of EEG signal in CSP feature extraction（right）

在用小波包分析噪聲較大的腦電信號(hào)數(shù)據(jù)時(shí)，依然可以得到很明顯的左右手分類效果，而此時(shí)CSP的特征提取效果變得較為不理想。見圖3，CSP特征提取的腦電信號(hào)平均功率在圖上的差異不是很明顯。對(duì)Emotiv epoc＋系統(tǒng)測(cè)得的數(shù)據(jù)分類結(jié)果見表2，小波包分析結(jié)果準(zhǔn)確率均在98%以上，而CSP特征提取的準(zhǔn)確率在65%左右?？梢娦〔ò治鲈趯?duì)非純凈的腦電信號(hào)處理時(shí)也能有很好的效果。

表2 Emotiv epoc＋系統(tǒng)測(cè)試數(shù)據(jù)分類結(jié)果Table 2 Classification results of test data by Emotiv epoc＋system

5 結(jié)論

CSP算法和小波包分析方法作為近來使用率最高的兩種腦電信號(hào)處理方法，各有優(yōu)劣。本研究對(duì)關(guān)于腦電信號(hào)的CSP算法和小波包分析法進(jìn)行了比較，實(shí)驗(yàn)結(jié)果表明小波包分析法相對(duì)于CSP算法有更好的效果和普適性。經(jīng)過分析原因如下：CSP算法雖然不需要預(yù)先選擇被試者特異性的頻帶，但是需要大量的電極和對(duì)特定的腦電信號(hào)通道分析。本研究采用的競(jìng)賽數(shù)據(jù)具有59個(gè)通道，因此，利用CSP算法對(duì)競(jìng)賽數(shù)據(jù)進(jìn)行特征提取與分類時(shí)，效果尚可。但是由于Emotiv epoc＋系統(tǒng)只包含14個(gè)電極，因此對(duì)本實(shí)驗(yàn)采集的腦電信號(hào)，利用CSP算法進(jìn)行特征提取，效果并不理想。而小波包分析則不存在此類問題，因此，利用小波包分析對(duì)競(jìng)賽數(shù)據(jù)和測(cè)量數(shù)據(jù)進(jìn)行特征提取，得到的分類精確度均較高。小波包分析方法的主要缺點(diǎn)是很難獲得準(zhǔn)確的先驗(yàn)知識(shí)，在本次試驗(yàn)中由于選用10段交叉驗(yàn)證的方法，180組數(shù)據(jù)為訓(xùn)練集，20組數(shù)據(jù)為測(cè)試集，給出了較為充分的先驗(yàn)知識(shí)，因此小波包分析的效果較好。