非侵入式腦-機接口已經逐步成為當前研究的熱點，在精神障礙檢測、生理監測等多方面都有所應用。但是非侵入式腦-機接口所需的腦電信號容易受到眼電偽跡污染，會嚴重影響對腦電信號的解碼分析。對此，本文提出了一種結合頻率濾波器的改進型獨立成分分析算法，以相關系數和峰度雙重閾值為依據自動識別偽跡組件；利用眼電與腦電頻率的差異，通過頻率濾波器去除偽跡組件中的眼電信息，從而保留更多腦電信息。在公開數據集和本實驗室數據上的實驗結果表明，本文算法可以有效提升眼電偽跡去除效果，同時改善腦電信息損失，這有助于非侵入式腦-機接口的推廣。
引用本文： 宋昊, 徐頌, 劉國明, 劉京, 熊鵬. 基于獨立成分分析的非侵入式腦-機接口眼電偽跡自動去除算法. 生物醫學工程學雜志, 2022, 39(6): 1074-1081. doi: 10.7507/1001-5515.202111060 復制
引言
腦-機接口（brain computer interface，BCI）技術是近些年的一大研究熱點，該技術的應用為許多視聽功能障礙患者帶來了希望，也為癱瘓患者帶來了重新運動的可能。醫學上為了獲取更加精確的腦電（electroencephalography，EEG）信號，通常采用侵入式BCI，將電極植入特定腦區。這種方式雖然采集的腦電信號更加干凈有效，但是存在一定的手術風險，不利于BCI技術的推廣。非侵入式BCI的使用在一定程度上解決了這個問題。目前，由于安全、便捷、無創等原因，非侵入式BCI已被廣泛應用于精神障礙預測[1]、藥物作用評估[2]、生理狀態監測[3]、神經系統疾病診斷[4]與康復[5]等領域。但同時，非侵入性BCI容易受到生理偽跡干擾的缺點也暴露了出來，最常見的有心電偽跡、肌電偽跡和眼電（electrooculographic，EOG）偽跡等[6]。其中，對腦電信號影響最大的是眼電偽跡，眨眼時產生的眼電信號電勢峰值可比腦電信號大十倍以上，并且單個峰持續時間長達200～400 ms，它產生的電勢可以蔓延到所有電極[7]。眼電偽跡的存在，將導致腦電信號的嚴重失真，極大地影響對腦電信號的解碼分析，甚至導致對疾病的診斷出現偏差和誤判。因此，如何在腦電信號中有效去除眼電偽跡是一個十分重要而且必要的問題。
通常眼電偽跡的去除使用數據丟棄法和偽跡消除法[8]。數據丟棄法最為簡單，即人工檢查采集到的腦電信號，然后丟棄受到眼電偽跡污染的片段，這將導致EEG信息的丟失。與使用受到EOG污染的EEG相比，丟失大量EEG信息可能存在更大的問題。并且，在非侵入式BCI中該類方法難以應用。偽跡消除法即識別和去除腦電信號中的偽跡，并保留原始信號的神經信息，主要分為兩類。第一類是采用回歸和濾波的方法，Fatourechi等[9]指出，眼電偽跡主要存在于低頻波段，可以通過高通濾波的方式予以去除；基于自回歸模型和自適應濾波的方法也已經應用于EOG偽跡去除，此類方法較為依賴參考通道的數據質量[10-11]。第二類是通過將腦電信號和偽跡信號分離或者分解到其他域進行去除[11]。此類方法較多，較為著名的是小波變換（wavelet transform，WT），將信號從時域映射到小波域，通過消除噪聲的小波系數消除噪聲[12]；使用最為廣泛的是獨立成分分析（independent component analysis，ICA），該方法主要原理是假設所有導聯的信號都來自于特定的源，觀測信號經過線性變換可以得到若干個源，通過去除EOG的源消除EOG偽跡；作為一種全數據驅動、高度自適應的非線性非平穩信號分析技術，經驗模態分解（empirical mode decomposition，EMD）和針對噪聲敏感問題改進的完備經驗模態分解（ensemble empirical mode decomposition，EEMD）也得到了廣泛應用[13]。
現階段更為主流的方法是將以上方法中的兩種或三種結合使用。樣本熵-完備經驗模態分解（sample entropy complete ensemble empirical mode decomposition with adaptive noise，SE-CEEMDAN）算法可以避免回歸算法中的雙向污染問題，但是同時引入了ICA、WT與樣本熵，大大增加了計算的復雜度[14]。基于ICA的算法中，去除包含偽跡的整個獨立成分（independent component，IC）的做法可能會造成某些混合在偽跡IC中的部分重要EEG信息的丟失[15-18]。Pion-Tonachini等[19]利用在線遞歸ICA等方法針對在線去除眼電偽跡做了改進。Castellanos等[18]提出的小波ICA方法（wavelet independent component analysis，wICA）可以減少信息損失，并且效果好過整個IC去除的方法，但是小波變換很大程度上取決于母波的選擇，而母波通常是包含正弦和余弦的函數，因此大多數時候不適用于非線性的EEG信號[7]。
針對上述問題，考慮到眼電偽跡主要存在于低頻波段[20]，為了更多地保留EEG信息，更好地適應非侵入式BCI，本文提出了一種結合ICA與頻率濾波器消除眼電偽跡的方法。首先使用ICA算法對獲取的原始數據進行分解以獲得不同的IC，之后采用相關系數與峰度雙重閾值判定的方式自動識別偽跡IC[21-22]，然后使用高通濾波器對偽跡IC進行濾波處理，最后通過逆ICA重構出觀測信號。
1 基于獨立成分分析的眼電偽跡自動去除算法
1.1 算法整體設計
本文提出的方法整體流程圖如圖1所示。主要處理步驟可分為三部分：首先對原始數據預處理，然后通過ICA算法分離并正確識別偽跡IC，最后通過帶通濾波器濾除偽跡成分并重構以獲得干凈的EEG。
圖1
算法流程圖
Figure1.
The flow chart of algorithm
1.2 獨立成分分析及眼電偽跡識別
ICA的提出最早是為了解決盲源分離問題，通常指可以根據高階統計量從線性混合中恢復統計獨立成分的一類算法，目前被認為是用于識別和去除腦電信號中偽跡魯棒性較強的算法，其數學表達式如式(1)所示：
 |
其中，A稱為混合矩陣，
表示源信號，
表示觀測信號。
對經過預處理（0.5 Hz高通濾波，設置參考電極為M1、M2）的原始信號使用Infomax ICA進行分解得到若干IC，引入皮爾遜相關系數、z分數閾值和峰度閾值以正確識別偽跡IC。由于EOG是EEG強度的幾十到幾百倍，偽跡IC與EOG通道信號的相關性要遠高于非偽跡IC與EOG通道信號的相關性，因此可以使用皮爾遜相關性對二者進行區分。z分數使得閾值設置具有自適應性，降低了直接使用相關性閾值判定導致誤判的可能性。z分數閾值公式為：
 |
其中，X為原始數據，
為平均數，S 為標準差。
通過以上計算，可以得到每個IC與EOG通道數據的相關性，并且通過z分數自動標記出候選偽跡IC。如圖2所示，橫軸表示IC數量，縱軸表示相關性；標紅部分為滿足閾值條件的IC，即候選偽跡IC。
圖2
ICs與垂直眼電（上）通道、水平眼電通道（下）數據的相關性得分
Figure2.
The correlation score of ICs with the data of vertical EOG (VEOG) (upper) and horizontal EOG (HEOG) (lower) channels
通過z分數閾值檢測可能仍然會對某些與EOG相似的IC產生誤判，信號的高階統計可以反映更多低階統計所不具備的信息。四階統計量峰度是具有簡單計算和理論的信號峰值的度量。通常，高斯信號的峰度為0，亞高斯信號峰度為負值，而EOG等超高斯信號峰度為正值，因此可以用峰度對候選偽跡IC做最終判定[22-24]。針對不同的數據集，峰度的閾值可以通過多次驗證確定，本文中峰度閾值為1.5[25]。若已知x是一個標量隨機變量，則峰度計算公式為：
 |
其中，
表示變量x的n階中心距，
為隨機變量x的均值。
1.3 眼電偽跡去除
EEG失真的主要原因是受到了EOG的污染，通過ICA分解可以獲得含有EEG信息的偽跡IC。去除偽跡IC中的EOG成分即可保留其EEG信息，從而減少重要腦電信息的丟失。由于腦電信號的主要頻帶為0.1～40 Hz，眼電信號的主要頻帶在7 Hz以下[20]，濾除存在于偽跡IC中的眼電成分即可保留絕大部分腦電信息。
2 實驗數據及評價指標
2.1 實驗數據
2.1.1 公開數據集
選用2016年Klados和Bamidis公開發表的數據集[26]（https://data.mendeley.com/datasets/wb6yvr725d/1）。該數據采集實驗使用了19通道電極帽，電極排布符合國際10-20規則，以200 Hz的采樣頻率記錄了27名受試者的閉眼腦電數據，每段數據時長28 s。該數據集還提供了垂直眼電和水平眼電的數據，并且根據式(4)將眼電數據添加到了腦電信號中：
 |
其中，
表示受到眼電偽跡污染的腦電信號，
表示無眼電偽跡的腦電信號，
表示垂直眼電信號，
表示水平眼電信號，
表示垂直眼電污染系數，
表示水平眼電污染系數，i表示受試者編號，j表示通道編號。
2.1.2 本實驗室數據
本實驗室數據來自河北省數字醫療工程重點實驗室，該實驗室位于河北大學新校區C座1號樓一層。使用40導便攜式腦電信號放大器（Grael，Compumedics Ltd.，澳大利亞）以及數據采集和在線處理軟件CURRY8 X（Compumedics Ltd.，澳大利亞）采集了12名（9名男生，3名女生）年齡19～25周歲在校大學生執行目標識別任務時的EEG數據。電極排布依據國際10-20標準，包含32個EEG導聯和2個EOG導聯（垂直EOG通道VEOG和水平EOG通道HEOG），數據采樣頻率為512 Hz。所有受試人員均由招募而來，視力正常或矯正視力正常，均為右利手，且均簽署了知情同意書。同時，該數據采集實驗通過了河北大學附屬醫院倫理委員會審查。
2.2 評價方法
有研究對EOG偽跡去除前后的數據進行了定量分析，主要評判內容為：腦電有效信息丟失程度與眼電偽跡去除效果[27]。根據數據的不同，本文對公開數據集和本驗室數據采用了不同的評價指標。
2.2.1 公開數據集評價
相關系數（correlation coefficient，CC）是在時域上判斷信號失真的一個重要指標，表征了兩個變量之間的相關程度。通過計算去除偽跡之后數據與純凈腦電之間的相關性，可以有效量化腦電信號丟失程度，相關系數越大，表明腦電信號保留越完整。均方根誤差（root mean square error，RMSE）表征了眼電偽跡去除效果，該值越小，表示眼電偽跡去除越完全。
 '/> |
 '/> |
其中，
表示EOG偽跡去除之后的EEG信號，x 表示純凈的EEG信號。
2.2.2 本實驗室數據評價
由于本實驗室數據無法獲得純凈的EEG，無法計算RSME和CC，因此采用幅度平方相干性（magnitude squared coherence，MSC）作為指標。MSC是定量評價偽跡去除前后頻域信息變化的常用度量[28-31]，通過相干性得分，可以有效評估在不同頻率上信息的保留程度。主要包含眼電信息的低頻段相干性得分越低，表明眼電信息去除越多，其余頻段相干性得分越高，表明腦電信息保留越完善。計算方法如式(7)所示：
 |
其中，
是信號x與信號y在頻率 
下的互譜密度，
和 
分別指信號x和信號y在頻率 
下的自譜密度。
3 結果與討論
3.1 公開數據集結果分析
本文提出的方法在公開數據集中的35組數據上進行了應用，每組數據都選取全部時長。圖3展示了公開數據集中數據1的FP1通道部分數據在EOG偽跡去除前后與純凈EEG的對比結果，藍線表示受到EOG偽跡污染的EEG，橙線表示經本文算法處理之后的EEG，綠線表示純凈EEG。如圖3所示，經過本文算法的處理，由EOG干擾產生的EEG波形畸變得到了很好的校正，說明該算法可以實現單個EOG偽跡的良好去除。
圖3
數據1 FP1通道純凈EEG、受污染的EEG和偽跡去除后EEG的對比
Figure3.
Comparison of the artifact-free, the contaminated and the cleaned EEG signals of dataset1, channel FP1
如圖4所示，由于在數據采集過程中受試者多次產生了眨眼活動，原始數據受到EOG干擾而多通道多處呈現出嚴重的波形畸變。經過本文算法的處理，所有通道的多處形變都得到了校正，說明本算法可適用于多通道多處EOG偽跡的自動化消除。
圖4
公開數據集眼電偽跡去除效果
Figure4.
The effect of EOG artifact removal in open datasets
為了評估該算法的有效性和穩定性，眼電偽跡去除之后，對每組數據都計算了相關系數與均方根誤差，并與SE-CEEMDAN、小波增強ICA進行了對比。如表1所示，受到眼電偽跡污染的EEG信號與純凈腦電信號之間的均方根誤差很大、相關系數很小，通過三種算法處理，實現了不同程度的改善。從結果看，本文算法平均RSME優于SE-CEEMDAN和小波增強ICA，即本文算法能夠去除更多的偽跡信息；平均CC明顯高于SE-CEEMDAN，即本文算法能保留更多的EEG信息。該結果表明，本文算法處理后的數據更加接近純凈EEG。
表1
不同算法性能對比（均值±標準差）
Table1.
Performance comparison of different methods (mean±SD)
3.2 本實驗室數據分析
本實驗室數據可以通過分析偽跡去除前后的MSC來定量分析EEG成分的保留程度。另外，由于EOG偽跡本身存在高振幅、低頻率的特點，可以通過觀察波形圖定性分析EOG偽跡去除效果。圖5展示了本文方法處理前后的波形圖。從圖中可以看出，采集到的原始數據由于受到EOG的影響而多處產生嚴重形變。經過本文算法處理后，在13.5 s附近和19 s附近嚴重的形變消失，21 s附近較弱的形變也得到了很好的校正，而其他沒有EOG偽跡的部分未受到影響。說明本文的算法對于EOG偽跡具有較高的識別度，并且能夠很好地去除不同幅度的EOG偽跡。
圖5
本實驗室數據眼電偽跡去除效果
Figure5.
The effect of removal of EOG artifacts from our laboratory data
人工檢查發現本實驗室獲取的12組數據均含有EOG偽跡，對每名受試者的數據都進行了EOG偽跡去除處理。由于眼睛距離額葉最近，因此額葉通道受到EOG影響更大，此處計算了EOG偽跡去除前后額葉通道FP1的MSC。由于MSC是一個表征信號在不同頻率上差異的量，而EOG偽跡主要存在于低頻波段，將其去除后與原始信號在低頻段差異較大，因此低頻波段MSC得分應該較低；相反，在高頻波段MSC得分應該較高。如圖6所示，低頻部分本文算法MSC得分較小波增強ICA更低，表明本文算法去除了更多的眼電偽跡，而高頻部分本文算法MSC得分較小波增強ICA更高，表明本文算法對高頻部分EEG信息的影響更小。
圖6
本實驗室數據經不同算法處理后的MSC結果
Figure6.
MSC results of the laboratory data processed by different algorithms
4 結論
基于非侵入式BCI快速發展卻容易受到EOG偽跡影響的現狀，針對現有眼電偽跡去除算法存在EEG信息丟失的問題，本著自動、高效的原則，本文提出了一種基于ICA的眼電偽跡自動去除算法。該算法通過z分數與峰度雙重閾值自動判斷偽跡IC，采用頻率濾波器去除EOG偽跡，實現了EOG偽跡的自動選擇和去除，提高了數據處理效率。在公開數據集和本實驗室數據上分別進行了測試，公開數據集結果顯示，相比于SE-CEEMDAN和小波增強ICA算法，本文提出的方法具有更小的RSME和更大的CC，表明該算法能夠更加有效地去除眼電偽跡，同時能夠保留更多的EEG信息；本實驗室數據分析結果表明，本文所提出的算法對不同幅度的EOG偽跡都實現了良好的去除，又幾乎不影響高頻EEG信息。綜上所述，該算法實現了EOG偽跡更大程度的自動去除，又大大減少了EEG信息的丟失，可能為非侵入式BCI的推廣起到一定助力作用。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：宋昊主要負責數據收集與分析、算法程序設計以及論文編寫；徐頌、劉國明主要負責實驗素材收集；劉京主要負責實驗指導；熊鵬主要負責實驗設計、數據分析指導以及論文審閱修訂。
倫理聲明：本研究通過了河北大學附屬醫院倫理委員會的審批（批文編號：HDFY-LL-2019-035）。
引言
腦-機接口（brain computer interface，BCI）技術是近些年的一大研究熱點，該技術的應用為許多視聽功能障礙患者帶來了希望，也為癱瘓患者帶來了重新運動的可能。醫學上為了獲取更加精確的腦電（electroencephalography，EEG）信號，通常采用侵入式BCI，將電極植入特定腦區。這種方式雖然采集的腦電信號更加干凈有效，但是存在一定的手術風險，不利于BCI技術的推廣。非侵入式BCI的使用在一定程度上解決了這個問題。目前，由于安全、便捷、無創等原因，非侵入式BCI已被廣泛應用于精神障礙預測[1]、藥物作用評估[2]、生理狀態監測[3]、神經系統疾病診斷[4]與康復[5]等領域。但同時，非侵入性BCI容易受到生理偽跡干擾的缺點也暴露了出來，最常見的有心電偽跡、肌電偽跡和眼電（electrooculographic，EOG）偽跡等[6]。其中，對腦電信號影響最大的是眼電偽跡，眨眼時產生的眼電信號電勢峰值可比腦電信號大十倍以上，并且單個峰持續時間長達200～400 ms，它產生的電勢可以蔓延到所有電極[7]。眼電偽跡的存在，將導致腦電信號的嚴重失真，極大地影響對腦電信號的解碼分析，甚至導致對疾病的診斷出現偏差和誤判。因此，如何在腦電信號中有效去除眼電偽跡是一個十分重要而且必要的問題。
通常眼電偽跡的去除使用數據丟棄法和偽跡消除法[8]。數據丟棄法最為簡單，即人工檢查采集到的腦電信號，然后丟棄受到眼電偽跡污染的片段，這將導致EEG信息的丟失。與使用受到EOG污染的EEG相比，丟失大量EEG信息可能存在更大的問題。并且，在非侵入式BCI中該類方法難以應用。偽跡消除法即識別和去除腦電信號中的偽跡，并保留原始信號的神經信息，主要分為兩類。第一類是采用回歸和濾波的方法，Fatourechi等[9]指出，眼電偽跡主要存在于低頻波段，可以通過高通濾波的方式予以去除；基于自回歸模型和自適應濾波的方法也已經應用于EOG偽跡去除，此類方法較為依賴參考通道的數據質量[10-11]。第二類是通過將腦電信號和偽跡信號分離或者分解到其他域進行去除[11]。此類方法較多，較為著名的是小波變換（wavelet transform，WT），將信號從時域映射到小波域，通過消除噪聲的小波系數消除噪聲[12]；使用最為廣泛的是獨立成分分析（independent component analysis，ICA），該方法主要原理是假設所有導聯的信號都來自于特定的源，觀測信號經過線性變換可以得到若干個源，通過去除EOG的源消除EOG偽跡；作為一種全數據驅動、高度自適應的非線性非平穩信號分析技術，經驗模態分解（empirical mode decomposition，EMD）和針對噪聲敏感問題改進的完備經驗模態分解（ensemble empirical mode decomposition，EEMD）也得到了廣泛應用[13]。
現階段更為主流的方法是將以上方法中的兩種或三種結合使用。樣本熵-完備經驗模態分解（sample entropy complete ensemble empirical mode decomposition with adaptive noise，SE-CEEMDAN）算法可以避免回歸算法中的雙向污染問題，但是同時引入了ICA、WT與樣本熵，大大增加了計算的復雜度[14]。基于ICA的算法中，去除包含偽跡的整個獨立成分（independent component，IC）的做法可能會造成某些混合在偽跡IC中的部分重要EEG信息的丟失[15-18]。Pion-Tonachini等[19]利用在線遞歸ICA等方法針對在線去除眼電偽跡做了改進。Castellanos等[18]提出的小波ICA方法（wavelet independent component analysis，wICA）可以減少信息損失，并且效果好過整個IC去除的方法，但是小波變換很大程度上取決于母波的選擇，而母波通常是包含正弦和余弦的函數，因此大多數時候不適用于非線性的EEG信號[7]。
針對上述問題，考慮到眼電偽跡主要存在于低頻波段[20]，為了更多地保留EEG信息，更好地適應非侵入式BCI，本文提出了一種結合ICA與頻率濾波器消除眼電偽跡的方法。首先使用ICA算法對獲取的原始數據進行分解以獲得不同的IC，之后采用相關系數與峰度雙重閾值判定的方式自動識別偽跡IC[21-22]，然后使用高通濾波器對偽跡IC進行濾波處理，最后通過逆ICA重構出觀測信號。
1 基于獨立成分分析的眼電偽跡自動去除算法
1.1 算法整體設計
本文提出的方法整體流程圖如圖1所示。主要處理步驟可分為三部分：首先對原始數據預處理，然后通過ICA算法分離并正確識別偽跡IC，最后通過帶通濾波器濾除偽跡成分并重構以獲得干凈的EEG。
圖1
算法流程圖
Figure1.
The flow chart of algorithm
1.2 獨立成分分析及眼電偽跡識別
ICA的提出最早是為了解決盲源分離問題，通常指可以根據高階統計量從線性混合中恢復統計獨立成分的一類算法，目前被認為是用于識別和去除腦電信號中偽跡魯棒性較強的算法，其數學表達式如式(1)所示：
|
其中，A稱為混合矩陣， 表示源信號， 表示觀測信號。
對經過預處理（0.5 Hz高通濾波，設置參考電極為M1、M2）的原始信號使用Infomax ICA進行分解得到若干IC，引入皮爾遜相關系數、z分數閾值和峰度閾值以正確識別偽跡IC。由于EOG是EEG強度的幾十到幾百倍，偽跡IC與EOG通道信號的相關性要遠高于非偽跡IC與EOG通道信號的相關性，因此可以使用皮爾遜相關性對二者進行區分。z分數使得閾值設置具有自適應性，降低了直接使用相關性閾值判定導致誤判的可能性。z分數閾值公式為：
|
其中，X為原始數據， 為平均數，S 為標準差。
通過以上計算，可以得到每個IC與EOG通道數據的相關性，并且通過z分數自動標記出候選偽跡IC。如圖2所示，橫軸表示IC數量，縱軸表示相關性；標紅部分為滿足閾值條件的IC，即候選偽跡IC。
圖2
ICs與垂直眼電（上）通道、水平眼電通道（下）數據的相關性得分
Figure2.
The correlation score of ICs with the data of vertical EOG (VEOG) (upper) and horizontal EOG (HEOG) (lower) channels
通過z分數閾值檢測可能仍然會對某些與EOG相似的IC產生誤判，信號的高階統計可以反映更多低階統計所不具備的信息。四階統計量峰度是具有簡單計算和理論的信號峰值的度量。通常，高斯信號的峰度為0，亞高斯信號峰度為負值，而EOG等超高斯信號峰度為正值，因此可以用峰度對候選偽跡IC做最終判定[22-24]。針對不同的數據集，峰度的閾值可以通過多次驗證確定，本文中峰度閾值為1.5[25]。若已知x是一個標量隨機變量，則峰度計算公式為：
|
其中， 表示變量x的n階中心距， 為隨機變量x的均值。
1.3 眼電偽跡去除
EEG失真的主要原因是受到了EOG的污染，通過ICA分解可以獲得含有EEG信息的偽跡IC。去除偽跡IC中的EOG成分即可保留其EEG信息，從而減少重要腦電信息的丟失。由于腦電信號的主要頻帶為0.1～40 Hz，眼電信號的主要頻帶在7 Hz以下[20]，濾除存在于偽跡IC中的眼電成分即可保留絕大部分腦電信息。
2 實驗數據及評價指標
2.1 實驗數據
2.1.1 公開數據集
選用2016年Klados和Bamidis公開發表的數據集[26]（https://data.mendeley.com/datasets/wb6yvr725d/1）。該數據采集實驗使用了19通道電極帽，電極排布符合國際10-20規則，以200 Hz的采樣頻率記錄了27名受試者的閉眼腦電數據，每段數據時長28 s。該數據集還提供了垂直眼電和水平眼電的數據，并且根據式(4)將眼電數據添加到了腦電信號中：
|
其中， 表示受到眼電偽跡污染的腦電信號， 表示無眼電偽跡的腦電信號， 表示垂直眼電信號， 表示水平眼電信號， 表示垂直眼電污染系數， 表示水平眼電污染系數，i表示受試者編號，j表示通道編號。
2.1.2 本實驗室數據
本實驗室數據來自河北省數字醫療工程重點實驗室，該實驗室位于河北大學新校區C座1號樓一層。使用40導便攜式腦電信號放大器（Grael，Compumedics Ltd.，澳大利亞）以及數據采集和在線處理軟件CURRY8 X（Compumedics Ltd.，澳大利亞）采集了12名（9名男生，3名女生）年齡19～25周歲在校大學生執行目標識別任務時的EEG數據。電極排布依據國際10-20標準，包含32個EEG導聯和2個EOG導聯（垂直EOG通道VEOG和水平EOG通道HEOG），數據采樣頻率為512 Hz。所有受試人員均由招募而來，視力正常或矯正視力正常，均為右利手，且均簽署了知情同意書。同時，該數據采集實驗通過了河北大學附屬醫院倫理委員會審查。
2.2 評價方法
有研究對EOG偽跡去除前后的數據進行了定量分析，主要評判內容為：腦電有效信息丟失程度與眼電偽跡去除效果[27]。根據數據的不同，本文對公開數據集和本驗室數據采用了不同的評價指標。
2.2.1 公開數據集評價
相關系數（correlation coefficient，CC）是在時域上判斷信號失真的一個重要指標，表征了兩個變量之間的相關程度。通過計算去除偽跡之后數據與純凈腦電之間的相關性，可以有效量化腦電信號丟失程度，相關系數越大，表明腦電信號保留越完整。均方根誤差（root mean square error，RMSE）表征了眼電偽跡去除效果，該值越小，表示眼電偽跡去除越完全。
| '/> |
| '/> |
其中， 表示EOG偽跡去除之后的EEG信號，x 表示純凈的EEG信號。
2.2.2 本實驗室數據評價
由于本實驗室數據無法獲得純凈的EEG，無法計算RSME和CC，因此采用幅度平方相干性（magnitude squared coherence，MSC）作為指標。MSC是定量評價偽跡去除前后頻域信息變化的常用度量[28-31]，通過相干性得分，可以有效評估在不同頻率上信息的保留程度。主要包含眼電信息的低頻段相干性得分越低，表明眼電信息去除越多，其余頻段相干性得分越高，表明腦電信息保留越完善。計算方法如式(7)所示：
|
其中， 是信號x與信號y在頻率 下的互譜密度， 和 分別指信號x和信號y在頻率 下的自譜密度。
3 結果與討論
3.1 公開數據集結果分析
本文提出的方法在公開數據集中的35組數據上進行了應用，每組數據都選取全部時長。圖3展示了公開數據集中數據1的FP1通道部分數據在EOG偽跡去除前后與純凈EEG的對比結果，藍線表示受到EOG偽跡污染的EEG，橙線表示經本文算法處理之后的EEG，綠線表示純凈EEG。如圖3所示，經過本文算法的處理，由EOG干擾產生的EEG波形畸變得到了很好的校正，說明該算法可以實現單個EOG偽跡的良好去除。
圖3
數據1 FP1通道純凈EEG、受污染的EEG和偽跡去除后EEG的對比
Figure3.
Comparison of the artifact-free, the contaminated and the cleaned EEG signals of dataset1, channel FP1
如圖4所示，由于在數據采集過程中受試者多次產生了眨眼活動，原始數據受到EOG干擾而多通道多處呈現出嚴重的波形畸變。經過本文算法的處理，所有通道的多處形變都得到了校正，說明本算法可適用于多通道多處EOG偽跡的自動化消除。
圖4
公開數據集眼電偽跡去除效果
Figure4.
The effect of EOG artifact removal in open datasets
為了評估該算法的有效性和穩定性，眼電偽跡去除之后，對每組數據都計算了相關系數與均方根誤差，并與SE-CEEMDAN、小波增強ICA進行了對比。如表1所示，受到眼電偽跡污染的EEG信號與純凈腦電信號之間的均方根誤差很大、相關系數很小，通過三種算法處理，實現了不同程度的改善。從結果看，本文算法平均RSME優于SE-CEEMDAN和小波增強ICA，即本文算法能夠去除更多的偽跡信息；平均CC明顯高于SE-CEEMDAN，即本文算法能保留更多的EEG信息。該結果表明，本文算法處理后的數據更加接近純凈EEG。
表1
不同算法性能對比（均值±標準差）
Table1.
Performance comparison of different methods (mean±SD)
3.2 本實驗室數據分析
本實驗室數據可以通過分析偽跡去除前后的MSC來定量分析EEG成分的保留程度。另外，由于EOG偽跡本身存在高振幅、低頻率的特點，可以通過觀察波形圖定性分析EOG偽跡去除效果。圖5展示了本文方法處理前后的波形圖。從圖中可以看出，采集到的原始數據由于受到EOG的影響而多處產生嚴重形變。經過本文算法處理后，在13.5 s附近和19 s附近嚴重的形變消失，21 s附近較弱的形變也得到了很好的校正，而其他沒有EOG偽跡的部分未受到影響。說明本文的算法對于EOG偽跡具有較高的識別度，并且能夠很好地去除不同幅度的EOG偽跡。
圖5
本實驗室數據眼電偽跡去除效果
Figure5.
The effect of removal of EOG artifacts from our laboratory data
人工檢查發現本實驗室獲取的12組數據均含有EOG偽跡，對每名受試者的數據都進行了EOG偽跡去除處理。由于眼睛距離額葉最近，因此額葉通道受到EOG影響更大，此處計算了EOG偽跡去除前后額葉通道FP1的MSC。由于MSC是一個表征信號在不同頻率上差異的量，而EOG偽跡主要存在于低頻波段，將其去除后與原始信號在低頻段差異較大，因此低頻波段MSC得分應該較低；相反，在高頻波段MSC得分應該較高。如圖6所示，低頻部分本文算法MSC得分較小波增強ICA更低，表明本文算法去除了更多的眼電偽跡，而高頻部分本文算法MSC得分較小波增強ICA更高，表明本文算法對高頻部分EEG信息的影響更小。
圖6
本實驗室數據經不同算法處理后的MSC結果
Figure6.
MSC results of the laboratory data processed by different algorithms
4 結論
基于非侵入式BCI快速發展卻容易受到EOG偽跡影響的現狀，針對現有眼電偽跡去除算法存在EEG信息丟失的問題，本著自動、高效的原則，本文提出了一種基于ICA的眼電偽跡自動去除算法。該算法通過z分數與峰度雙重閾值自動判斷偽跡IC，采用頻率濾波器去除EOG偽跡，實現了EOG偽跡的自動選擇和去除，提高了數據處理效率。在公開數據集和本實驗室數據上分別進行了測試，公開數據集結果顯示，相比于SE-CEEMDAN和小波增強ICA算法，本文提出的方法具有更小的RSME和更大的CC，表明該算法能夠更加有效地去除眼電偽跡，同時能夠保留更多的EEG信息；本實驗室數據分析結果表明，本文所提出的算法對不同幅度的EOG偽跡都實現了良好的去除，又幾乎不影響高頻EEG信息。綜上所述，該算法實現了EOG偽跡更大程度的自動去除，又大大減少了EEG信息的丟失，可能為非侵入式BCI的推廣起到一定助力作用。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：宋昊主要負責數據收集與分析、算法程序設計以及論文編寫；徐頌、劉國明主要負責實驗素材收集；劉京主要負責實驗指導；熊鵬主要負責實驗設計、數據分析指導以及論文審閱修訂。
倫理聲明：本研究通過了河北大學附屬醫院倫理委員會的審批（批文編號：HDFY-LL-2019-035）。
