专利名称：一种无需设定阈值的脑电信号中眼电伪迹自动去除的方法
技术领域：
本发明涉及生物信息技术领域，特别涉及脑电信号(Electooenc-ephalogram,EEG)的预处理技术。具体是涉及一种无需人工设定阈值的脑电信号中眼电伪迹(Ocular Artifact, 0A)的自动去除技术。
背景技术：
脑电信号是一种微弱的电生理信号，其中包含了大量的生理与疾病信息，在临床医学、认知心理学等研究领域中发挥着重要作用。头皮表面记录脑电的方式是目前脑电信号研究中所采用的最主要的信号采集方式，该方式记录技术简单、易操作，但信号易受噪声的干扰。脑电干扰给脑电信号的解释与分析造成了很大的困难。因此，有效识别并去除脑电信号中的伪迹，得到纯净的脑电信号，对进一步研究大脑的真实活动具有重要意义。影响脑电信号的伪迹主要包括眼电伪迹、心电伪迹、肌电伪迹、50Hz工频干扰等。相对于其它伪迹成分而言，眼电伪迹是使EEG信号失真的一种最主要的干扰成分。眼电伪迹的幅值较大，最大可达到IOOmV,并且眼电信号源距离测量电极较近，对脑电信号造成了很大的干扰。眼电的产生是由于跨膜间存在静息电位，角膜侧为正，巩膜侧为负，这可以看作是有一个电偶极子存在。眼球转动过程中，该偶极子改变了眼睛周围的电场，从而影响头皮电场，眼球运动方向一侧的电极记录到眼电伪迹为正，眼球运动相反方向一侧的电极记录得到的眼电伪迹为负，由于大脑是电的良导体，眼电信号从前额处向后传播，遍历整个头部，放置在头部用于采集脑电信号的电极可以采集到这些信号，造成对脑电信号的干扰。而且，眼电伪迹的频带与脑电信号的频带相混叠，很难使用传统的滤波方法进行去除，因此，脑电信号预处理阶段最重要的工作就是如何有效地去除眼电伪迹。从脑电信号中有效地去除眼电伪迹一直是生物医学信号处理领域的研究热点，其中一些主要的研究成果有I实验控制法这种方法要求研究人员给予受试者一定的实验指导，例如要求受试者在实验过程中始终保持闭眼状态或尽量避免眨眼和眼球转动以减少眼电伪迹的干扰，这样做虽然可以在一定程度上抑制干扰，但是前者会引入α波的干扰，后者在长时间的实验过程中极易引起疲劳，并且不易做到。2伪迹剔除法它的主要思想是将受到眼电伪迹影响脑电信号剔除以获得较纯净的脑电信号，但是这样做的同时剔除了大量的脑电信号，造成了信息的遗漏，并且人工剔除干扰信号的劳动量过大，主观性强，而且不适应在线脑机接口系统的需求。3回归方法研究表明眼电伪迹在测量电极处的衰减程度与电极到眼部的距离有关，此方法利用这一原理通过计算得到眼电信号在每个头皮电极处的衰减因子，然后将采集到的每导脑电信号减去乘以衰减因子的眼电信号，从而得到伪迹去除后的数据。它的模型为IiliC,{ε) = !]‘:(;、(ε) - ηΕΟΟ ( I )
式中，EOG表不测量得到的眼电信号；EEGS( ε )表不原始脑电信号；EEG( ε )表不伪迹去除后的脑电信号；ε表示相应的测量通道；η表示比例系数，即眼电信号在第ε导处的衰减因子。H值常用回归分析的方法进行估计。这种伪迹去除方法的优点是直观易懂，物理意义明确，缺点是没有考虑到脑电信号的双向传播性，即实际测得的眼电信号中也包含一些脑电成分。因此，回归方法去除伪迹的同时不可避免的去除某些脑电成分。4小波变换法小波变换(Wavelet Transform, WT)是数学领域的一个重要分支，近年来越来越受到各个学科领域的极大重视。由于其在时频两域都具有表征信号局部特征的能力和多分辨率分析的特点，因此被誉为“数学显微镜”。小波变换的基本思想是将原始信号通过伸缩和平移后，分解为一系列具有不同时间分辨率、不同频率特性的子带信号，这些子带信号具有良好的时域、频域等局部特征，进而可以实现对信号在时域和频域的局部化分析，从而克服了傅里叶分析在处理非平稳信号和复杂图像时所存在的局限性。但是使用小波变换法去除伪迹要求信号和噪声的频带不能混叠。而脑电信号与眼电伪迹的频带相混叠，于是研究者开始将传统的小波变换法与其它去除伪迹的方法相结合，取长补短，期望获得更好的效果。
5主成分分析法(Principal Components Analysis, PCA):主成分分析是一种多通道的信号处理方法，其基本思想是利用正交原理将一组相关变量转化成为另一组相互独立的变量(即主成分)，分解出来的主成分按大小排序，通常前面几个主成分就包括了信号80%的信息，因此PCA也被用来数据降维，去除冗余。PCA可以把多导联的脑电信号分解为相互正交的成分，去掉噪声和伪迹成分，再对其余成分进行重构，就可以得到较为纯净的脑电信号。PCA的主要缺点在于I信号和噪音波形相似时，难以从EEG中去除伪迹；2它只对正交的信号源比较有效，在去除伪迹时只能去除与脑电信号相正交的干扰成分；3它只适合于分解低阶相关信号，难以分解高阶的相关信号。6独立分量分析(Independent Component Analysis, I CA)独立分量分析理论的发展可追溯到上个世纪80年代初期，法国学者J. Hernah和C. Juttne等人首次提出了ICA分析的基本概念，最初被用来解决“鸡尾酒会”问题，随着近年来在ICA方面研究兴趣的增加，该理论已被广泛应用在图像处理、天线阵列、控制科学、生物医学信号处理等领域，其研究热潮方兴未艾。ICA是实现盲源分离(Blind Source Separation, BSS)最主要的方法之一,是一种基于高阶统计量的分析方法，处理对象为非高斯信号。理论上可以认为脑电信号、眼电伪迹、心电伪迹、肌电伪迹以及其他干扰源所产生的干扰信号都是由相互独立的源产生的，而且可以将观测信号视为是伪迹与脑电信号的线性组合，并且伪迹源的数目通常比头皮上测量的脑电通道数要少，因而独立成分分析方法适用于分离脑电信号中的各种伪迹。在实际应用中，ICA方法可以很好的分离频带相重叠的信号，弥补了小波分析去除伪迹的不足。但是ICA方法也有缺陷在分离出信号和伪迹后，需要人工识别并筛选出伪迹成分再进行去除，费时费力，且具有很大的主观性。
发明内容
本发明针对传统脑电信号中眼电伪迹去除的过程中需要人工对眼电伪迹进行识另|J、费时费力、工作量大的缺陷，提出了一种眼电伪迹自动去除的方法。其主要特征包括首先利用独立分量分析对含有眼电伪迹的脑电信号进行分解，得到若干个独立分量，然后求取每个独立分量的峭度、序列renyi熵和样本熵作为区分脑电信号和眼电伪迹成分的特征，进而采用k均值聚类分析——一种无监督的分类算法进行特征分类，通过将含有眼电伪迹的独立分量置零重构出纯净脑电信号。具体步骤如下I独立分量分解本发明采用1999年芬兰学者Aapo Hyvarinen等人提出的基于负熵的固定点ICA算法，简称fastica算法，进行独立分量分解。设η通道观测信号x (t) = [Xl (t)，X2 (t)，…，xn(t)],经过 ICA 分解得到 m 个独立分量 y (O = Iiy1 (t)，y2 (t),…，ym(t)]。 2特征提取经过ICA分解得到m个独立分量y (t) = Hy1 (t)，Yi (t)，…，ym(t)],分别求取它们的峭度、序列renyi熵和样本熵作为分类特征。2. I峭度峭度是时间序列的四阶累计量，可以表征信号的波动情况。眼电伪迹的峭度值要显著大于脑电信号，因此，利用该参数作为识别眼电伪迹的特征。2. 2 renyi熵renyi熵是一种反映信号随机性的参数,眼电伪迹的renyi熵值明显小于脑电成分，因而使用其作为区分眼电伪迹和脑电成分的特征。2.3样本熵样本熵是一种测量时间序列复杂度的方法，而且眼电伪迹的样本熵值比脑电成分低，因此使用样本熵作为特征区分眼电伪迹和脑电成分。3特征分类本发明使用k均值聚类分析的方法进行分类，这是一种动态的、无监督的分类方法。区分眼电伪迹和脑电成分显然是一个两分类问题，因此先随机取两个样本作为聚类中心，计算每个样本到这两个中心的距离，把样本归入最近的聚类中心，然后用本类所有样本的均值代替原有的聚类中心的值，这样就得到了新的聚类中心，再重新计算各样本到新的聚类中心的距离，重新分类、修改中心点，直到新的聚类中心与上一次的中心点一致时，算法停止，然后再按照最小距离原则将独立分量分为两类，算法流程如图2所示。4信号的重构将含有眼电伪迹的独立分量置零，其它独立分量不变，对信号进行重构，得到纯净的脑电成分。具体如下一种无需人工设定阈值的脑电信号中眼电伪迹自动去除的方法，其特征在于首先，对采集的脑电信号进行独立分量分解，得到m个独立分量；其次，对于每一个独立分量，求取它们的峭度、序列renyi熵和样本熵作为特征向量；然后，使用k均值聚类分析——一种无监督的分类算法自动识别出含有眼电伪迹的独立分量；最后，将含有眼电伪迹的独立分量置零，其余独立分量不变，对信号进行重构，得到纯净的脑电信号；包括以下步骤I)独立分量分解使用fastICA算法对信号进行分解，它的基本模型表达式为x(/) = As[l)(I)式中，χα)=[χια)，χ2α)，···，χηα)]τ e RnxN是η通道的观测向量，为测量得到，I彡t彡N，N表示采样点总数，S (t) e rXN为源信号，混合矩阵RnXm，独立分量分析的目标是求解分离矩阵W，通过w从观测信号x(t)中恢复出未知的源信号s (t)，源信号s (t)使用下式估计得到⑵式中，y⑴是源信号s⑴的估计，y (t) = [Yl (t)，Yi⑴，…，ym (t) ]T e RmXN是求得的m个独立分量，其中Ii (t)是第i个独立分量，I^i 为整数，t表示第t个采样点，分离矩阵w e Rmxn ；2)特征提取2. I峭度对于每一个独立分量Ji (t)，求取它的峭度；k = m4-3m( 3)
式中，k为独立分量yi(t)的峭度值，I彡i彡m，i为整数；me为yi⑴的β阶中心矩，Hie=EKyi (O-Iii1) 0}, β ={2，4}, Ε( ·)表示对括号里的表达式求均值，Hi1为Yi (t)的均值；2. 2 renyi熵对于每一个独立分量Ji (t),其计算公式为H (y；{>)) = logΣΣ('( ~^ (Λ))(4)
y—oc iv /=1 h=i式中，G(g.2cr2)=^-L.^ e—F , Θ =Yi(J-)-Yi (h);方差 σ =0. 2 0. 3, α =1. 8, N 为独
立分量Ii (t)的采样点个数，yi (j)与Ii (h)是独立分量Ii⑴在采样点j和h的采样值；2. 3样本熵对于每一个独立分量Ji (t)，其样本熵的计算过程如下第一步设独立分量Yi⑴有N个采样点，t表示第t个采样点，I彡t彡N，t为整数；第二步从独立分量Ji (t)中提取b维矢量，b=2,即
_8] Yb{t) = VyMyl^m(5)式中，yi(t)和7七+1)是独立分量Yi⑴在采样点t和t+Ι的采样值，I彡t彡N-I，t为整数；第三步使用第二步的方法得到两个b维矢量为Yb (t)和YJt1)，定义它们之间的距离为d[γ ( = max[|χ· {t + c)~y, {ti + c)|](6)式中，I彡 Lt1SN-I, t 关整数，Yb(t)和 YbU1)为公式(5)中 Yb(t)在采样点t和的值％ (t+c)和Yi (ti+c)为独立分量Yi (t)在采样点t+c和h+c的采样值，C e {O, 1}，求取差值的最大值，即为两个b维矢量的距离；第四步设定一个阈值r, r= (O. I "O. 25) StcKyi (t)), StcKyi (t))为 Yi (t)的标准差；对于每一个t值,统计d[Yb (t), YbU1)]小于r的数目，记为T, I ^ t, ^ N-1, t幸t1；t、&为整数，&可以取N-2个值，因此，总的距离数为N-2，计算T与总的距离数N-2的比
值，记为O·)，即
^ (O = ^η)第五步计算B1; (/·)的均值，即^(/-) = 77-7^^00C8)
—丨，1第六步再从同一个独立分量Ji⑴中提取Id1维矢量，b1=3,即V、(/) = ^,(0^,(/ + 1)^,(/ + 2)] (9)式中，yi (t)、Yi (t+1)和yi (t+2)是独立分量Ii (t)在采 样点t、t+1和t+2的采样值，I彡t彡N-2，t为整数；第七步使用第六步的方法得到两个Id1维矢量为&的和1；1@1)，定义它们之间的距离为d[7&i (t) Jbi^)] = max^[|>>; (t + e)-yt (t, + e)|](10)式中，I^ t, ti ^ N-2, t 幸 t1； t、ti 为整数，Yi (O和K (G)为公式(9)中八 O)在
采样点t、&的采样值；yi (t+e)和yi (t+e)为独立分量Ji (t)在采样点t+e和tje的采样值，e e {0，1，2}，求取差值的最大值，即为两个匕维矢量的距离；第八步在第四步中已计算得到阈值r，对于每一个t值，统计<& {t)Jh (^1)]小于r的数目，记为T1,1彡t, h ( N-2, t Φ t1； U1为整数，&可以取N-3个值，因此总的距离数为N-3，计算T1与总的距离数N-3的比值，记力W (/)即B； (r) = -^-(11)
Jy—j第九步计算皮1M的均值Bb^ {r) = ~^~MBHr) (12)
W-Z t=i第十步样本熵的计算公式为SiimpKii = - In[5*1 (r)/ Bh (/·)] (13)3)特征分类经过上面的计算，对于每个独立分量yi(t) {I ( i Sm}，可计算得三个特征量，将其依序构成如下特征向量，即
_2] f = [fu f2t /3；Γ(14)式中，分别表示yi(t)的峭度、renyi熵和样本熵；分别针对每一个特征进行一次k均值聚类分析，共进行三次分类，以f 为例，给出计算步骤3. I给定样本种类，选定初始聚类中心；显然，这是一个两分类问题，初始聚类中心选为 fn, f12 两个分量，即 Center1 (O) =fn, Center2 (O) =f12 ； center, (O)表示第 q 类初始聚类中心，q=l或2,第P次迭代计算的聚类中心用Centerq(P)表示；3. 2将样本按最小距离原则分为两个类，即在第P次迭代时，若I Ai-Center1 (P) I I < I I ^i-Center2 (P) I I，则 e C1 (p),否则 e c2 (p) ; | | · | | 表示
对其中的两个向量求取它们之间的欧式距离；cq(p)表示在第P次迭代计算中第q类的聚类域，fn e Cq(P)表示将样本fu划分到第q类，q=l或2 ;3. 3重新计算聚类中心；以本类所有样本的均值代替原有的聚类中心值；即centerq (夕+ O = 士 Σ ，Ii1为第q类样本总数；
nI3. 4若两次迭代计算的聚类中心相同，即Centerq(P-I)=Centerq(P),则计算完毕， 否则令P=P+I，返回3. 2 ；4信号的重构眼电伪迹的峭度值明显高于脑电信号，而renyi熵和样本熵的数值低于脑电信号，因此，综合三次分类的结果，将所有独立分量分为C1和C2两类，设C1类表示含有眼电伪迹的独立分量，C2类表示脑电独立分量；将C1类中所有成分置零，C2不变，得到去除眼电伪迹后的独立分量f，利用式(15)进行重构，即可得到纯净的脑电成分i ；X= H' 1J'(15)式中，w—1为第I部分求得的分离矩阵w的伪逆矩阵。本方法无需人工设定阈值就可以达到自动识别并去除眼电伪迹的目的，改善了以往方法中设定阈值时需要研究人员具备一定的先验知识、主观性强的不足
图I脑电信号中眼电伪迹自动去除的算法流程2 K均值聚类分析算法流程3国际10-20系统标准电极放置4 C3导联眼电伪迹去除效果5 Cz导联眼电伪迹去除效果6 C4导联眼电伪迹去除效果图
具体实施例方式本发明提出一种无需人工设定阈值的脑电信号中眼电伪迹的自动去除方法，完整的算法流程包括下面四(I. 2. 3. 4)个部分。其中，第I部分是现有方法，本发明申请的特征包括三(2. 3. 4)个部分完整的算法流程如图I所示。I独立分量分解本发明采用fastICA算法，该算法有如下优点①收敛速度比批处理和自适应处理都快；②负熵作为高斯性度量的效果优于累计量；③采用牛顿迭代法，收敛有保证；它的基本模型表达式为jc(t) = As{t)(2)
式中，χα)=[χια)，χ2α)，···，χηα)]τ e RnxN是η通道的观测向量，为测量得到，1≤t≤ N，N表示采样点总数，s (t) e RmXNS源信号，混合矩阵A e Rnxm，独立分量分析的
目标是求解分离矩阵w，通过w从观测信号X (t)中恢复出未知的源信号s (t)，源信号s (t)使用下式估计得到
权利要求
1.一种无需人工设定阈值的脑电信号中眼电伪迹自动去除的方法，其特征在于首先，对采集的脑电信号进行独立分量分解，得到m个独立分量；其次，对于每一个独立分量，求取它们的峭度、序列renyi熵和样本熵作为特征向量；然后，使用k均值聚类分析——一种无监督的分类算法自动识别出含有眼电伪迹的独立分量；最后，将含有眼电伪迹的独立分量置零，其余独立分量不变，对信号进行重构，得到纯净的脑电信号；包括以下步骤 1)独立分量分解 使用fastICA算法对信号进行分解，它的基本模型表达式为
全文摘要
本发明提出一种无需人工设定阈值的脑电信号中眼电伪迹自动去除的方法，属于生物信息技术领域，主要应用于脑电信号的预处理过程中。具体包括对采集的含有眼电伪迹的脑电信号进行独立分量分解；并且求取每个独立分量的峭度、序列renyi熵和样本熵作为特征向量，进而使用k均值聚类分析的方法自动识别出含有眼电伪迹的独立分量，并将其置零，其余分量不变，对信号进行重构，得到纯净的脑电信号。本发明解决了传统的眼电伪迹去除过程中需要人工对伪迹进行识别、费时费力、工作量大的问题，并且本方法无需人工设定阈值就可以实现自动识别并去除眼电伪迹的目的，弥补了以往方法中设定阈值时需要研究人员具备一定的先验知识、主观性强的不足。
文档编号A61B5/0476GK102835955SQ20121033119
公开日2012年12月26日 申请日期2012年9月8日 优先权日2012年9月8日
发明者李明爱, 李骧, 张译帆, 崔燕, 乔俊飞, 杨金福, 郝冬梅 申请人:北京工业大学