嚴 娜,喬曉艷,李 鵬
(山西大學電子信息技術系,山西太原 030006)
視覺誘發腦電是指對視覺系統施加適當的刺激時所引起的大腦電位變化.其幅值約為0.5 μ V ~50 μ V,頻率范圍0.5 Hz~30 Hz,信噪比為0~10 dB,屬于非平穩和非高斯的微弱信號.誘發腦電反映了大腦在決策和判定過程中的認知功能,成為臨床和實驗室腦功能測試的通常選擇,也被大量應用于腦機接口(BCI)系統中作為反映大腦活動的特征信號.近年來,基于P300的BCI系統仍在不斷研究中[1,2],2008年瑞士聯邦洛桑理工大學BCI小組使用圖片輪換的方案,要求受試者默數目標刺激圖片的數目獲取視覺誘發特征信號,并成功構建腦機接口系統[1].同時,視覺誘發特征信號在測謊技術[3]、犯罪意識檢測、病態人格測試中均得到廣泛應用[4].視覺誘發腦電的特征提取方法主要包括相干平均法、功率譜分析、AAR模型法、獨立分量分析、支持向量機、高階譜分析等.疊加平均技術能有效提取誘發腦電特征,但需要多次刺激疊加,容易引起受試者神經系統疲勞,而且該方法忽略了每次刺激之間誘發電位的變異,不能反應誘發腦電的逐次變化,從而影響特征提取準確性;功率譜分析法利用了信號頻域的能量信息,但損失了時域信息;AAR模型法更適合分析平穩隨機信號,而EEG信號是高度非平穩的信號.此外,AAR模型法對偽跡很敏感;獨立分量分析能從單次(或少次)刺激中提取出誘發腦電,但各個獨立分量所蘊含的物理意義有待進一步研究;支持向量機只能提取特征,不能提取信號,因而丟失了部分信息;小波變換具有多尺度、相對帶寬恒定特點,適當選擇基本小波可在時域和頻域都具有表征信號局部特性的能力[5],但是小波分析不能有效獲取信號的相位信息,也不能抑制加性高斯噪聲.高階累積量定義的高階譜分析可以從更高階的概率結構來表征信號,完全抑制了信號中的高斯噪聲,彌補二階統計量的缺陷[6].但高階譜分析的數據來自信號時域信息,丟失了有用的頻域信息,無法進行時頻域分析.為了更有效、更精確并能從少次刺激中獲得視覺誘發腦電特征,本文采用基于小波-雙譜分析相結合的方法提取誘發腦電.
相干平均是腦電信號研究中提取腦誘發電位的傳統方法,其利用對齊疊加平均的原理,可獲得倍信噪比增益[7,8].
設各次測量信號xi(t)、誘發電位si(t)和自發腦電(可視為噪聲)ni(t)的關系表示如下



相干平均之后信噪比為

式中:σ2表示各次測量信號方差,σ2/N表示N次疊加平均誘發響應方差.可見,信噪比提高了倍.
小波變換是一種時頻分析方法,適用于非平穩信號處理,其變焦距特性,容易將特征間的差異突出表現.小波分析描述非平穩信號是在時間-尺度平面上,將信號分解成一系列小波函數的疊加,小波窗口大小隨頻率改變.在低頻段,時間分辨率較低而頻率分辨率較高;在高頻段,時間分辨率較高而頻率分辨率較低.由于小波變換的多分辨率特點,適合提取非平穩的腦電信號特征.
在多分辨分析中,尺度函數為φ(t),對應的小波函數為 ψ(t),它們滿足二尺度差分方程

式中:h0(k)對應正交低通濾波器系數;h1(k)對應正交高通濾波器系數.
由 φj,k和 ψj,k各自的正交性,h0(k),h1(k)可由下式求得

基于小波變換的分析技術主要包括小波分解和小波重構兩部分.設原始信號為f(n),如果aj(k),dj(k)是多分辨率分析中的離散逼近系數和離散細節系數,h0(k),h1(k)是滿足尺度方程的兩個正交濾波器,則小波分解算法可表示為

令j由0逐漸增大,可得到多分辨率的逐級實現[9].
小波重構算法是其分解算法的逆過程,也可由濾波器組實現,即

由于小波函數具有不唯一性,選擇不同的基本小波會產生不同的分析結果.Daubechies系列小波是工程上應用最廣泛、最成熟的緊支集正交實小波函數族,簡稱dbN小波系(其中N為小波序號這一系列),小波共有49個(db1~db49);其特點是:支集長度L=2N,消失矩階數p=N;隨著序號N的增大,時間局部性變差,但頻域局部性變好.其它類似的緊支集正交小波,例如coif小波,Sym小波特性均不如db小波,因此,通過比較分析,選擇db7小波作為腦電信號分析的小波函數.
雙譜即三階累積量譜,是具有幅值和相位的復值譜.在高階譜分析中,雙譜具有高階譜的所有特性,并且它的階次最低,計算量也是最小的,雙譜可以較好地反映信號的特征信息.
設x(n)為零均值、三階實平穩隨機序列,其三階相關函數為[10]

其雙譜表達式為

采用非參數直接雙譜估計,先計算觀測序列的FFT,再求頻域相關,算法具體步驟描述如下[7]:
設x(0),x(1),…,x(N-1)為一組測試數據,fs為采樣頻率.在雙譜域內,若 ω1和 ω2軸的頻率采樣點數為N0,則頻率抽樣間隔為Δ0=fs/N0.
1)將所給數據分成K段,每段M個觀測樣本,且每段數據之間允許重疊,對每段數據減去該段的均值,使每段成為零均值序列.
2)計算每段的DFT系數

3)根據DFT系數,分別求出每段數據的雙譜估計,即

式中:λ1,λ2分別對應 ω1,ω2軸的 DFT變換后的點,L1表示平滑點數且N0和L1應選擇為滿足的值.
4)對各段雙譜估計的結果進行統計平均

采用Donchin提供的視覺誘發Oddball腦-機接口實驗范式,該數據被BCI2005競賽選為標準數據.實驗范式構成如圖1所示,實驗中,受試者需注視此6×6字符矩陣中某個字符,目標是通過視覺刺激誘發腦電特征,識別受試者正在注視的字符,從而完成拼寫單詞的任務.
對需要拼寫的單詞中每個字母而言,字符矩陣以 2.5 s為周期顯示,其中,前400 ms字符矩陣沒有被點亮,即無顯示,接下來字符矩陣中的行和列以5.7 Hz的頻率隨機加亮,即每次加亮100 ms,加亮完成后,矩陣變黑75 ms.6×6字符矩陣共需進行12次行、列加亮,其中會有兩次包含目標字符(字符所在的特定行和特定列),這些包含目標字符的刺激引起的誘發腦電響應與不包含目標刺激的誘發腦電響應是不同的.每一個字符重復15次試驗,以保證獲得足夠強度的誘發電位.數據通過0.1 Hz~60 Hz的帶通濾波和 240 Hz頻率采樣后,以Matlab數據格式*.mat提供.有兩位受試者A和B參與了實驗,各有一個包含85個字符的訓練文件和包含100個字符的測試文件.所有數據文件以單精度存儲,每個文件內都包含64導聯的誘發腦電數據.

圖1 視覺誘發腦控拼寫器界面Fig.1 The interface of visual evoked brain control speller
由于誘發腦電信號的微弱性,在腦電采集過程中常常會存在各種干擾和噪聲,同時在誘發腦電測量過程中,不可避免地會記錄到自發腦電、眼電、肌電等多種電生理信號的偽跡以及由于基線漂移引起的記錄誤差,通常在提取特征之前,需要對采集的腦電信號進行預處理,以便抑制和消除測量信號中的干擾和噪聲,而又不損失腦電信息.由于BCI CompetionⅢ提供的數據是未經處理的原始數據,本文對腦電信號的預處理包括了導聯選擇,共平均參考,帶通濾波處理,偽跡去除,疊加平均等幾個部分.首先,對數據進行共平均參考(CAR),獲得一個相對理想的參考電極,選擇10個最優導聯,通過 0.5 Hz~30 Hz的帶通濾波器,使信號更加平滑,按照每次刺激起始點后1s進行數據分割,并將數據重排和降采樣,在以上預處理的基礎上,對數據進行15次疊加平均;然后,按照2.2中描述的多分辨率分析原理,選用db7母函數進行6層小波分解,并對第6層的細節系數進行重構濾波并白化處理,去掉相關性.最后,采用2.3中描述的雙譜分析方法分別對靶刺激和非靶刺激下重構的誘發腦電細節信號進行提取特征.
實驗數據采用BCI2000腦機接口系統和64導聯的腦電放大器記錄得到,本文僅使用Fz,Fc1,Fc2,C3,Cz,C4,Pz,PO7,PO8,FP1共10個導聯,因為在這些導聯上P300響應較為明顯.
首先,對原始標準實驗的腦電數據按照上述預處理方法進行處理.因訓練集中共有85個字符,本文以Cz導聯上的第一個字符為例,按照每次刺激起始點后1s內的數據進行處理,圖2是共平均參考后的結果,圖3顯示的是15次疊加平均的結果.經過以上預處理之后,信噪比可提高30%,該腦電數據用于后續的特征提取中.

圖2 Cz導聯上使用共平均參考進行重參考處理Fig.2 Re-reference handling used CAR on Cz channel

圖3 Cz導聯上原始腦電信號15次疊加平均Fig.3 15 times average?stack of?the original EEG on Cz channel
其次,為了便于比較,本文對記錄的視覺誘發腦電信號進行雙譜分析,同樣以Cz導聯為例,圖4和圖5分別顯士出了對第一個字符的靶刺激與非靶刺激進行雙譜分析的結果.由圖中可知不論腦功能處于何種狀態,腦電信號各分量均發生平方相位耦合現象,說明腦電是一種非常復雜的非高斯隨機過程,因此采用高階譜技術處理腦電信號,能獲得更多有關腦功能狀態的有用信息.

圖4 雙譜處理三維圖(靶刺激)Fig.4 Three-dimensional map of bispectrum processing(target stimulus)

圖5 雙譜處理三維圖(非靶刺激)Fig.5 Three-dimensional map of bispectrum processing(non-target stimulus)
最后,采用本文提出的小波-雙譜分析方法,即以db7小波母函數進行6層小波分解和重構,并對細節系數進行白化處理,最后再進行雙譜估計,本文處理的數據長度為60,將其數據分成6段,重疊率為50%,圖6和圖7分別顯示了對第一個字符的靶刺激和非靶刺激特征提取的結果.

圖6 小波-雙譜處理三維圖(靶刺激)Fig.6 Three-dimensional map of Wavelet-bispectrum processing(target stimulus)

圖7 小波-雙譜處理三維圖(非靶刺激)Fig.7 Three-dimensional map of wavelet-bispectrum processing(non-target stimuli)
由圖6,圖7可知:腦功能狀態的變化在雙譜中直接反映在腦電分量發生平方相位耦合的差異上,在大腦處于注意力高度集中時,雙譜中的三階能量幾乎集中在α波段的9 Hz的頻率分量上,且這一頻率分量的有序性迅速增強,非靶刺激的腦電誘發信號的雙譜在 α波段和θ波段均出現明顯的非線性相位耦合現象,有幾個比較明顯的譜峰,但在 9 Hz分量(μ節律)上幅值遠遠低于靶刺激在該分量的幅值.基于高階累積量的雙譜分析方法,能將腦電中反映大腦狀態的信息提取出來,較好地保留了注意相關和非注意相關時腦功能狀態的差異性.
表1給出了靶刺激和非靶刺激在 μ節律處的幅值,從表中可以明顯看出,小波-雙譜分析方法與雙譜分析方法相比,提取的腦電特征更加明顯,靶刺激相對非靶刺激在μ節律上腦電信號幅值相差近16倍.

表1 靶刺激與非靶刺激在 9 Hz頻率上腦電幅值相位耦合現象,Tab.1 The EEG amplitude of target stimulus and non-target stimulus in the frequency 9 Hz
采用現代信號處理技術分析非平穩和非線性的腦電信號并快速有效地提取特征進而對腦電信號進行自動識別與分類是腦-機接口技術研究的一個重要環節.由于通常的腦電雙譜分析僅僅利用了腦電的時域數據,丟失了頻域中大量有用的信息.因此,本文采用時頻分析中的小波變換方法對視覺誘發腦電數據進行小波分解和重構,再對重構的細節系數進行白化和雙譜分析,而且應用小波分解多分辨率特性,可以更有效地提取非平穩的腦電特征,并且仿真試驗驗證了該算法具有很好的泛化能力,可以用于腦機接口系統中.
[1]Piccione F,Giorgi F,Tonin P,et al.P300-based brain computer interface:Reliability and performance in healthy and paralysed participants[J].Clinical Neurophysiology,2006,117(3):531-537.
[2]Lebedev M A,Nicolelis M A.Brain-machine interfaces:Past,present and future[J].Trends in Neurosciences,2006,29(9):536-546.
[3]Vahid Abootalebi.A new approach for EEG feature extraction in P300-based lie detection[J].Computer Methods and Programs in Biomedicine,2008,10-14.
[4]Lebedev MA,Nicolelis MA.Brain-machine interfaces:past,present and future[J].Trends in Neurosciences,2006,29(9):536-546.
[5]徐寧壽,張建華,曹正才.小波變換在視覺誘發腦電信號提取中的應用[J].北京工業大學學報,2000,26(4):9-14.
Xu Ningchun,Zhang Jianhua,Cao Zhengcai.Application of wavelet transform technology in VEP signal extract[J].Journal of Beijing Polytechnic University,2000,26(4):9-14.(in Chinese)
[6]王群,樂建威,金頌揚,等.腦電信號高階譜分析技術的研究[J].中國醫療器械雜志,2009,33(2):79-82.
Wang Qun,Le Jianwei,Jin Songyang,et al.The study of EEG higher oder spectral analysis technology[J].Chinese Journal of Medical Instrumentation,2009,33(2):79-82.(in Chinese)
[7]黃日輝,李霆,阜艷,等.誘發腦電提取方法的研究進展[J].現代電子技術,2008,22(285):139-144.
Huang Rihui,Li Ting,Fu Yan,et al.Review of methods for extracting evoked potential[J].Modern Electronics Technique,2008,22(285):139-144.(in Chinese)
[8]王洪濤.視覺誘發電位腦機接口關鍵技術研究[J].重慶文理學院學報(自然科學版),2010,29(1):69-74.
Wang Hongtao.Researches of key techniques on brain-computer interface of visual evoked potential[J].Journal of Chongqing University of Arts and Science(Natural Science Edition),2010,29(1):69-74.(in Chinese)
[9]Eduardo Bayro-corrochano.The theory and use of the quaternion wavelet transform[J].Journal of Mathematical Imaging and Vision,2006,24:19-35.
[10]Huang Liyu,Wang Weirong,Singare S.Bispectrum quantification analysis of EEG and artificial neural network may classify ischemic states[J].Lecture Notes in Computer Science,2006,4233:533-542.