基于样本熵的运动想象分类方法

基于样本熵的运动想象分类方法

1基于前脑电特征的运动想象分类方法脑机接口（bci）是近10年来发展起来的一种人机接口方式。它不是依赖于大脑的正常输出路径（即外部神经和肌肉），而是基于脑的计算和控制技术。作为一种直接投入大脑和外部交流的手段，脑机接口在脑认知科学、军事、医疗恢复等领域具有广阔的应用前景。已成为神经工程研究的热点课题。然而,目前的脑机接口技术远不能达到心想事成的地步,即不是人的所有的思维活动都能被正确地识别,我们只能寻找那些具有明显脑电特征的思维活动进行分类,用来与外界交流和进行控制.基于左右手运动想象的脑机接口,因其属于内源性诱发脑电,有明显的统计特征,无须大量训练,而且不依赖人的其他感官系统,无须额外的刺激装置,因而受到研究者们更多的关注.这一领域的先驱者——奥地利Graz大学的Pfurtscheller等人——已经将运动想象脑机接口和功能电刺激相结合,重建了一名脊髓损伤患者的手的运动功能,使其能够完成日常基本活动.对于基于左右手想象意识任务的脑机接口,其实现的关键技术在于如何进行有效特征提取和正确分类.目前特征提取方法主要是依据左右手想象时EEG不同频段的能量变化特性,即左右手运动想象时,会在大脑对侧感觉运动皮层的EEG信号α、β波段出现事件相关去同步(ERD)现象,在大脑同侧出现事件相关同步(ERS)现象.基于此原理,诸如频段能量特征、AR模型系数、小波系数都被用于特征向量加以分类.但实践表明,这种基于频带能量的特征提取方法对分类正确率的进一步提高很难再有大的突破,因此,寻找左右手运动想象的其它有效特征,也许会对此类脑机接口的性能有大的改善.基于此,本研究首次提出了采用非线性动力学参数——样本熵——作为脑电特征进行运动想象分类的思想.分析了左右手运动想象时感觉运动皮层的脑电信号样本熵动态变化规律,证明了样本熵作为特征进行分类的可行性,并阐明了其变化规律的神经电生理意义.在此基础上,利用Fisher线性分类器进行了左右手运动想象动态分类,得到了满意的结果,平均最大分类正确率达到87.8%.本文还提出了一种样本熵的快速算法,其计算量小、速度快,满足BCI实时系统要求,具有较高的实用价值.2样本熵原理principlesyner2.1样本熵与1.2模型拟合样本熵是在Pincus近似熵概念的基础上,于2000年由Richman提出的一种改进的复杂度测量方法.它不但具备近似熵所有的优点,而且避免了近似熵中统计量的不一致性.本研究之所以选用样本熵作为特征对左右手运动想象进行分类,是考虑样本熵比其他非线性参数具备以下明显优点:①样本熵的物理意义明确,能够衡量时间序列中产生新模式概率的大小.样本熵值越大,产生新模式的概率越大,序列越复杂.这在脑电分析中,正好可以用来衡量大脑皮层被激活时脑电活动的复杂性变化.②计算所需数据短,算法运算量小.只需要很短的数据(100～5000点)便可得到稳健的样本熵估计值,这对于脑机接口的实时应用是非常重要的.其他一些非线性动力学参数(如相关维数、李亚普诺夫指数等)经研究虽然也能提取脑电信号有效的非线性特征,但由于对其计算往往需要大量数据才能得到稳定估计值(一般最少在几万点以上),而且算法复杂、速度慢,从而失去了在脑机接口中的应用价值.③具有较好的抗噪和抗干扰能力.脑机接口中采集EEG时,不可避免地受到各种伪迹和工频干扰,样本熵对这些干扰都有较好的抗噪能力,从而增强了脑机接口的实用性.④对确定信号和随机信号都适用.这一点对于分析本身就是确定性成分和随机成分组合的脑电信号非常有利.⑤样本熵具备相对一致性,比近似熵更符合理论性.在近似熵算法中,数据和其自身之间的比较会产生误差,而样本熵不对自身匹配值计数,因为熵是新信息产生率的测度,比较数据和其自身的匹配毫无意义,所以样本熵统计量可以维持相对一致性而近似熵做不到.样本熵的定义算法描述如下:对于给定的N点时间序列{u(i)},其样本熵通过如下步骤计算得到(其中m是预先选定的模式维数,r是预先选定的相似容限).①将序列{u(i)}按顺序组成m维矢量,即Xm(i)=[u(i),u(i+1),⋯,u(i+m-1)](i=1~Ν-m)(1)Xm(i)=[u(i),u(i+1),⋯,u(i+m−1)](i=1~N−m)(1)②对每一个i值计算矢量Xm(i)与其余矢量Xm(j)之间的距离:d[Xm(i),Xm(j)]=maxk=0~m-1︳u(i+k)-u(j+k)︳(i,j=1~Ν-m;i≠j)(2)③给定一个相似容限r值(r>0),对每一个i值统计d[Xm(i),Xm(j)]<r,(i,j=1～N-m;i≠j)的数目,以及此数目与总的距离数目N-m-1的比值,记作Bmi(r),即:Bmi(r)={d[Xm(i),Xm(j)]＜r的数目(i≠j)}/(Ν-m-1)(3)④求Bmi(r)对所有i值的平均值,记做Bm(r),即Bm(r)=1Ν-mΝ-m∑i=1Bmi(r)(4)⑤增加维数为m+1,构造一个m+1维矢量,即Xm+1(i)=[u(i),u(i+1),⋯,u(i+m)](i=1~Ν-m)(5)⑥对每一个i值计算矢量Xm+1(i)与其余矢量Xm+1(j)之间的距离:d[Xm+1(i),Xm+1(j)]=maxk=0~m︳u(i+k)-u(j+k)︳(i,j=1~Ν-m;i≠j)(6)⑦对每一个i值统计d[Xm+1(i),Xm+1(j)]<r,(i,j=1～N-m;i≠j)的数目,以及此数目与总的距离数目N-m-1的比值,记作Bm+1i(r),即:Bm+1i(r)={d[Xm+1(i),Xm+1(j)]＜r的数目(i≠j)}/(Ν-m-1)(7)⑧求Bm+1i(r)对所有i值的平均值,记做Bm+1(r),即Bm+1(r)=1Ν-mΝ-m∑i=1Bm+1i(r)(8)理论上次序列的样本熵为:SampEn(m,r)=limΝ→∞{-ln[Bm+1(r)/Bm(r)]}(9)当N为有限值时,得到序列长度为N时的样本熵的估计值,记作:样本熵的值显然与m、r的取值有关,根据经验,通常取m=2,r=0.1～0.25SD(u),其中SD(u)表示原始时间序列的标准差.2.2dij,nn的二值距离矩阵为了将样本熵用于运动想象的实时分类,本文在分析样本熵定义的基础上,参考洪波等人近似熵的快速算法,引入了二值距离阵的概念,提出了一种计算时间序列样本熵的实用快速算法,实验证明,此快速算法减少了定义中的冗余计算,大大提高了计算速度,在脑机接口实时系统中有很大应用前景.此快速算法的计算步骤如下:第一步:对N点序列,定义其N×N的二值距离矩阵为D,D的第i行第j列的元素记作dij:dij={1,︳u(i)-u(j)︳＜r0,︳u(i)-u(j)︳≥r,(i,j=1~Ν)(11)第二步:利用矩阵D中的元素,可以方便地计算B2i(r)和B3i(r){B2i(r)=1Ν-3Ν-1∑j=1dij∩d(i+1)(j+1)(i,j=1~Ν-2;i≠j)B3i(r)=1Ν-3Ν-2∑j=1dij∩d(i+1)(j+1)∩d(i+2)(j+2)(i,j=1~Ν-2;i≠j)(12)实际计算时,取m=2,按照行递增的顺序,每两行或每三行的矩阵元素斜对角的方向组合进行“与”运算,把每一行斜线的求“与”的结果累加后除以N-m-1,即可得到B2i(r)和B3i(r).其中,两行和三行的求“与”过程可以放在同一个循环中进行,进一步简化运算量.第三步:计算B2(r)和B3(r){B2(r)=1Ν-2Ν-2∑i=1B2iB3(r)=1Ν-2Ν-2∑i=1B3i(13)第四步:计算样本熵SampEn(2,r,Ν)=-ln[B3(r)/B2(r)](14)3单次实验流程受试对象为3名在校大学生,均为男性,右利手,无神经系统病史,年龄在23～27岁之间.受试者坐在一个舒适的高靠背椅子上,全身肌肉处于放松状态,不产生任何的肌肉紧张与实际运动.要求受试者根据屏幕提示进行相应的左右手抓握动作想象.单次实验流程如图1所示.前2秒,屏幕显示固定的十字叉,受试者此阶段处于思维放松状态;在第2秒开始时,产生一个提示音,提示受试者集中注意力,等待想象指令的出现;第3～4.25秒时,屏幕出现左右箭头提示,受试者根据提示内容,分别想象左右手抓握动作.受试者的想象思维活动一直持续到第8秒,直到下一个十字叉出现时结束.为了增加脑电信号的空间分辨率和定位性能,本研究采用双极同心圆环拉普拉斯电极采集大脑C3、Cz、C4所在位置头皮的拉普拉斯电位(LaplacianElectroencephalography,LEEG),采样率为128Hz.每个受试者按上述试验流程执行左右手想象各300次,去除大的眼电肌电干扰数据后,随机抽取50%的数据作为测试集,其余作为训练集.4左、右运动想象下的eeg样本熵erd时程为了分析左右手运动想象时,大脑感觉运动皮层EEG信号的复杂度变化,本研究参考了Pfurtscheller的思想,即采用α、β频段能量量化ERD/ERS随时间变化,通过计算左右手运动想象时的EEG样本熵相对于参考时间段的百分比来量化复杂度ERD时程.主要步骤为:首先从单次实验流程的第一秒开始,以1s长度为滑动时间窗(即128点),计算该采样点处前1秒时间C3、C4脑电信号的样本熵.每次移动窗口一个采样点,直到计算出最后一秒样本熵为止,从而得到在单次试验流程中C3、C4脑电信号样本熵的时间序列.然后将同一运动想象任务的所有试验的样本熵序列叠加平均,得到单个受试者的平均样本熵序列.再计算各时间点处平均样本熵相对于参考时段的百分比(参考时段为0～2s),从而得到单个受试者的样本熵ERD时程.最后将不同受试者的样本熵ERD曲线叠加平均,得到最终的统计样本熵ERD序列.因为样本熵的物理意义为序列的复杂度,所以EEG样本熵ERD时程的本质表征了EEG复杂度随时间的变化规律.图2为经过上述步骤计算得到的左右手运动想象样本熵ERD时程曲线.从图中可以看出,左右手运动想象时,C3、C4导联的EEG样本熵有明显的区别.从第3秒以后,当左手运动想象时,大脑同侧C3位置的EEG信号的样本熵明显降低,即复杂度明显降低,同时大脑对侧C4位置的EEG信号的样本熵略有上升,即有复杂度升高的趋势;与之相反,当进行右手运动想象时,同侧C4位置的EEG复杂度明显降低,而对侧C3位置的EEG复杂度明显升高.试验结果与α、β频段能量法所表征的ERD曲线规律正好相反.结合神经解剖学与神经生理学知识可知,当想象左右手运动时,大脑对侧的感觉运动皮层的手功能区被激活,该区域的代谢和血流量增加,大量神经元的理化活动增强,其同时进行信息加工导致了脑电信号独立性增强,同步化程度降低,进而导致了α、β频段能量降低(ERD),脑电复杂度反而升高的现象.与此同时,大脑同侧感觉运动皮层的对应区域相对处于静息或惰性状态,神经元活动被抑制而表现出宏观脑电信号的一致性,因而导致α、β频段能量升高(ERS)和脑电复杂度的降低.所以,从以上分析可知,用样本熵所表征的EEG复杂度特征和经典的能量特征本质上具有相同的生理基础,从而也证明了利用样本熵作为EEG复杂度特征进行左右手运动想象任务的分类是完全可行的.5资源提取和分类fortandclassgo5.1维特征分区本研究采用上文提出的快速算法直接计算出单次实验流程C3、C4导联的EEG信号样本熵的时间序列SampEn_C3(t)和SampEn_C4(t),并将其组合成二维时变特征向量Ft=(SampEn_C3(t),SampEn_C4(t)),用于后续运动想象任务的实时连续分类.每一时刻样本熵的值均以该时刻前1秒脑电信号进行计算,相似容限r取0.25SD(u).5.2fps的动态连续分类目前,脑机接口的分类方法很多,诸如Mahalanobis距离、神经网络、支持向量机等,本研究采用了Fisher线性判别式实现对左右手运动想象意识任务动态连续分类,取得了较好的效果.首先分别计算训练样本和测试样本的时变特征向量,由训练样本特征向量得到每一个采样点处Fisher判别式的权值系数W(t)和阈值ω0,然后对测试样本的每个采样点特征向量计算Fisher判别式距离,判断其类别属性,从而实现动态连续分类.线性判别函数为:g(x)=WΤx+ω0(15)其中x为被测试特征向量,当g(x)>0时,决策x属于第一类(左手想象);当g(x)<0时,决策x属于第二类(右手想象);当g(x)=0时,决策x不属于任何一类.利用Fisher准则,计算W的公式为:W=(SL+SR)-1(uL-uR)(16)其中SL、SR分别为左右手想象训练集特征向量的类内聚散度矩阵,uL、uR分别为左右手运动想象训练集特征向量的均值.ω0的取值根据经验用下式进行估计:ω0=-WΤ(uL+uR)2(17)6模型分类结果使用二维时变EEG样本熵作为特征向量,并用Fisher线性判别式加以分类,得到每个受试者测试集随时间变化的连续分类结果,最后对不同受试者分类结果进行叠加平均得到平均分类正确率曲线,如图3所示.在t=4.8s附近,平均分类正确率达到最大87.8%.实验结果表明,样本熵能很好地表征左右手运动想象时的EEG复杂度特征,对其分类取得了较好的效果.7实验结论与讨论本文以样本熵作为特征向量对左右手运动想象任务进行了动态分类研究,取得了较好的效果,为脑机接口的特征提取提供了一种新的思路.实验表明,样本熵可以表征EEG的复杂度变化,有明确的生理意义,而且在左右手运动想象时有明显的区别,有利于BCI分类正确率的提高,是一个很有实用价值的特征参数.另外,本文还对EEG信号的互样本熵做了分析,将Cz导联