腦電中眼電偽跡的自動識別與去除

腦電中眼電偽跡的自動識別與去除李明愛;劉帆【摘要】目的為改善腦電中眼電偽跡的去除效果,基于腦電的非平穩(wěn)性和模糊特點,提出一種將離散小波變換與二階盲辨識相結合，并以模糊熵為眼電偽跡判別準則的眼電偽跡去除方法.方法首先,采用離散小波變換對含噪的腦電信號進行多分辯分析,獲得平穩(wěn)性更好的多尺度小波系數(shù);進而，選擇同層的小波系數(shù)構成小波系數(shù)矩陣，并基于二階盲辨識對其盲源分離,得到源信號的估計;進一步以模糊熵為判別依據(jù)，實現(xiàn)眼電偽跡的自動判別與剔除.實驗數(shù)據(jù)采用BCICompetitionIV公開數(shù)據(jù)庫，使用信噪比、相關系數(shù)及均方根誤差等常用偽跡判別指標進行衡量.結果本文方法相對于常用的眼電偽跡去除方法在多個性能指標上均取得最大值.結論本文提出的眼電偽跡去除方法，實現(xiàn)了眼電偽跡的自動精確判斷與剔除，并表現(xiàn)出很好的穩(wěn)定性.【期刊名稱】《北京生物醫(yī)學工程》【年(卷)，期】2018(037)006【總頁數(shù)】7頁(P559-565)【關鍵詞】腦電信號;眼電偽跡;離散小波變換;二階盲辨識;模糊熵【作者】李明愛瀏帆【作者單位】北京工業(yè)大學信息學部北京100124;計算智能與智能系統(tǒng)北京市重點實驗室北京100124;北京工業(yè)大學信息學部北京100124【正文語種】中文【中圖分類】R318.040引言腦電信號(electroencephalogram,EEG)是大腦中的電生理活動在大腦皮質(zhì)或頭皮表面神經(jīng)細胞的整體反映，其蘊含豐富的生理、心理及病理信息，在神經(jīng)科學、臨床醫(yī)學、康復工程等研究領域的作用日益顯著。但EEG信號微弱，在采集過程中不可避免地受到眼電(electrooculogram,EOG)、心電(electrocardiogram,ECG)、肌電(electromyogram,EMG)等噪聲干擾，特別是眼電信號的巨大幅值使得眼電偽跡(oculararitifial,OA)對EEG信號產(chǎn)生嚴重影響。如何有效去除眼電偽跡，獲得純凈的腦電信號具有重要的研究價值[1-2]。腦電信號的偽跡去除方法主要有兩類，小波變換法[3](wavelettransform,WT)和盲源分離法[4-5](blindsourcesseparation,BSS)。小波變換法利用WT具有良好的時頻特性和自適應特點，將含眼電偽跡的EEG源信號進行多尺度分解，通過將與眼電偽跡頻帶吻合的部分多尺度信號置零實現(xiàn)偽跡去除。由于EEG和EOG存在頻帶重疊現(xiàn)象，致使偽跡去除效果并不理想。盲源分離法利用BSS將含噪聲的EEG觀測信號分解成多個源信號分量，通過判斷并剔除與眼電偽跡相似或相關的部分信號分量，還原出純凈的腦電信號。常用的BSS包括：(1)主成分分析(principlecomponentanalysis,PCA)[6]O該方法假設源信號是非相關的，基于協(xié)方差矩陣對觀測信號進行正交分解，得到一組線性無關的變量，即主成分。當眼電信號與腦電信號成分具有相當?shù)姆禃r，不能徹底去除眼跡成分，無法得到理想的去噪效果。(2)獨立成分分析(independentcomponentanalysis,ICA)[7]OICA假設源信號是統(tǒng)計獨立的，且至多一個服從高斯分布，基于信號的二階及高階統(tǒng)計信息進行信號分解，從混合信號中恢復出原始信號的各統(tǒng)計獨立分量，計算量較大，且EEG難以滿足假設條件，去噪效果欠佳。FastICA是一種典型的ICA的快速實現(xiàn)算法，計算效率得到很好提升，但易受噪聲干擾，且存在迭代運算不收斂現(xiàn)象。⑶二階盲辨識(secondorderblindidentification,SOBI)［8-9］。該方法假設源信號為不相關的平穩(wěn)信號，利用觀測信號的二階統(tǒng)計量實現(xiàn)源信號的盲分離，不僅減少了計算量，而且增強了方法的魯棒性。腦電信號的統(tǒng)計特性與各種BSS方法假設條件的吻合程度是影響基于BSS的偽跡去除方法性能的重要因素。為了進一步改善眼跡去除效果，利用各種方法的優(yōu)勢，將兩種及以上去噪方法相結合是很好的選擇［10］。文獻［11］提出離散小波變換(discretewavelettransform,DWT)和FastICA相結合的偽跡去除方法，即DWICA，它將EEG信號進行小波變換，獲得比原始信號超高斯性更強、峰度更大的小波系數(shù)，進而對其進行獨立成分分析，使收斂速度和抗噪能力同時得以提高。影響眼跡去除效果的另一重要因素是偽跡判斷準則?；赪T的眼跡去除方法中，頻率范圍是衡量眼電偽跡成分的標準，由于不同受試者的眼電頻率范圍存在差異，且無法精確估算，與多尺度信號間頻率匹配時還會有誤差，導致眼電偽跡剔除不準確?；贐SS的眼跡去除方法中，通過計算眼電參考信號與信號分量間的相關系數(shù)或夾角余弦進行眼電偽跡判斷。相關系數(shù)與夾角余弦準則均具有較強的分辨能力，且計算速度快，但要求腦電信號的數(shù)據(jù)序列足夠長。文獻［12］根據(jù)眼電偽跡相對腦電信號的復雜度更低，因而具有更小的熵值，采用樣本熵來識別眼電偽跡。然而，樣本熵的相似度量公式為Heaviside二值函數(shù)，其突變性使得樣本熵值缺乏連續(xù)性，對閾值的取值非常敏感，其微弱變化，可能導致樣本熵值的突變，造成偽跡的錯誤判斷。本文提出將離散小波變換和二階盲辨識算法相結合的腦電中眼電偽跡的去除方法。利用離散小波變換對含眼電噪聲的非平穩(wěn)腦電信號進行多尺度分析，以放大突變的眼電偽跡信號；對同一尺度下具有更好平穩(wěn)性的小波系數(shù)［13］進行二階盲辨識，并引入模糊熵作為眼電偽跡的判別依據(jù)，實現(xiàn)眼電偽跡的自動判別與去除。1基本原理1.1二階盲辨識二階盲辨識是一種基于二階統(tǒng)計量的盲源分離方法。設信號線性混合模型如式(1)所示：x(t)=As(t)+n(t)(1)式中：s(t)=[s1(t),s2(t),...，sn(t)]T代表n個未知源信號；T代表矩陣的轉(zhuǎn)置；x(t)=[x1(t),x2(t),...,xm(t)]T為通過傳感器獲得的m個觀測信號；AuRmn表示線性混合矩陣，n(t)=[n1(t),n2(t),...,nm(t)]T表示噪聲信號。假設源信號s(t)為多變量平穩(wěn)信號，混合矩陣A為列滿秩矩陣，n(t)為獨立于源信號的零均值加性白噪聲。則二階盲辨識方法通過對一批協(xié)方差矩陣進行聯(lián)合近似對角化來求解分離矩陣，以獲得源信號s(t)的一個估計y(t)，主要步驟如下[8,14]:計算白化矩陣W,將觀測信號x(t)白化處理，使得z(t)的協(xié)方差矩陣是單位矩陣，除去各分量之間的相關性。z(t)=Wx(t)(2)對于固定延時TU｛Tj|j=12...,k｝，計算白化數(shù)據(jù)的采樣協(xié)方差矩陣R(t):R(T)=E[z(t+T)zT(t)]=ARz(T)AT⑶對所有的R(Tj),采用聯(lián)合近似對角化法，求得正交歸一化矩陣V:VTR(Tj)V=Dj(4)式中：｛Dj｝表示一組對角矩陣。計算混合矩陣A=W#V，求得去除相關后的源信號的估計y(t):y(t)=A-1x(t)=VTWx(t)⑸1.2離散小波變換小波變換是時間、頻率的局部化分析方法，通過伸縮平移實現(xiàn)對信號的逐步細化，反映信號在時域和頻域的局部特征。因具有多分辨分析和對信號的自適應特點，小波變換特別適合于處理EEG這類非平穩(wěn)信號。在實際中，為方便計算，減少運算量，通常將處理的信號離散化，此即為離散小波變換［15］。在多分辨分析基礎上，Mallat提出離散小波快速算法，即Mallat算法［16］。利用Mallat算法可以根據(jù)不同尺度進行多級分解，把源信號分解成細節(jié)分量和逼近分量，如式(6)所示。(6)式中：an,k、dn,k分別代表第n層的逼近分量和細節(jié)分量。Mallat塔式重構算法是分解算法的逆過程，如式(7)所示，當n逐級遞減至1時，完成對信號的重構。⑺1.3模糊熵模糊熵(fuzzyentropy，F(xiàn)E)［17］是一種時間序列復雜度的測量方法。熵值越小，序列越簡單；熵值越大，序列越復雜。針對眼電偽跡頻率低、復雜性低，其模糊熵值相對腦電信號更小的特點，通過計算FE可實現(xiàn)偽跡判別，計算式如式(8)。(8)式中：em(n,r)由相似度計算得到；m為嵌入維數(shù)；n為指數(shù)邊界的梯度；r為相似容限。具體過程詳見文獻[16]。2基于二階盲辨識和模糊熵的眼電偽跡自動去除方法本文提出一種眼電偽跡自動去除方法，將二階盲辨識方法與離散小波變換相結合，實現(xiàn)測量信號在小波域的盲源分離，并采用模糊熵進行眼電偽跡判別，具體步驟如下。2.1基于sym8小波基函數(shù)，對原始信號進行離散小波變換分解假設x(t)=[x1(t),x2(t),...,xm(t)]T表示通過電極采集的m導信號，xi(t)為第i導腦電信號，i=12...,m-1,xm(t)為一導參考眼電信號。采用sym8小波母函數(shù)，對m導測量信號分別進行3層離散小波分解，并將得到的逼近分量和細節(jié)分量構成小波系數(shù)向量Cx,i(t)=[Ai3Di3Di2Di,1](i=12...,m),進而構造小波系數(shù)矩陣Cx(t)=[Cx,1(t),Cx,2(t),...,Cx,m(t)]T。2.2基于SOBI算法的小波系數(shù)矩陣盲源分離采用二階盲辨識算法，對小波系數(shù)矩陣進行盲源分離，求得其源信號的估計，記為Cs(t)=[Cs,1(t),Cs,2(t),...,Cs,m(t)]T。2.3基于模糊熵的眼電偽跡判別與去除計算各個源估計信號Cs,i(t)的模糊熵值FE(Cs,i(t))(i=12...,m),對其進行排序，其中具有最小模糊熵值的源估計信號即自動判定為眼電偽跡。將其置零，其余源估計信號保持不變，以達到偽跡剔除的目的。2.4基于SOBI投影逆變換假設眼電偽跡去除后的源估計信號記為則基于SOBI算法對其投影變換后的信號為2.5小波的逆變換基于式(7)所示Mallat塔式重構算法，對信號進行小波逆變換，還原出偽跡去除后的頭皮腦電信號xc(t)。3實驗研究與結果分析3.1實驗數(shù)據(jù)3.1.1數(shù)據(jù)來源實驗數(shù)據(jù)來自公開的BCICompetitionIV的DataSets2b數(shù)據(jù)庫。采用Ag/AgCl電極，電極位置服從國際標準的10~20導聯(lián)系統(tǒng)分布［18］。數(shù)據(jù)記錄9名右利手、視覺良好的健康受試者。每個受試者的實驗數(shù)據(jù)采集包含3個階段，如圖1所示。(1)空載階段，包括2min睜眼和1min閉眼。(2)眼電采集階段80s，包括眨眼、轉(zhuǎn)動眼球、上下眼動和左右眼動各15s，受試者根據(jù)屏幕上的提示，依次進行。每種眼電偽跡采集完成時，隨后休息5s進入下一階段。(3)腦電采集階段，采集左、右手運動想像EEG各60次共計120組數(shù)據(jù)。選擇C3、Cz和C4,3個電極位置記錄腦電數(shù)據(jù)(Cz為C3和C4的參考電極)，采樣頻率為250Hz,并將采集的信號通過0.5~100Hz的帶通濾波及50Hz的陷波濾波器濾除工頻干擾。采集實驗時序如圖1所示。圖1實驗采集時序Figure1Thetimingdiagramofcollectionexperiment3.1.2實驗數(shù)據(jù)構造為了精確衡量腦電信號中眼電偽跡的去除效果，這里將根據(jù)腦電信號和眼電信號的雙向激活擴散特性［19］，利用純凈的眼電信號和腦電信號構造含噪聲的眼電信號和腦電信號。模型表達式如下：EOGraw(t)=EOGclean(t)+n1EEGclean(t)(9)EEGraw(t)=EEGclean(t)+n2EOGclean(t)+WGN(10)式中:EOGraw(t)、EEGraw(t)分別代表實際采集到的含噪聲的眼電信號和腦電信號;EOGclean(t).EEGclean(t)代表純凈的眼電信號和腦電信號;WGN為高斯白噪聲；系數(shù)n1和n2為干擾影響系數(shù)，n1u(0,0.5］,n2u(0,0.2］,干擾系數(shù)均隨機產(chǎn)生，以保證實驗結果的可信度。同時，為模擬肌電信號、心電信號等其他噪聲對腦電信號的干擾，在C3，Cz，C4三導腦電信號中加入10dBW的高斯白噪聲。3.2眼電偽跡去除性能評價指標3.2.1均方根誤差均方根誤差(rootmeansquarederror，RMSE)是一種衡量偽跡去除效果最常用的性能指標。假設xclean(t)和xc(t)分別表示某一導聯(lián)上的純凈腦電信號和偽跡去除后的腦電信號，則均方根誤差計算式為：(11)式中:N為采樣點個數(shù)。RMSE的值越小，說明偽跡去除后的EEG越接近純凈的EEG，即去除效果越好。3.2.2信噪比第二種常用的偽跡去除效果評價指標為信噪比(signalnoiserate，SNR)，其計算式為：(12)xclean(t)和xc(t)式中標量含義同上。SNR的值越大，意味著xclean(t)和xc(t)越接近，說明去噪效果越好。3.2.3相關系數(shù)相關系數(shù)(correlationcoefficient,CC)是用以反映變量之間相關關系密切程度的統(tǒng)計指標，是用于評價偽跡去除效果的又一常用指標，計算式為：(13)式中：分別表示信號xclean(t)和xc(t)的均值。CC數(shù)值越大，表示xclean(t)和xc(t)之間相似度越高，即去噪效果越理想。3.3實驗結果與分析3.3.1眼電偽跡去除的實驗研究基于數(shù)據(jù)庫中某受試者(subject1)的水平眼電(horizontalEOG,H-EOG)和運動想像腦電實驗數(shù)據(jù)，將本文方法與常用眼電偽跡去除方法進行對比實驗研究，并使用多種性能指標進行性能評價［9］。為克服實驗數(shù)據(jù)的偶然性對實驗結果的影響，依式(9)和式(10)構造含噪眼電和運動想像腦電的實驗數(shù)據(jù)時，相互干擾系數(shù)n1、n2分別隨機產(chǎn)生10組。針對每組干擾系數(shù)，從120組運動想像腦電數(shù)據(jù)中任選一組與水平眼電分別截取等長的750個采樣點，構造含噪的運動想像腦電和眼電信號，并構成4導源信號進行實驗研究。采用WT、PCA、ICA等多種方法去除眼電偽跡，計算偽跡去除后的腦電信號與純凈腦電信號的10組實驗均方根誤差的平均值如表1所示，WT方法去除水平眼電偽跡時效果最差，這是由于WT方法無法針對EEG信號與EOG信號的頻帶混疊問題進行偽跡有效去除；基于PCA、ICA、SOBI的水平眼電偽跡去除效果得到明顯改善，且ICA和SOBI的效果較為接近并優(yōu)于PCA,這說明腦電源信號和觀測信號能夠較好地滿足盲源分離的原理和條件，但眼電偽跡去除效果還取決于腦電源信號的統(tǒng)計特點和分布特性；基于DWICA的去除效果進一步得到提高，這是由于小波系數(shù)相對于EEG源信號超高斯性更強；本文方法(DWSOFE)去除水平眼電偽跡效果最佳，在C3、Cz、C4中任意導聯(lián)上均方根誤差都取得相對最小值，這主要是由于該方法利用小波變換對非平穩(wěn)的MI-EEG信號進行多分辨分析，獲得了具有較好平穩(wěn)性的小波系數(shù)，使之更加吻合SOBI假設條件，有效增強了對于EEG信號的適應性。表1多種方法去除H-EOG后的平均RMSE比較(單位：^V)Table1AverageRMSEcomparisonofmultiplemethodsafterH-EOGremoval(unit:pV)導聯(lián)WTPCAICASOBIDWICADWSOFEC320.4111.797.247.367.176.40Cz26.0413.978.298.958.237.41C421.2912.007.607.076.976.08多種方法用于去除水平眼電在各導聯(lián)上的平均信噪比和相關系數(shù)如圖2和圖3所示。從圖2可知，WT方法獲得的平均信噪比最低，ICA和PCA方法獲得的平均信噪比依次得到極大提高，SOBI和DWICA的效果進一步得以改善，且效果相當，而本文方法得到最高平均信噪比。另外，從圖3可見，WT方法在三導聯(lián)上的平均相關系數(shù)指標均取得最小值，PCA和SOBI得到的平均相關系數(shù)依次得以明顯提升，ICA和DWICA以微弱優(yōu)勢得以進一步提高，DWSOFE獲得最高平均相關系數(shù)。充分展現(xiàn)了DWSOFE在以相關系數(shù)為偽跡去除性能指標時的絕對優(yōu)勢和更好的結果穩(wěn)定性。綜上可見，采用多種指標衡量水平眼電偽跡去除效果時，不同方法的表現(xiàn)并不具有完全一致性。例如，當采用RMSE和CC性能指標時，ICA效果好于PCA，而當采用SNR指標時，PCA卻優(yōu)于ICA，這與腦電源信號的統(tǒng)計特點以及分布特性密切相關。但對3種性能指標而言，本文方法均獲得最佳效果，足以證明本文方法去除水平眼電偽跡的有效性和穩(wěn)定性。圖2多種方法去除H-EOG后的平均SNR比較Figure2AverageSNRcomparisonofmultiplemethodsafterH-EOGremoval3.3.2眼電偽跡判別方法對偽跡去除效果影響的實驗研究表2眼電偽跡判別方法對偽跡去除效果的影響Table2ComparisonofOAremovaleffectfordifferentmethods評價指標RMSE/pVSNR/dBCC判別方法C3CzC4C3CzC4C3CzC4相關系數(shù)9.17610.6499.11511.75812.89412.3270.9220.9140.924夾角余弦7.0427.8306.74513.23614.50613.8830.9640.9690.969樣本熵7.0387.8266.74013.24014.51213.8890.9650.9690.970模糊熵6.4007.4106.08014.09014.79014.7400.9720.9740.976圖3多種方法去除H-EOG后的平均CC比較Figure3AverageCCcomparisonofmultiplemethodsafterH-EOGremoval腦電信號與眼電信號的復雜度不同，且腦電信號具有模糊特性，為此，本文通過計算模糊熵值來達到偽跡辨別的目的。采用3.3.1節(jié)中構造的含噪聲的眼電數(shù)據(jù)和運動想像腦電數(shù)據(jù)，依據(jù)本文提出的SOBI和WT相結合的偽跡去除方法，僅改變眼電偽跡的判別準則，實現(xiàn)模糊熵和相關系數(shù)、夾角余弦、樣本熵等常用偽跡判別準則的對比研究，基于3種評價指標的平均實驗結果如表2所示。表2表明，以相關系數(shù)進行眼電偽跡效果相對較差，夾角余弦與樣本熵的實驗結果有明顯改善，且二者效果較為接近，而模糊熵則對3種評價指標均取得最好的效果。究其原因，在模糊熵的相似度量計算式中采用指數(shù)函數(shù)，使熵值變化更加平穩(wěn)平滑，不易造成偽跡判別錯誤，因而有利于提高偽跡識別的準確性。4討論與結論本文提出一種二階盲辨識結合離散小波變換的眼電偽跡去除方法，并引入模糊熵實現(xiàn)偽跡自動判別的同時，有效增強了偽跡判別的準確性和偽跡去除的精度?；趪H公開腦電和眼電數(shù)據(jù)庫進行實驗研究，結果表明，采用均方根誤差、信噪比和相關系數(shù)等多項性能指標時，本文方法在去除水平眼電偽跡時均有明顯優(yōu)勢，對腦電信號的在線預處理具有重要意義。在未來的工作中，將加強本文方法用于去除多個受試者垂直眼電、眼睛轉(zhuǎn)動等多種眼電偽跡的適用性研究，并從統(tǒng)計學的角度驗證本文方法的顯著性。參考文獻【相關文獻】BurgerC,vandenHeeverDJ.RemovalofEOGartefactsbycombiningwaveletneuralnetworkandindependentcomponentanalysis[J].BiomedicalSignalProcessingandControl,2015,15:67-79.QnarS,AcirN.AnovelsystemforautomaticremovalofocularartefactsinEEGbyusingoutlierdetectionmethodsandindependentcomponentanalysis[J].ExpertSystemswithApplications,2017,68:36-44.KrishnaveniV,JayaramanS,AravindS,etal.AutomaticidentificationandremovalofocularartifactsfromEEGusingwavelettransform[J].MeasurementScienceReview,2008,6(4):45-57.OosugiN,KitajoK,HasegawaN,etal.AnewmethodforquantifyingtheperformanceofEEGblindsourceseparationalgorithmsbyreferencingasimultaneouslyrecordedECoGsignal[J].NeuralNetworkstheofficialjournaloftheInternationalNeuralNetworkSociety,2017,93:1-6.FitzgibbonSP,PowersDM,PopeKJ,etal.RemovalofEEGnoiseandartifactusingblindsourceseparation[J].JournalofClinicalNeurophysiology:officialpublicationoftheAmericanElectroencephalographicSociety,2007,24(3):232-243.CasarottoS,BianchiAM,CeruttiS,etal.Principalcomponentanalysisforreductionofocularartefactsinevent-relatedpotentialsofnormalanddyslexicchildren[J].ClinicalNeurophysiology,2004,115(3):609-619.KongW,ZhouZ,HuS,etal.AutomaticanddirectidentificationofblinkcomponentsfromscalpEEG[J].Sensors,2013,13(8):10783-10801.BelouchraniA,Abed-MeraimK,CardosoJF,etal.Ablindsourceseparationtechniqueusingsecond-orderstatistics[J].IEEETransactionsonSignalProcessing,1997,45(2):434-444.BridwellDA,RachakondaS,SilvaRF,etal.Spatiospectraldecompositionofmultisubjecteeg:evaluatingblindsourceseparationalgorithmsonrealandrealisticsimulateddata[J].BrainTopography,2016,31(1):47-61.UriguenJA,Garcia-ZapirainB.EEGartifactremoval-state-of-the-artandguidelines[J].JournalofNeuralEngineering,2015,12(3):031001.LiM,CuiY,YangJ.AutomaticremovalofocularartifactfromEEGwithDWTandICAmethod[J].AppliedMathematics&InformationSciences,2013,7(2):809-816.MowlaMR,NgSC,ZilanyMSA,etal.Artifacts-matchedblindsourceseparationandwavelettransformformultichannelEEGdenoising[J].BiomedicalSignalProcessingandControl,2015,22:111-118.黎志勇，李寧.基于小波的非平穩(wěn)時間序列預測方法研究[J].計算機工程與應用,2014,50(10):38-43.LiZY,LiN.Researchonnon-stationarytimeseriespredictionbasedonwavelet[J].ComputerEngineeringandApplications,2014,50(10):38-43.BrewickPT,SmythAW.Increasingtheefficiencyandefficac