胡春海 信思旭 劉永紅

(燕山大學河北省測試計量技術及儀器重點實驗室 秦皇島 066004)

?

單導聯(lián)腦電信號P300少試次提取算法研究①

胡春海②信思旭 劉永紅

(燕山大學河北省測試計量技術及儀器重點實驗室 秦皇島 066004)

研究了腦電(EEG)信號中的正向電位P300的提取，為了克服P300信號提取中對電極數(shù)量和試驗次數(shù)的苛刻要求，提出了一種單導聯(lián)腦電信號P300少試次提取新方法，該方法首先對原始信號進行小波變換去噪處理，利用分離盲源的特征矩陣聯(lián)合近似對角化(JADE)算法對得到的觀測信號進行分解；其次為避免少試次提取中可能存在的突發(fā)性錯誤，依據(jù)皮爾遜相關系數(shù)(PCC)和列向標準差構(gòu)建層次分析模型，對P300分量進行最優(yōu)提取；最后利用時域能量熵的統(tǒng)計分布特性和小波變換原理對重構(gòu)后的P300信號進行時頻域的補償。試驗結(jié)果表明，該方法只需對單導聯(lián)3試次的數(shù)據(jù)進行處理，并能夠有效地提取P300成分。

腦電信號處理， P300， 時域能量熵， 小波變換， 層次分析模型

0 引 言

腦-機接口(brain computer interface, BCI)是在人腦和計算機或其他電子設備之間建立的不依賴于常規(guī)大腦信息輸出通路的對外信息交流和控制的通信系統(tǒng)。作為事件相關電位(event related potential, ERP)中重要的內(nèi)源性成分，P300腦電信號(EEG倍號)被廣泛應用于BCI系統(tǒng)中[1-3]。P300是受試者受到發(fā)生概率較小的外界事件刺激時在頭頂中縫區(qū)域內(nèi)誘發(fā)的潛伏期為300ms左右的正向電位，不需要初始訓練，刺激物與被試者的利害關系及被試者的情緒都會在P300上有所反映。然而P300成分十分微弱，實際采集到的EEG信號伴隨著大量的神經(jīng)源性噪聲和非神經(jīng)源性干擾，因此如何從EEG數(shù)據(jù)中提取出P300成分對認知生物科學和電子與信息科學等研究具有十分重要的意義[4]。

目前針對P300特征提取技術有很多方法，包括相干平均、參數(shù)模型、自適應濾波、小波分析[5]、神經(jīng)網(wǎng)絡[6]以及主成分分析[7, 8]等，均能夠較好的提取P300成分，并達到可觀的分類效果，但是其大量的導聯(lián)數(shù)目增加了BCI系統(tǒng)操作的繁瑣程度，冗雜的EEG數(shù)據(jù)也會嚴重加大算法的計算負荷，任一導聯(lián)數(shù)據(jù)采集的好壞均會對BCI系統(tǒng)的識別結(jié)果產(chǎn)生影響；另外，由于EEG信號低信噪比和非平穩(wěn)隨機的特點，即使對于同一被試，誘發(fā)的P300也可能存在一定的區(qū)別，多次平均忽略了試驗之間的差異性，多次刺激容易引起被試者的神經(jīng)性疲勞，從而影響波形，這些均很難滿足BCI系統(tǒng)在線應用的要求。針對上述問題，本文提出一種單導聯(lián)腦電信號P300少試次提取的新方法，將小波變換與盲源分離(blind source separation, BSS)相結(jié)合，利用層次分析模型提取源信號估計分量中的P300成分，并結(jié)合時域能量熵的統(tǒng)計分布特征對重構(gòu)后的信號進行時域和頻域的補償。

1 試驗描述及數(shù)據(jù)來源

試驗所用數(shù)據(jù)來自于BCI Competition II中的Data set IIb，該數(shù)據(jù)由基于經(jīng)典Oddball范式的拼寫器試驗產(chǎn)生，使用標準化的64導聯(lián)電極帽采集，并通過BCI2000系統(tǒng)記錄，記錄部位按照國際10-20導聯(lián)系統(tǒng)分布。試驗中，被試者面對一個6×6字符矩陣，并集中注意力觀察要拼寫的字符。對于每一個待選的字符，矩陣的行和列都會隨機高亮12次，其中必定有一行和一列的高亮包含了所需要的字符，這樣即構(gòu)成了概率分別為0.17和0.83的偏差刺激和標準刺激。被試者對偏差刺激所在行和列的高亮次數(shù)進行計數(shù)，便構(gòu)成了具有隨機性和突發(fā)性的靶刺激，該刺激誘發(fā)的腦電信號中便含有強度較大的P300成分。

P300的頭皮分布相對集中在中線部位，其波幅在頂后部Pz處和中央頂部Cz處較大。大量試驗表明中央頂部記錄的信號效果更為明顯，因此只選擇電極Cz處的EEG數(shù)據(jù)進行分析處理。對單導聯(lián)采集到的數(shù)據(jù)分段整理，試驗共拼寫42個字符，每個字符進行15次試驗，共得到1260次靶刺激和6300次非靶刺激數(shù)據(jù)段，每個數(shù)據(jù)段以刺激時刻0開始，到1000ms結(jié)束，因采樣頻率是240Hz，則每一次刺激的數(shù)據(jù)長度為240。

2 單導聯(lián)腦電信號P300少試次提取算法

提取單導聯(lián)3試次腦電信號中P300成分的過程如圖1所示。具體提取步驟如下：

步驟1：對3試次的原始信號分別進行小波變換處理，濾掉0～7.5Hz之外的噪聲干擾。

步驟2：利用盲源分離的特征矩陣聯(lián)合近似對角化(joint approximate diagonalization of eigenmatrices, JADE)算法對觀測信號進行分解，得到3個源信號的估計分量。

步驟3：依據(jù)相關系數(shù)和列項標準差構(gòu)建層次分析模型，對P300分量進行最優(yōu)提取，并結(jié)合混合矩陣重構(gòu)P300信號。

步驟4：依據(jù)時域能量熵的統(tǒng)計分布特性對重構(gòu)后的P300信號進行熵值調(diào)整，并利用小波變換作平滑處理，進行時頻域的信號補償。

圖1 單導聯(lián)腦電信號P300少試次提取算法流程

2.1 基于小波變換的頻域濾波

P300的頻譜主要集中分布在低頻區(qū)域，尤其是2～8Hz之間的EEG信號對P300的提取結(jié)果影響較大，為了濾除其頻率之外的噪聲干擾，選取小波變換[9]作為P300腦電信號的頻域濾波方法。具體步驟為：

(1) 母小波采用Daubechies4小波，對每個EEG信號進行4層小波分解。

(2) 將分解得到的不同尺度下的高頻細節(jié)分量置零處理。

(3) 進行多層小波重構(gòu)，得到頻帶為0～7.5Hz的小波去噪后的腦電信號。

如圖2、圖3所示，隨機選取3組電極Cz處的靶刺激信號作為原始EEG信號，對其進行上述小波變換處理，得到頻域濾波后的3組觀測信號。

圖2 原始EEG信號

圖3 頻域濾波后的EEG信號

P300腦電信號與外界刺激之間具有良好的鎖相關系，其波形不會因重復試驗而發(fā)生很大改變，通?？刹捎?0～100次累加平均來提高P300成分的信噪比。為便于比較分析，隨機選取40組經(jīng)過頻域濾波的靶刺激信號進行累加平均，將其作為標準信號，如圖4所示。將上述3組經(jīng)過頻域濾波的靶刺激信號進行累加平均，將其作為參考信號，如圖5所示。

圖4 標準信號-40試次累加平均

圖5 參考信號-3試次累加平均

2.2 基于盲源分離JADE算法的特征提取

經(jīng)過小波變換子帶重組的多通道輸出信號中，P300成分得到顯著增強，非目標信號得到有效削弱，可以認為小波變換后的觀測信號由3個信源組成，分別為P300、自發(fā)EEG以及噪聲干擾，因此應用JADE算法對3試次信號進行盲源分離，能夠有效估計出相對較強的P300成分[10]。對前述經(jīng)過頻域濾波后的3組靶刺激信號進行JADE盲源分離，得到3個源信號的估計分量，可以通過如下兩種方法來提取其中的P300分量：

(1) 相關系數(shù)最大化

樣本變量X、Y間的皮爾遜相關系數(shù)(Pearsom correlation coefficient, PCC)定義為

(1)

將得到的3個源信號的估計分量分別和靶刺激標準信號作相關分析，相關系數(shù)最大的即為P300分量。但是，由于觀測信號均為單導聯(lián)單次試，對非靶刺激來說，任一個觀測信號含有大量的P300成分，就會導致分離得到的其中一個估計分量的相關系數(shù)較大，導致判斷錯誤。

(2) 列項標準差最小化

混合矩陣的列向標準差定義為[11]

(2)

盲源分離中分離矩陣W的逆矩陣H的列反映了源信號在各觀測信號中的變化情況，P300成分相對穩(wěn)定，變化差異較小，其對應的混合矩陣H的列向標準差較小。對H的3個列向量依次作列向分析，標準差最小的即為P300分量。但是，由于自發(fā)EEG和噪聲干擾的隨機性和不穩(wěn)定性，其對應的列向標準差最小的情況也是可能存在的。

針對上述問題，本文利用層次分析法，構(gòu)建了一個具有豐富層次結(jié)構(gòu)的分析評價模型，將從源信號的估計分量中有效、準確地提取P300分量轉(zhuǎn)化為一個最優(yōu)化選擇問題，其實質(zhì)就是依據(jù)相關系數(shù)及列項標準差，選擇一個最有可能是P300分量的部分。根據(jù)層次分析原理[12]，最高層次A即為P300分量的相似程度；第2層次B包含2個一級評價指標，列向標準差B1和相關系數(shù)B2；對于第3層次C，為便于計算，將列向標準差小、中、大分別定義為3個二級評價指標C1、C2、C3，將相關系數(shù)大、中、小分別定義為3個二級評價指標C4、C5、C6，依據(jù)表1中的1～9比較尺度構(gòu)造各層次的判斷矩陣，如表2～4所示。

表1 判斷矩陣標度及其含義

表2 一級評價指標判斷矩陣

表3 列向標準差層次下二級指標判斷矩陣

表4 相關系數(shù)層次下二級指標判斷矩陣

求解各判斷矩陣的特征向量，對一級評價指標、二級評價指標的權(quán)重進行計算，得到列向標準差B1的權(quán)重為0.333，相關系數(shù)B2的權(quán)重為0.667，C1～C6的權(quán)重分別為0.635、0.287、0.078、0.582、0.348、0.070，且判斷矩陣均滿足一致性準則。依據(jù)以下步驟對P300分量進行最優(yōu)提?。?/p>

(1) 依據(jù)式(2)，對混合矩陣H的3個列向量做列向分析，得到標準差序列σ1、σ2、σ3，按數(shù)值大小分別對其重新賦值為相應的權(quán)重C1～C3；

(2) 依據(jù)式(1)，將3個源信號的估計分量分別和靶刺激標準信號作相關分析，得到相關系數(shù)序列ρ1、ρ2、ρ3，按數(shù)值大小分別對其重新賦值為相應的權(quán)重C4-C6。

(3) 依據(jù)式

Wmax=max{B1×σi+B2×ρi}， i=1,2,3

(3)

權(quán)重最大值Wmax所對應的估計分量即為所要提取的P300分量。式中，權(quán)重最大值可取為0.600、0.484、0.444或0.414。根據(jù)上述方法選擇出P300分量后，將其余代表自發(fā)腦電和噪聲干擾的分量置零，利用混合矩陣對其進行重構(gòu)，完成對腦電信號P300成分的特征提取。

2.3 基于時域能量熵和小波變換的時頻補償

根據(jù)信息論中香濃熵的概念，定義時域能量熵

Ei=-log(Pi/P)

(4)

式中，Ei為第i段信號的時域能量熵，Pi為第i段信號的能量，P為信號的總能量。時域能量熵的大小描述了信號能量在時域上的分布情況：能量熵大，對應的時域段能量在信號總能量中的比重??；反之，在信號總能量中的比重大。

對前100，200，…，1200試次頻域濾波后的靶刺激和非靶刺激EEG信號分別進行累加平均，把得到的12個信號分別分割成20段，每段50ms，包含12個數(shù)據(jù)點，進行時域能量熵分析，如圖6、圖7所示。

圖6 靶刺激時域能量熵值的統(tǒng)計分布

圖7 非靶刺激時域能量熵值的統(tǒng)計分布

從圖中可以看出，靶刺激和非靶刺激時域能量熵的分布具有統(tǒng)計特性，而且差異性顯著，靶刺激信號在300ms左右的時域段能量較為集中，非靶刺激信號在整個時域上的能量分布較為平均。分別對20段的時域能量熵值進行處理，每段靶刺激和非靶刺激的各12個數(shù)據(jù)中，去掉2個最大和最小熵值，求出剩余的8個數(shù)據(jù)中的極大和極小熵值，這樣對于每一個段即產(chǎn)生靶刺激最大、靶刺激最小、非靶刺激最大和非靶刺激最小4個閾值，分別選擇其中的最大和最小熵值作為門限閾值，如表5所示。

表5 時域能量熵門限閾值分布表

如圖8、圖9所示，利用時域能量熵的統(tǒng)計分布特征對重構(gòu)后的信號進行補償，具體步驟如下：

(1) 對每一個重構(gòu)后的腦電信號進行時域能量熵處理，選擇門限閾值之外距離最遠的6段作為調(diào)整段。

圖8 經(jīng)過時頻補償后的時域能量熵值分布

圖9 經(jīng)過時頻補償后的時域分布

(2) 對其中較遠的5個調(diào)整段重新賦值為門限閾值，根據(jù)總能量P不變的等式求出第6個調(diào)整段調(diào)整之后的熵值。

(3) 利用Daubechies8小波，按照前述小波變換的頻域濾波方法，對調(diào)整之后的信號作平滑處理，完成第1次補償。

(4) 若補償后信號某一段的幅值超過400μV，則將該段的數(shù)據(jù)還原為調(diào)整前的大小，再次做平滑處理，完成第2次補償。

3 試驗結(jié)果及數(shù)據(jù)分析

3.1 靶刺激單導聯(lián)P300少試次提取相關程度分析

由于P300信號的潛伏期為300ms左右，因此這里只對前500ms內(nèi)的信號作相關分析處理。上述參考信號、經(jīng)過盲源分離JADE算法提取并重構(gòu)后的信號、第1次補償后的信號、第2次補償后的信號與標準信號的皮爾遜相關系數(shù)(PCC)依次為0.490、0.670、0.663和0.887，由此可見，經(jīng)過本文算法提取的靶刺激P300成分與標準信號高度相關。這里第1次時頻補償后相關系數(shù)較低主要是由于較大的波峰或波谷導致。

為進一步說明本文算法的有效性，在MATLAB R2010b Win32平臺下，隨機選取100組靶刺激信號作分析處理，并對重構(gòu)后PCC系數(shù)小于0.7的信號作時頻補償。如圖10、表6所示，依據(jù)本文算法，在導聯(lián)數(shù)量和試驗次數(shù)均較少的情況下仍能獲得較高的相似程度，相比于累加平均、盲源分離等方法有很大程度的提高。

圖10 提取的P300成分與標準信號相關系數(shù)試驗次數(shù)的統(tǒng)計分布對比

類別PCC>0．7(%)PCC>0．5(%)PCC>0．3(%)參考信號275476盲源分離重構(gòu)后的信號346084第1次補償后的信號496981第2次補償后的信號577989

3.2 靶刺激、非靶刺激單導聯(lián)P300少試次提取效果對比

如圖11所示，本文算法所提取的靶刺激、非靶刺激信號中的P300成分，相比于3試次累加平均和經(jīng)過盲源分離重構(gòu)后得到的結(jié)果，靶刺激波形相對于非靶刺激波形的P300成分更加明顯，更加接近于標準信號。

圖11 靶刺激、非靶刺激信號提取P300成分的效果對比

4 結(jié) 論

P300腦電信號作為事件相關電位(ERP)中重要的內(nèi)源性成分被廣泛應用于腦-機接口(BCI)系統(tǒng)中，但是臨床上使用的ERP檢測儀器均采用累加平均技術，忽略了不同試次之間的差異性和多試次引起的神經(jīng)性疲勞，造成操作繁瑣、數(shù)據(jù)冗雜等問題。本文提出的單導聯(lián)腦電信號P300少次提取算法，在采用小波變換與盲源分離JADE算法作濾波處理的同時，利用層次分析模型提取其中的P300分量，有效避免了少次提取中可能存在的突發(fā)性錯誤對提取結(jié)果造成的影響，并采用時域能量熵的統(tǒng)計分布特性對重構(gòu)后的信號進行時頻補償，進一步強化P300特征成分。本文方法在單通道、3試次的情況下，仍能夠有效地提取EEG中的P300成分，而且有效克服了多導聯(lián)、多試次產(chǎn)生的缺陷，為基于P300腦電信號的BCI系統(tǒng)在線應用打下基礎。

[1] Tu Y, Hung Y S, Hu L, et al. An automated and fast approach to detect single-trial visual evoked potentials with application to brain-computer interface.ClinicalNeurophysiology, 2014, 125(12): 2372-2383

[2] 郭凱, 裴為華, 王宇等. 基于干電極的頭帶式射頻無線腦-機接口系統(tǒng). 高技術通訊, 2012, 22(2): 211-216

[3] 程龍龍, 邱爽, 徐瑞等. 基于極能差與共空間模式的腦電信號特征增強研究. 高技術通訊, 2013, 23(9): 980-987

[4] Postelnicu C C, Talaba D. P300-based brain-neuronal computer interaction for spelling applications.IEEETransactionsonBiomedicalEngineering, 2013, 60(2): 534-543

[5] Vareka L, Mautner P. Off-line analysis of the P300 event-related potential using discrete wavelet transform. In: Proceedings of the 36th International Conference on Telecommunications and Signal Processing, Rome, Italy, 2013. 569-572

[6] Mousavi S A, Arshad M R H M, Mohamed H H, et al. P300 detection in electroencephalographic signals for brain-computer interface systems: a neural networks approach.ComputerEngineeringandNetworking, 2014, 277: 355-363

[7] Elsawy A S, Eldawlatly S, Taher M, et al. A principal component analysis ensemble classifier for P300 speller applications. In: Proceedings of the 8th International Symposium on Image and Signal Processing and Analysis, Trieste, Italy, 2013. 444-449

[8] 艾玲梅, 李營, 馬苗. 基于EMD和PCA的P300分類算法. 計算機工程, 2010, 36(5): 182-184

[9] Wessel J R, Ullsperger M. Selection of independent components representing event-related brain potentials: a data-driven approach for greater objectivity.NeuroImage, 2011, 54(3): 2105-2115

[10] 黃璐, 王宏. 單導聯(lián)腦電信號P300分類研究. 儀器儀表學報, 2014, 35(4): 814-819

[11] 陳洪波, 李蓓蕾, 陳真誠. 基于ICA的腦電信號P300少次自動提取. 電子學報, 2012, 40(6): 1257-1262

[12] Ishizaka A, Labib A. Review of the main developments in the analytic hierarchy process.ExpertSystemswithApplications, 2011, 38(11): 14336-14345

Study on extraction of P300 from monopolar lead EEG signal using several trials

Hu Chunhai， Xin Sixu， Liu Yonghong

(Measurement Technology and Instrumentation Key Lab of Hebei Province,Yanshan University, Qinhuangdao 066004)

The extraction of the positive potential P300 from electroencephalography (EEG) signals was studied. To vercome P300 extraction’s harsh requirements in electrode number and trial times, a new method for extraction of P300 from monopolar lead EEG signals using several trials was proposed. The method removes the noises in original signals by using wavelet transform and uses the Joint Approximate Diagonalization of Eigenmatrices (JADE), algorithm for blind source separation to decompose the observed signals, then, performs the optimal P300 extraction by using the analytic hierarchy process model built according to the pearson correlation coefficient (PCC) and the column standard deviation to avoid the possible sudden errors in P300 extraction in several trials, and finally, uses the statistical properties of temporal energy entropy and the wavelet transform principle to compensate reconstructed P300 signals in the time and frequency domain. As shown by the experimental results, only three trials of monopolar lead data need to be processed with the proposed method, and the P300 ingredient was extracted more effectively.

EEG processing, P300, temporal energy entropy, wavelet transform, analytic hierarchy process model

10.3772/j.issn.1002-0470.2016.04.006

①國家自然科學基金(61201110)資助項目。

2016-03-14)

②男，1966年生，博士，教授；研究方向：視覺檢測技術，信息融合，腦機接口等；聯(lián)系人，E-mail: fred-hu@ysu.edu.cn(