李 芳 邱天爽 馬 征

（大連理工大學(xué)電子與信息工程學(xué)院，大連 116024）

引言

腦機(jī)接口（BCI）是指建立一種不依賴于人體外周神經(jīng)和肌肉組織而在人與周圍環(huán)境間進(jìn)行信息交流與控制的新型通道。近年來，隨著世界人口的增長和老齡化的加劇，帕金森病、癲癇、脊髓嚴(yán)重?fù)p傷及糖尿病引起的外周神經(jīng)疾病患者的數(shù)量也以驚人的速度增加。這些疾病給病人造成了嚴(yán)重的感覺運(yùn)動(dòng)功能損傷，病人如同閉鎖在自己的世界里，無法同外界進(jìn)行交流。腦機(jī)接口為這些患者找到了一種可能的解決方案。為了提高腦電信號(hào)分類正確率，需要改善腦電信號(hào)的特征提取和機(jī)器學(xué)習(xí)方法。

在特征提取方面，由于腦電信號(hào)是一種典型的非平穩(wěn)信號(hào)，因此在方法的選擇上受到了一定的局限，而有“數(shù)字顯微鏡”之稱的小波變換（wavelet transform，WT）是適合處理非平穩(wěn)信號(hào)的一種方法［1］；目前在BCI中所使用的EEG信號(hào)，如SCP（皮層慢電位）或P300（事件相關(guān)電位）等往往不是相互獨(dú)立的，從腦科學(xué)的角度來看，相關(guān)信息的結(jié)合將會(huì)提高腦電信號(hào)的識(shí)別性能，這意味著時(shí)頻特征的結(jié)合會(huì)更有利，而小波分析正是非常有效的時(shí)頻分析方法。在實(shí)現(xiàn)基于小波分解系數(shù)及子帶系數(shù)均值特征算法的基礎(chǔ)上［2］，發(fā)現(xiàn)采用更多的信息可以得到更好的分類性能，在此啟發(fā)下，本研究嘗試在特征組合中加入小波熵這一普適特征量［3-4］。結(jié)果表明，小波系數(shù)、子帶系數(shù)均值和小波熵的組合特征量能更準(zhǔn)確地表示含有瞬變成分的信號(hào)的特征信息，是腦電信號(hào)分析和處理的有效途經(jīng)。

在特征分類方面，BCI中使用的分類方法主要有人工神經(jīng)網(wǎng)絡(luò)［5］、支持向量機(jī)、K均值聚類、遺傳算法等等。其中常用的BP網(wǎng)絡(luò)雖然實(shí)現(xiàn)方便，但訓(xùn)練過程很費(fèi)時(shí)，因此針對(duì)少量的有標(biāo)注樣本和大量的無標(biāo)注樣本，半監(jiān)督方法是最適合的?，F(xiàn)有的半監(jiān)督算法包括：自訓(xùn)練、合作訓(xùn)練、期望最大化算法和直推式向量機(jī)等；在經(jīng)過一系列改進(jìn)之后，目前比較流行的應(yīng)用于腦電信號(hào)分類的方法是基于自訓(xùn)練半監(jiān)督的支持向量機(jī)［6］，本研究針對(duì)皮層慢電位數(shù)據(jù)對(duì)上述方法進(jìn)行改進(jìn)，依據(jù)此原理還演繹出基于自訓(xùn)練半監(jiān)督的K均值聚類法，實(shí)驗(yàn)數(shù)據(jù)的分析結(jié)果證明了其快速的收斂性和分類的準(zhǔn)確性。

1 基于小波系數(shù)及子帶系數(shù)均值與小波熵結(jié)合的特征提取方法

1.1 特征提取算法原理

1.1.1 腦電信號(hào)的小波分解

在選定了小波函數(shù)和分解層數(shù)之后即可進(jìn)行小波分解。如以Daubechies類db4小波對(duì)信號(hào)進(jìn)行6層分解為例，輸出的分解結(jié)構(gòu)包含小波分解系數(shù)向量C和系數(shù)長度向量L。系數(shù)向量C結(jié)構(gòu)為：C=［cA6，cD6，cD5，cD4，cD3，cD2，cD1］，其中cA6為逼近系數(shù)，cD6，cD5，…，cD1為小波系數(shù)。

1.1.2 基于分解系數(shù)的特征

以Fisher距離指標(biāo)作為特征的可分離性度量［2］，式中Sb為類間離散度矩陣，Sw為類內(nèi)離散度矩陣，J值越大代表該特征的類別可分離性越好。各子帶可以分別按不同的比例，從中選取J值大的逼近系數(shù)、小波系數(shù)。具體比例值應(yīng)依據(jù)人不同的思維和實(shí)驗(yàn)范例來確定。在小波分解C結(jié)構(gòu)中，通過關(guān)鍵系數(shù)選擇得到C′=［c1，c2，…，cn］，即特征系數(shù)向量，n為系數(shù)向量的維數(shù)。

1.1.3 基于小波熵的特征

一個(gè)隨機(jī)信號(hào)，如果是由完全無序的過程生成的，那么其在各頻段的幅度和能量近似相同，即信號(hào)的概率分布接近無序，熵值理論上接近最大值。如果用小波對(duì)該信號(hào)進(jìn)行多尺度分析，當(dāng)該隨機(jī)信號(hào)的不同尺度的系數(shù)矢量越接近時(shí)，信號(hào)的熵值越大，信號(hào)信息量越小，相反則信息量越大。因此將熵值作為衡量信號(hào)信息量大小的特征值。

將尺度s上的小波系數(shù)矢量記為：Ws=（ws1，ws2，…，wsn），其中下標(biāo)s表示尺度參數(shù)，矢量元素wsi（i=1，2，…，n）為小波系數(shù)。也可以用Rn空間中矢量的范數(shù)來衡量系數(shù)矢量的接近程度，即：

根據(jù)小波能量的定義，用小波系數(shù)矢量的范數(shù)序列來衡量信號(hào)在各尺度的接近程度，即用各尺度的能量來衡量各頻段的信號(hào)成分的接近程度。這里對(duì)能量序列E1，E2，…，EK歸一化，通過歸一化的能量序列的分布來分析信號(hào)的結(jié)構(gòu)和復(fù)雜程度，過程如下：

假設(shè)對(duì)信號(hào)在K個(gè)尺度上進(jìn)行分解，令尺度上的小波系數(shù)矢量為Wj=（wj1，wj2，…，wjn），尺度j的能量定義為

能量序列的分布定義為各尺度的歸一化能量

1.1.4 基于子帶系數(shù)均值的特征

EEG經(jīng)多尺度小波分解后，各小波子帶系數(shù)均值分別代表不同頻率段的時(shí)間均值信息，與信號(hào)整體均值相比，更能細(xì)致地刻畫均值信息。小波子帶系數(shù)均值為

式中，ni為子帶系數(shù)長度，i為子帶序號(hào)。求出每個(gè)子帶的系數(shù)均值，并基于Fisher距離判別準(zhǔn)則，按一定比例，從中選取部分系數(shù)均值作為特征，得到M′=｛m1，m2，…，mq｝特征均值向量，q為特征向量的維數(shù)。

1.1.5 3種特征量的形成與結(jié)合

1）對(duì)每個(gè)訓(xùn)練樣本，用指定的小波函數(shù)進(jìn)行指定層次的分解，得到系數(shù)向量C。

2）各子帶按不同的比例，從C中選取部分具有較高可分離能力的逼近系數(shù)、小波系數(shù)，得到C′，可分離能力指標(biāo)基于上述的Fisher準(zhǔn)則。基于同樣的過程可以從小波各子帶系數(shù)均值中選擇部分系數(shù)均值得到M′。

3）對(duì)每個(gè)訓(xùn)練樣本，用指定的小波函數(shù)進(jìn)行不同尺度（j=1，2，…，k）的分解。用每次分解后得到的小波系數(shù)來計(jì)算該尺度的能量，最終進(jìn)行歸一化，對(duì)歸一化的能量序列進(jìn)行小波熵的計(jì)算，得到小波熵的特征序列N′。

4）按照上述方法，得到C′、M′和N′，假設(shè)EEG信號(hào)的采樣通道為CH1～CHL，則結(jié)合所有待分析的EEG通道的C′、M′和N′形成特征向量，L為所分析的EEG通道數(shù)。將特征向量其歸一化，形成最終的特征向量：X=｛C′1，M′1，N′1，C′2，M′2，N′2，…，C′L，M′L，N′L｝。

2 基于自訓(xùn)練半監(jiān)督的支持向量機(jī)和K均值聚類分類算法

2.1 改進(jìn)的基于自訓(xùn)練半監(jiān)督的支持向量機(jī)

2.1.1 算法原理

標(biāo)準(zhǔn)SVM分類器定義如下（針對(duì)兩分類問題）［6］：

yi（wTxi+b）≥1-ξi，ξi≥0，i=1，2，…，N（6）式中，xi∈Rn是一個(gè)訓(xùn)練樣本（特征向量），yi∈｛-1，1｝是xi的標(biāo)注，C＞0是歸一化常量。

自訓(xùn)練半監(jiān)督算法即是以少量訓(xùn)練樣本（有標(biāo)注）訓(xùn)練SVM，并以循環(huán)的方式不斷更新訓(xùn)練樣本（添加無標(biāo)注的樣本）及SVM，最終得到一個(gè)優(yōu)化后的支持向量機(jī)，并對(duì)測(cè)試樣本進(jìn)行最終分類。該方法的改進(jìn)之處在于每次循環(huán)更新訓(xùn)練樣本時(shí)，以等量的有標(biāo)注訓(xùn)練樣本加入到新建的訓(xùn)練樣本中，這樣充分利用了已有的訓(xùn)練樣本，同時(shí)也減少了時(shí)間的消耗。

2.1.2 算法實(shí)現(xiàn)的具體步驟

已知條件：一個(gè)小的訓(xùn)練數(shù)據(jù)集FI，含有N1個(gè)樣本｛xi，i=1，…，N1｝，其標(biāo)注為［y0（1），…，y0（N1）］；一個(gè)大的測(cè)試樣本集FT，含有N2個(gè)樣本｛xN1+i，i=1，…，N2｝。

步驟1：用FI的少量樣本訓(xùn)練一個(gè)SVM，并對(duì)FT進(jìn)行分類。此時(shí)SVM的參數(shù)為w（1）∈Rn，ξ（1）∈RN1，b（1），得到的預(yù)測(cè)標(biāo)注為［y（1）（1），…，y（1）（N2）］；

步驟2：k次循環(huán)（k≥2）遵循步驟2.1～2.3

步驟2.1：定義一個(gè)新的訓(xùn)練樣本集FN=FI+FT+（k-1）×△F′I，F(xiàn)T的標(biāo)注由第k-1次循環(huán)預(yù)測(cè)得到，△F′I為每次循環(huán)增加的訓(xùn)練集中已標(biāo)注的樣本增量；

步驟2.2：用FN訓(xùn)練一個(gè)SVM，并再次對(duì)FT進(jìn)行分類。這時(shí)的SVM參數(shù)表示為w（k）∈Rn，ξ（k）∈RN1+N2，b（k），預(yù)測(cè)標(biāo)注更新為［y（k）（1），…，y（k）（（N2）］；

步驟2.3：計(jì)算目標(biāo)函數(shù)值

步驟3：（循環(huán)終止步驟）假定給一個(gè)正常量δ0，如果滿足下式

或者已標(biāo)注的訓(xùn)練樣本全部訓(xùn)練過，則算法在第k次循環(huán)后停止，測(cè)試樣本的最終分類標(biāo)注為［y（k）（1），…，y（k）（N2）］；否則進(jìn)行第k+1次循環(huán)。

2.2 改進(jìn)的基于自訓(xùn)練半監(jiān)督的K均值聚類

2.2.1 改進(jìn)算法的具體步驟

在了解K均值聚類算法［8］的基礎(chǔ)上，確定訓(xùn)練樣本的中心，并在自訓(xùn)練半監(jiān)督的過程中不斷修改中心值，直到達(dá)到收斂為止。算法流程如下：

1）將兩類訓(xùn)練樣本分別求均值構(gòu)成初始聚類中心；

2）取出一小部分訓(xùn)練樣本（該實(shí)驗(yàn)取了前L個(gè)，L即一個(gè)通道內(nèi)的特征量維數(shù)），聚類中心為2（根據(jù)類別數(shù)決定），K均值算法內(nèi)最大循環(huán)次數(shù)為Max=1 000，容錯(cuò)誤差在10-8左右；

3）初始的訓(xùn)練及測(cè)試：用L個(gè)訓(xùn)練樣本構(gòu)造均值聚類中心，獲得新的聚類中心；以新的聚類中心來對(duì)測(cè)試樣本進(jìn)行劃分T0、T1；

4）循環(huán)：將上次劃分的測(cè)試樣本T0、T1分別加入到兩類訓(xùn)練樣本中，并重新確定聚類中心，若達(dá)到最大循環(huán)次數(shù)或滿足容錯(cuò)誤差，則以此中心重新對(duì)測(cè)試樣本進(jìn)行分類；

5）判斷算法是否達(dá)到收斂（收斂的標(biāo)準(zhǔn)為新的聚類中心與兩類訓(xùn)練樣本的原始聚類中心重合），若收斂，則停止循環(huán)并計(jì)算最終分類結(jié)果；否則繼續(xù)4）。

3 實(shí)驗(yàn)結(jié)果及討論

3.1 P300信號(hào)的特征提取處理結(jié)果

該實(shí)驗(yàn)主要是驗(yàn)證改進(jìn)的特征提取算法的有效性?；贠ddball模式，設(shè)計(jì)了靶刺激與非靶刺激的實(shí)驗(yàn)來產(chǎn)生P300信號(hào)，并在本實(shí)驗(yàn)室用Neuroscan采集了大量的數(shù)據(jù)，之后分別對(duì)3個(gè)被試的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行預(yù)處理來獲得4種組合方式的訓(xùn)練樣本和測(cè)試樣本。為了驗(yàn)證改進(jìn)的特征提取算法的優(yōu)越性，分別應(yīng)用文獻(xiàn)［2］中的小波系數(shù)及子帶均值法和所提出的結(jié)合小波熵的特征提取算法來獲得最終的特征量，輸入到BP神經(jīng)網(wǎng)絡(luò)得到測(cè)試結(jié)果［5］。

從表1中可以看出，利用改進(jìn)的3種特征量結(jié)合的算法，最優(yōu)分類正確率可達(dá)到95.90％，而且相比于小波分解系數(shù)和子帶系數(shù)均值法，分類效果基本上都有所提高，最多提高2.27％。這說明小波熵的加入，補(bǔ)充了信號(hào)信息量大小這一特征，豐富了單一的不同頻段上的時(shí)域特征，從而提高了信號(hào)的識(shí)別率。為了進(jìn)一步提高分類正確率，還利用了改進(jìn)的自訓(xùn)練半監(jiān)督分類器，下文主要利用BCI2003競賽數(shù)據(jù)來驗(yàn)證算法的有效性。

表1 實(shí)驗(yàn)數(shù)據(jù)的分析結(jié)果以及兩種特征提取算法的比較Tab.1 Data analysis results and contrast of two feature extraction methods

3.2 BCI2003競賽數(shù)據(jù)處理結(jié)果

3.2.1 數(shù)據(jù)描述

數(shù)據(jù)取自一個(gè)健康的被試，為BCI2003競賽Data Set Ia的標(biāo)準(zhǔn)數(shù)據(jù)。被試按要求對(duì)計(jì)算機(jī)屏幕上的一個(gè)光標(biāo)進(jìn)行向上或向下的移動(dòng)，同時(shí)他的皮層電位被記錄下來。在記錄的同時(shí)將其皮層慢電位（Cz-乳突）通過視覺的方式反饋給被試［9］。正皮層電位將會(huì)使屏幕光標(biāo)向下移動(dòng)，而負(fù)電位將使光標(biāo)向上移動(dòng)。

刺激間隔為6s。對(duì)于每次試驗(yàn)，從開始后0.5s到試驗(yàn)結(jié)束，將會(huì)在屏幕的上部或下部高亮顯示一個(gè)目標(biāo)，被試的任務(wù)就是根據(jù)這一目標(biāo)所表明的正向或負(fù)向含義進(jìn)行相應(yīng)的活動(dòng)。視覺反饋的呈現(xiàn)時(shí)間為第2s到第5.5s。本數(shù)據(jù)集將每次試驗(yàn)僅取這3.5s的數(shù)據(jù)用于訓(xùn)練和檢驗(yàn)。采樣頻率為256Hz，因此每個(gè)導(dǎo)聯(lián)將得到896個(gè)樣本點(diǎn)。

EEG記錄來自6個(gè)位置：A1，A2，C3f，C3p，C4f，C4p，參考電極為Cz。共有268個(gè)訓(xùn)練樣本，293個(gè)檢驗(yàn)樣本。訓(xùn)練樣本中前135個(gè)為光標(biāo)向下（類別“0”），后133個(gè)為光標(biāo)向上（類別“1”）。數(shù)據(jù)分析的目標(biāo)是依據(jù)這些原始數(shù)據(jù)將檢驗(yàn)樣本的類別“0”與“1”正確區(qū)分。

3.2.2 特征提取處理結(jié)果

以小波變換分解系數(shù)、子帶系數(shù)均值和小波熵作為初始特征向量，以Fisher距離指標(biāo)作為特征的可分離性度量，針對(duì)所采用的數(shù)據(jù)背景，各子帶系數(shù)的選取比例設(shè)為：30％、35％、23.5％、12.9％、5％、2％、0.22％，小波分解的不同尺度為j=1，2，…，6，最終獲得特征向量。圖1表示通道1獲取的特征量分布情況，每個(gè)通道共43維，系數(shù)分布取36維，系數(shù)均值分布取1維，小波熵為6維。為了更清晰地顯示兩類樣本的分類效果，如圖2所示。

3.2.3 分類算法處理結(jié)果及分析

1）改進(jìn)的基于自訓(xùn)練半監(jiān)督的支持向量機(jī)

文獻(xiàn)［6］中自訓(xùn)練半監(jiān)督SVM算法中的初始訓(xùn)練樣本是固定的，這是由于訓(xùn)練樣本數(shù)較少的關(guān)系，但針對(duì)本實(shí)驗(yàn)數(shù)據(jù)，為了充分利用已知標(biāo)注的訓(xùn)練樣本，算法中將不斷更新初始訓(xùn)練樣本。

圖1 訓(xùn)練樣本的通道一特征量分布Fig.1 The feature distribution of train data

圖2 訓(xùn)練樣本的通道一特征量分布Fig.2 Two-dimension feature distribution of train data

更新方法：初始樣本量為50（已標(biāo)注的），然后以20為遞增量，將上一輪的測(cè)試樣本的最終標(biāo)注及相應(yīng)的樣本加入初始樣本中訓(xùn)練SVM，更新次數(shù)根據(jù)測(cè)試樣本量的大小而定，遞增量可變。

分別以不同的特征量來驗(yàn)證該分類器的有效性，結(jié)果如下：

特征1：只以小波分解系數(shù)和子帶系數(shù)均值為特征量，每個(gè)通道37維。由于特征維數(shù)的大小會(huì)影響SVM的分類效果，所以對(duì)特征量進(jìn)行了均值化降維處理，得到的最終樣本維數(shù)為12。

特征2：以小波變換分解系數(shù)、子帶系數(shù)均值和小波熵作為特征向量，每個(gè)通道43維。降維處理后得到的特征量維數(shù)是18。

如圖3所示，對(duì)于特征一的分類結(jié)果，最高分類正確率為88.40％，而基于改進(jìn)的特征量的分類結(jié)果有了顯著的提高，最高達(dá)到91.13％。經(jīng)過計(jì)時(shí)，它們的收斂時(shí)間相差約5s，而且從曲線的趨勢(shì)來看，特征二能夠很快地達(dá)到一定的分類精度，有利于今后實(shí)現(xiàn)實(shí)時(shí)在線的工作。這里，特征量經(jīng)過降維處理后達(dá)到很好的分類效果，說明對(duì)于該分類器，過多的原始特征向量反而造成分類時(shí)的干擾。

圖3 改進(jìn)的基于半監(jiān)督的SVM的分類正確率曲線。（a）針對(duì)特征1；（b）針對(duì)特征2Fig.3 The classification accuracy using improved semi-supervised SVM classifier.（a）feature 1；（b）feature 2

2）基于自訓(xùn)練半監(jiān)督的K均值聚類

對(duì)由小波系數(shù)、子帶系數(shù)均值和小波熵組成的258維特征量進(jìn)行模式識(shí)別，分別采用K均值聚類法和改進(jìn)的基于自訓(xùn)練半監(jiān)督的K均值聚類法來分類測(cè)試樣本，實(shí)驗(yàn)結(jié)果如下：

①單純的K均值聚類：如果所有特征量都加入，分類正確率為87.03％；如果只加入43維的特征量（即一個(gè)通道的特征量），則分類正確率為72.01％

②K均值+半監(jiān)督：對(duì)于所有訓(xùn)練樣本，在加入43維特征量的情況下，最優(yōu)的結(jié)果是90.10％

分類器②比較穩(wěn)定，而且花費(fèi)的時(shí)間很短，達(dá)到的分類正確率也比較理想。

3）算法比較

分別用BP神經(jīng)網(wǎng)絡(luò)、基于半監(jiān)督的支持向量機(jī)和基于自訓(xùn)練半監(jiān)督的K均值聚類算法對(duì)特征一和特征二進(jìn)行訓(xùn)練和識(shí)別，結(jié)果見表2。

顯然，本研究提出的兩種基于自訓(xùn)練半監(jiān)督的改進(jìn)方法，相比于BP網(wǎng)絡(luò)，大大提高了收斂速度，有效地減少了訓(xùn)練時(shí)間，并且比較可靠穩(wěn)定。其中，改進(jìn)的支持向量機(jī)對(duì)特征量的維數(shù)要求較高（通常需要降維），但同時(shí)能很好的分辨特征量的不同組合；而改進(jìn)的K均值聚類法則對(duì)特征量的變化不太敏感，但需要大量的特征量來支持。

表2 不同特征量組合下三種分類方法的分類精度的比較Tab.2 The accuracy rates of three classification methods in different features

4 結(jié)論

實(shí)驗(yàn)結(jié)果表明，本研究提出的新的小波特征量組合法以及改進(jìn)的基于自訓(xùn)練半監(jiān)督SVM和K均值聚類法，能夠有效地從自發(fā)腦電中提取出皮層慢電位，使被試實(shí)現(xiàn)對(duì)光標(biāo)控制，不僅獲得較為理想的分類正確率，而且在訓(xùn)練時(shí)間和收斂速度上有所貢獻(xiàn)，為在實(shí)際系統(tǒng)中應(yīng)用奠定了基礎(chǔ)，為腦機(jī)接口的研究提供了一種更加實(shí)用、更加自然的控制方式。

［1］Qi Limei，Li Xiaofeng，Zhang Guozhu.Wavelet transform theory and its application in signal processing［J］.Journal of the University of Electronic Science and Technology of China，2008，37（3）：386-388.

［2］楊幫華.自發(fā)腦電腦機(jī)接口技術(shù)及腦電信號(hào)識(shí)別方法研究［D］.上海：上海交通大學(xué)，2007.

［3］He Zhengyou，Liu Zhigang，Qian Qingquan.Study on wavelet entropy theory and adaptability of its application in power system［J］.Power System Technology，2004，28（21）：17-21.

［4］He Zhengyou，Chen Xiaoqing，Luo Guoming.Wavelet entropy measure definition and its application for transmission line fault detection and identification（part I：Definition and methodology）［A］.In：International Conference on Power System Technology 2006［C］.Chongqing：IEEE，2006.1-6.

［5］周開利，康耀紅.神經(jīng)網(wǎng)絡(luò)模型及其MATLAB仿真程序設(shè)計(jì)［M］.北京：清華大學(xué)出版社，2005.

［6］Li Yuanqing，Li Huiqi，Guan Cuntai，et al.A self-training semisupervised support vector machine and its application in brain computer interface［A］.In：IEEE International Conference on Acoustics，Speech and Signal Processing［C］.Honolulu：IEEE，2007.385-388.

［7］李建勛，柯熙政，郭華.小波方差與小波熵在信號(hào)特征提取中的應(yīng)用［J］.西安理工大學(xué)學(xué)報(bào)，2007，23（4）：365-368.

［8］徐義峰，陳春明，徐云青.一種改進(jìn)的K-均值聚類算法［J］.計(jì)算機(jī)應(yīng)用與軟件.2008，25（3）：.275-277.

［9］Hinterberger T，Mellinger J，Birhaumer N，et al.The thought translation device：structure of a multimodal brain-computer communication system［A］.In：Proceedings of the 1st international IEEE EMBS Conference on Neural Engineering［C］.Capri Island：IEEE，2003.603-606.