基于局域均值分解与典型相关分析的眼电伪迹去除方法

李明爱，田晓霞，孙炎珺，杨金福

基于局域均值分解与典型相关分析的眼电伪迹去除方法

李明爱，田晓霞，孙炎珺，杨金福

（北京工业大学电子信息与控制工程学院，北京100124）

为消除眼电伪迹（ocular artifact，OA）对脑电信号（electroencephalography，EEG）造成的严重影响，提出一种基于局域均值分解法（local mean decomposition，LMD）与典型相关分析（canonical correlation analysis，CCA）的眼电伪迹自动去除方法，记为LMDC法.首先，利用LMD将每导脑电采集信号自适应地分解为一系列具有物理意义的乘积函数（production function，PF）分量，通过CCA去除PF分量之间的相关性，获得相应的典型变量；其次，计算每导脑电信号与多导眼电信号间的相关系数矩阵，实现眼迹成分的自动识别，将典型相关变量中对应眼迹成分的部分随机变量置零，其余随机变量不变，得到新的典型相关变量；最后，基于CCA逆变换将新的典型相关变量投影返回得到眼迹去除后的PF分量，并进一步重构出眼迹去除后的脑电信号.基于BCI竞赛数据库进行实验研究，结果表明:LMDC法相对其他常用方法获得了较好的眼迹去除效果，并对多位实验者和多种眼迹表现出较强的自适应性.

脑电信号；眼电伪迹；局域均值分解；典型相关分析；自适应

脑电信号（electroencephalography，EEG）是大脑皮层脑神经细胞电活动的总体反映，是人体大脑发出的自发性、节律性的重要生理信号，包含大量的生理和病理信息，在脑功能研究、神经系统疾病诊断及功能康复等多个方面发挥着重要的作用.然而，在EEG的采集过程中，不可避免地因眼睛或眼球的运动而产生眼电信号（electrooculography，EOG），从而导致 EEG发生严重畸变，形成眼电伪迹（ocular artifact，OA）.OA幅值较大，且频率范围与EEG频率有一定重叠，是EEG信号的最主要干扰之一，在很大程度上影响对脑电信号的进一步分析和应用［1-3］.

目前，去除眼电伪迹常用的方法主要包括以下几种:1）小波变换法（wavelet transform，WT）［4-5］.该方法是一种典型的时频分析方法，能在不同尺度下观察EEG非平稳信号的局部特征.但该方法受小波基函数和分解层数的主观选择的影响，缺乏自适应性；2）独立成分分析法 （independent component analysis，ICA）［6-7］.该方法基于非高斯性最大化原理，能有效去除信号的二阶及高阶统计相关性，但求解分离矩阵的计算量大，甚至会出现不收敛，实时性欠佳；3）典型相关分析法（canonical correlation analysis，CCA）［8-9］.该方法是一种研究2组变量间相关性的多元统计学方法，通过将多导EEG信号典型相关变量中与EOG干扰相关的部分分量置零，消除EOG对各导EEG信号的影响.该方法简单，计算速度快，适合EEG的在线处理，但去除效果有待提高.

为进一步提高眼迹去除方法的精度及对实验者、眼迹类型的自适应性，本文将局域均值分解法（local mean decomposition，LMD）和CCA相结合，提出一种自动去除脑电信号中眼电伪迹的方法，记为LMDC.利用脑机接口（brain computer interface，BCI）竞赛公共数据库提供的多位实验者的多种眼电数据进行实验研究，证明LMDC方法具有较强的去噪能力和自适应能力.

1　基本原理

1.1局域均值分解法

局部均值分解法是由Smith［10］于2005年提出的一种新的适合于非线性、非平稳信号的自适应时频分析方法.LMD根据复杂信号自身的特征自适应地将其从高频到低频逐级分解为一系列具有瞬时物理意义的乘积函数（production function，PF），每个PF分量都是一个纯调频信号和一个包络信号的乘积.目前，LMD已用于地震、脑电等信号分析及机械故障诊断等方面［11-12］.

给定信号x（t），对其进行局域均值分解的具体过程如下:

1）确定信号x（t）的所有局域极值点ai（i=1，2，…，M），M为极值点数.

2）依计算相邻局域极值点的平均值mi（1=1，2，…，M-1），将所有的相邻均值用折线连接并利用滑动平均法进行平滑处理，得到局域均值函数m11（t）.

3）依计算相邻局域均值点的包络估计值bi（i=1，2，…，M-1），将所有的相邻包络估计值用折线连接并利用滑动平均法进行平滑处理，得到包络估计函数n11（t）.

4）将平滑后的局域均值函数m11（t）从原信号x（t）中分离出来，得到信号

5）用h11（t）除以包络估计函数n11（t）得到调频信号

如果满足-1≤s11（t）≤1，并且它的包络估计值n12（t）=1，则s11（t）是一个纯调频信号，停止迭代；否则，将s11（t）作为新信号重复步骤1）～5），直到s1n（t）为纯调频信号.迭代过程为

式中

6）将上述每次迭代过程中产生的包络估计函数相乘得到第一个分量的包络信号

7）将纯调频信号s1n（t）乘以包络信号n1（t），得到原信号x（t）的第一个PF分量

Fp，1（t）包含了信号x（t）的最高频率成分.

8）将Fp，1（t）分量从信号x（t）中分离出来，获得新的原始信号y1（t），重复步骤1）～7），循环k次，直到yk（t）为单调函数.

综上所得，原信号x（t）分解为k个PF分量和一个单调函数yk（t），即

1.2典型相关分析

典型相关分析是由Hoteling于1936年提出的、研究2组随机变量间相互关系的一种多元统计分析方法.其目的是寻找2组随机变量各自的线性组合，使得线性组合后的2个随机变量的相关性达到最大.

假设2组随机变量分别为:X=［x1，x2，…，xl］T∈RRl×N和Y=［y1，y2，…，yl］T∈ RRl×N，l表示每组随机变量的个数，N为每个随机变量的样本点数.经过样本中心化处理后，2组随机变量分别记为:和对随机变量和分别进行线性组合，即

式中:u和v为典型相关变量；基向量矩阵Wx=［wx1，wx2，…，wxl］∈ RRl×l和Wy=［wy1，wy2，…，wyl］∈RRl×l应使u和v的相关系数ρ（u，v）=最大.

根据拉格朗日乘法，将相关系数ρ（u，v）的最大值问题转化为求的特征值问题.式中:是矩阵Cyx和的特征值，Wx和Wy则按特征值λ由大到小的顺序将其对应的特征向量排列而成.

2　基于LMDC方法去除眼电伪迹

近年来，许多学者将局域均值分解用于处理非线性、非平稳信号，同时，典型相关分析在信号的盲源分离中也发挥了很大的作用.因此，本文将局域均值分解和典型相关分析相结合（称为LMDC法），自动识别和去除脑电信号中的眼电伪迹.该方法能充分发挥两者的优点，更精细且尽可能多地去除眼迹的同时最大程度保证脑电信号的完整性.LMDC法具体过程如下:

1）基于局域均值分解对每导脑电信号进行自适应分解，得到相应的PF分量矩阵.假设第i导脑电采集信号为si∈R RN×1，i=1，2，…，p，p为脑电极个数，N为采样点数.利用LMD将脑电信号si自适应地从高频至低频逐级分解，得到k个PF分量，即Fp，ij∈R RN×1（j=1，2，…，k），并将其构成si的PF分量矩阵，记为 Fp，1=［Fp，i1，Fp，i2，…，Fp，ik］T∈R Rk×N，i=1，2，…，p.

2）对各导脑电信号的PF分量矩阵进行典型相关分析.

首先，对第i导脑电信号的PF分量矩阵Fp，i（i=1，2，…，p）进行去中心化处理，重新记为=；然后，利用CCA对进行典型相关分析，求得的典型相关变量2，…，p），式中Wx∈R Rk×k为基向量矩阵.从而，达到减少PF分量之间相关性的目的.

3）利用相关系数自动识别和去除眼电伪迹成分.

假设采集的参考眼电信号为Eo=［g1，g2，…，gq］T∈R Rq×N，q为眼电信号电极个数，N为采样点数.依据式（13）计算典型相关变量ui=［ui1，ui2，…uik］T（i=1，2，…，p）与参考眼电信号Eo的相关系数［13］，即

4）采用CCA逆变换将眼迹去除后典型变量进行投影变换.

通过CCA逆变换，将典型变量 ～ui投影返回到去除眼迹后的PF分量矩阵

5）将眼迹去除后的PF分量相加重构，获得去除眼电伪迹的脑电信号.

利用上述眼电伪迹去除后的PF分量矩阵～Fp，i，重构得到去除眼电伪迹后的第i导脑电信号～si，即

3　实验研究

基于公用的BCI竞赛数据库提供的脑电和眼电数据，将LMDC与其他常用眼迹去除方法进行实验对比研究，以验证本文方法的有效性.

3.1数据来源

实验数据来源于“BCI Competition 2008”竞赛数据库“two motor”，包括9位实验者在不同时段各进行5次采集实验所得45组实验数据.其中，实验者1的第2次实验及实验者5的第4次实验因采集设备故障未采集到有效数据.每次采集实验的时序如图1所示.首先，保持2 min睁眼和1 min闭眼状态；其次，1 min的眼球运动，包括眼球旋转、眨眼、眼球水平移动和眼球垂直移动各15 s；最后，想象左右手运动各60次.上述实验数据，眼电信号通过左眼左侧、眉心和右眼右侧3个位置获得；脑电信号则采用Ag/AgCl电极从国际标准10～20导联系统的C3、Cz和C4三个电极采集，采样频率为250 Hz，经过0.5～100 Hz的带通滤波和50 Hz的陷波滤波以去除工频干扰.

3.2性能评价指标

本文采用信噪比（signal noise rate，SNR）［14］和均方根误差（root mean squared error，RMSE）［15］2项指标评价眼电伪迹的去除效果，其计算式分别为

式中:c（n）为某电极上纯净的EEG；x（n）为该电极上含眼电伪迹的EEG进行眼电伪迹去除后的EEG；N为样本点数.SNR越大，或RMSE越小，说明去除眼迹后的脑电信号与纯净的脑电信号越接近，即眼电伪迹去除效果越好.

3.3LMDC去除眼迹的实验结果及分析

3.3.1基于实验者2眨眼眼电数据的眼迹去除实验

本部分将以实验者2的第1次实验采集的C3、Cz和C4三导联运动想象EEG信号为纯净的EEG信号，以左眼左侧、眉心和右眼右侧3个位置获得的5次三导联眨眼眼电数据为纯净的EOG信号，开展LMDC的眼迹去除实验研究.

1）含噪EEG和EOG信号的构造

为便于衡量LMDC去除眼迹的效果与性能，根据EEG和EOG的双向激活扩散性（bidirectionality），建立二者相互干扰的数学模型以获得真实的或含噪的脑电信号和眼电信号［16］.式中:Ee，r和Eo，r分别表示实际受到相互干扰的脑电和眼电信号；Ee，c和Eo，c分别表示纯净的脑电和眼电信号；随机矩阵对A∈［0，0.2］和B∈［0，1］分别为眼电对脑电信号、脑电对眼电信号的影响因子矩阵.为尽量减少随机矩阵取值对实验结果的影响，随机矩阵对A和B将随机产生10对，即共构造50组含眨眼EOG干扰的EEG数据及含EEG干扰的眨眼EOG数据.

2）部分中间实验结果

图2给出了利用LMD将某次含噪的C3、Cz和C4三导联EEG进行自适应分解所得到的PF分量及进一步利用CCA得到的去相关后的典型变量.其中，PFi表示第i个PF分量；Ui表示第i个典型相关变量.图3则展示了C3、Cz和C4三导联纯净的EEG、依式（17）构造的含眨眼眼电噪声的EEG及基于LMDC法去除眨眼眼迹后的EEG.可见，LMDC对实验者2的眨眼眼电伪迹取得了较理想的去除效果.

3）基于SNR和RMSE的实验结果

基于50组含噪EEG和眨眼EOG数据，将本文LMDC方法与其他常用眼迹去除方法进行比对实验研究，并用SNR和RMSE指标评价各方法的眼迹去除效果.图4给出了C3、Cz和C4三导联EEG进行眼迹去除时，SNR和RMSE指标的平均实验结果.可见，LMDC方法相比DWT、ICA及CCA几种相关伪迹去除传统方法呈现较明显的优势，在C3、Cz及C4导联上均取得了最小均方误差和最大信噪比.

3.3.2基于9位实验者4类眼电数据的眼迹去除实验

将9位实验者的眨眼眼电、水平眼电、垂直眼电和旋转眼电共4种、43组眼电数据作为纯净的EOG数据，仍以相同的竞赛数据库提供的运动想象脑电信号为纯净的EEG信号进行实验研究，多种方法的平均实验结果如表1所示.可见，LMDC方法仅在C3导去除旋转OA时的均方误差指标不及CCA，而其余情况下，针对眨眼眼电、水平眼电、垂直眼电和旋转眼电4种不同OA，相比其他方法SNR均有明显提高，RMSE均有显著下降，充分说明本文方法对多位实验者及多种OA表现出较强的自适应性.

表1　基于9位实验者4种OA数据的43次伪迹去除实验的平均结果Table 1　Average results of 43 times OA removal experiments with 4 kinds of OA sampled by nine subjects

4　结论

1）提出一种将局域均值分解和典型相关分析相结合的眼电伪迹去除方法，即LMDC法，能够实现对眼电伪迹的自动识别与去除.

2）基于国际BCI竞赛公用数据库，并根据眼电信号与脑电信号之间激活扩散的双向性构造混有眼电伪迹的脑电实验数据，对LMDC法的去噪能力进行了定量评估.相对其他传统方法，LMDC在信噪比和均方根误差2项性能指标下具有明显优势.

3）对9位实验者的4种眼电伪迹进行实验研究，LMDC法表现出较强的自适应能力.这将有利于脑电信号的进一步分析和实际应用.

［1］罗志增，曹铭.基于最大信噪比盲源分离的脑电信号伪迹滤波算法［J］.电子学报，2011，12（39）:2926-2931.

LUO Z Z，CAO M.An algorithm to filter artifacts in EEG based on blind source separation of maximum signal noise ration［J］.Acta Electronica Sin，2011，12（39）:2926-2931.（in Chinese）

［2］GE S，HAN M，HONG X.Erratum to a fully automatic ocular artifact removal from EEG based on fourth-order tensor method［J］.Biomedical Engineering Letters，2014，4（2）:195.

［3］LI M A，CUI Y，YANG J F.Automatic removal of ocular artifact from EEG with DWT and ICA method［J］.Applied Mathematics&Information Sciences，2013，7（2）:809-816.

［4］YUAN F L，LUO Z Z.The EEG De-noising research based on wavelet and Hilbert transform method［C］∥InternationalConferenceonComputerScienceand Electronics Engineering（ICCSEE）.Piscataway，NJ: IEEE，2012:361-365.

［5］徐伟业.一种快速离散小波变换算法及其在语音信号中的应用［J］.计算机工程与应用，2011，47（35）: 143-146.

XU WY.Fastdiscretewavelettransformandits application to speech signal processing［J］.Computer Engineering and Applications，2011，47（35）:143-146.（in Chinese）

［6］董洁，王涛，张爱桃.基于独立分量分析去除脑电中眨眼和水平扫视的伪迹［J］.航天医学与医学工程，2011，24（2）:122-127.

DONG J，WANG T，ZHANG A T.Removal of blind and saccade artifactinEEGrecordingswithindependent component analysis［J］.Space Medicine&Medical Engineering，2011，24（2）:122-127.（in Chinese）

［7］SOOMRO M H，BADRUDDIN N，YUSOFF M Z，et al.A method for automatic removal of eye blink artifacts from EEG based on EMD-ICA［C］∥ Colloquium on Signal Processing and Its Applications（CSPA），2013 IEEE 9th International.Piscataway，NJ:IEEE，2013:129-134.

［8］张莉，何传红，何为.典型相关分析去除脑电信号中眼电伪迹的研究［J］.计算机工程与应用，2009，45（31）:218-220.

ZHANG L，HE C H，HE W.Research of removing EOG artifacts fromEEGbasedonCCA［J］.Computer Engineering and Applications，2009，45（31）:218-220.（in Chinese）

［9］JIN X，QIU T，LIU W H.An ocular artifacts removal method based on canonical correlation analysis and twochannelEEGrecordings［J］∥ ChineseJournalof Biomedical Engneering（English Edition），2012，21（2）: 60-66.

［10］SMITH J S.The local mean decomposition and its application to EEG perception data［J］.J R Soc Interface，2005（2）:443-454.

［11］朱晓军，樊刘娟，吕士钦，等.LMD方法在脑电信号处理中的应用研究［J］.计算机科学，2012，39（2）: 273-275.

ZHU X J，FAN L J，LÜ S Q，et al.Application research of LMD method in EEG signal processing［J］.Computer Science，2012，39（2）:273-275.（in Chinese）

［12］杨斌，程军圣.基于LMD和主分量分析的齿轮损伤识别方法［J］.振动、测试与诊断，2013，33（5）:809-813.

YANG B，CHENG J S.Fault detection method for gear condition based on LMD method and principle component analysis［J］.Journal of Vibration，Measurement& Diagnosis，2013，33（5）:809-813.（in Chinese）

［13］赵欣欣，孔万增.基于ICA和相关性分析的脑电信号眨眼伪迹去除［J］.杭州电子科技大学学报，2013，33（1）:32-35.

ZHAO X X，KONG W Z.Removal of blink artifact in EEG recording based on FastICA and correlation analysis［J］.Journal of Hangzhou Dianzi University，2013，33（1）:32-35.（in Chinese）

［14］ZHAO C，QIU T.An automatic ocular artifacts removal method based on wavelet-enhanced canonical correlation analysis［C］∥Engineering in Medicine and Biology Society（EMBC），2011 Annual International Conference of the IEEE.Piscataway，NJ:IEEE，2011:4191-4194.

［15］王文波，张晓东，汪祥莉.基于独立成分分析和经验模态分解的混沌信号降噪［J］.物理学报，2013，62（5）:27-34.

WANG W B，ZHANG X D，WANG X L.Chaotic signal denoisingmethodbasedonindependentcomponent analysis and empirical mode decomposition［J］.Acta Phys Sin，2013，62（5）:27-34.（in Chinese）

［16］李明爱，梅意城，孙炎珺，等.眼电伪迹自动去除方法的研究与分析［J］.仪器仪表学报，2014，35（11）: 2515-2523.

LI M A，MEI Y C，SUN Y J，et al.Research and analysis of ocular artifact automatic removal［J］.Chinese Journal ofScientificInstrument，2014，35（11）: 2515-2523.（in chinese）

（责任编辑吕小红）

Ocular Artifact Removal Based on Local Mean Decomposition and Canonical Correlation Analysis

LI Ming蒺ai，TIAN Xiaoxia，SUN Yanjun，YANG Jinfu
（College of Electronic Information and Control Engineering，Beijing University of Technology，Beijing 100124，China）

Based on local mean decomposition（LMD）and canonical correlation analysis（CCA），an automatic removal method，denoted as LMDC，was proposed to eliminate the serious impact of ocular artifact（OA）from electroencephalography（EEG）.Each recorded EEG was decomposed into a series of physically meaningful production function（PF）components adaptively by LMD，and CCA was applied to eliminate the correlation among the PFs to get the corresponding canonical correlation variable.Then，the correlation coefficient matrix between each EEG and multi electrooculogram（EOG）was computed to recognize the OA component automatically.The random variables corresponding with OA components in the canonical correlation variable were set to zero，and the others remain unchanged to obtain a new canonical correlation variable.Finally，the inverse algorithm of CCA was utilized to project the new canonical correlation variable to the OA free PFs，and the OA removed EEG was reconstructed. Experimental research was conducted on a public brain computer interface（BCI）completion database. Experiment results show that LMDC has better performance than that of the other related methods，and has stronger adaptability for multi subjects and types of OA.

electroencephalography；ocular artifact；local mean decomposition；canonical correlation analysis；adaptivity

R 318

A

0254-0037（2016）06-0843-08

10.11936/bjutxb2015050009

2015-05-05

国家自然科学基金资助项目（81471770）；北京市自然科学基金资助项目（7132021）

李明爱（1966—），女，教授，主要从事脑机接口技术与智能控制方面的研究，E-mail:limingai@bjut.edu.cn