基于LMD和CSP的多域融合脑电信号分类方法

陈 舒 周 青

1(广州市公用事业技师学院信息传媒产业系 广东 广州 510100) 2(中山大学信息科技学院 广东 广州 510006)

0 引 言

脑-机接口(Brain Computer Interface, BCI)是一种将人脑与外部计算机等辅助设备直接联系起来的实时通信系统,BCI技术直接对脑电信号进行采集和分析,提取其中包含的运动意图信息并进一步将其解译为驱动外部设备的命令,能够帮助运动障碍患者恢复与外界的信息交互,提升其生活便利性[1-3]。脑电图(Electro-Encephalon-Graph, EEG)由于其低成本和非侵入等优点成为当前BCI系统中应用最为广泛的一种脑电信号。

EEG信号的产生机理复杂,是一种典型的非平稳、非线性和微弱的时变信号。如何从中提取反映患者不同思维活动的差异信息并将其与特定的控制命令联系起来是BCI技术的难点,也是当前研究的热点[4]。目前常用的方法可以分为时-频变换类方法和空间滤波类方法,其中典型的时-频变换类方法有短时傅里叶变换(Short Time Fourier Transformation, STFT)、小波变换ST和经验模态分解(Empirical Mode Decomposition, EMD)等。其中,FT将时域EEG信号转换到时-频域,能够较好地刻画信号的局部特征,但是受制于海森堡测不准原理,该方法无法同时获得较高的时间分辨率和频率分辨率[5];小波变换作为一种经典的非平稳、非线性信号分析方法,能够对信号进行不同“尺度”的细化分析,被广泛应用于EEG信号分类中,然而“小波基”函数的选择和分解层数的确定对小波方法的性能影响较大[6-8];EMD不需要预先设定“基函数”,直接从原始信号中提取本征模函数,能够同时获得高时间分辨率和频率分辨率,但是边缘效应的存在制约了该类方法的应用[9]。与时-频变换类方法不同的是,空间滤波类方法从提升EEG信号空间分辨率的角度出发,通过多通道的方式对不同脑域的特征信息进行提取和分析,典型的空间滤波类方法有独立主分量分析法,共空间模式法(Common Spatial Pattern, CSP)、空间稀疏共空间模式(Spatially Sparse Common Spatial Pattern, SSCSP)、滤波器组共空间模式(Filter Bank Common Spatial Pattern, FBCSP)等多种CSP变化方法[10-13]。该类方法相对于时-频变换类方法由于采用了更多观测通道,通常可以获得较高的分类结果,但是多通道的数据获取是以增加系统复杂度为代价的,而且该类方法丢失了频率信息。

上述分析表明,现有方法虽然能够提取出反映不同运动脑电信号中包含的差异信息,但是由于这些方法仅仅是从某个特定的变换域对EEG进行分析,提取的特征较为单一,没能包含不同维度的信息,具有一定的局限性,同时在低信噪比条件下会出现分类识别性能下降的问题[14-15]。

针对上述问题,本文提出一种基于局部均值分解(LMD)和CSP的多域融合脑电信号分类方法,首先利用LMD对EEG信号进行分析,得到各阶乘积分量(PF)的同时提取时-频域特征,然后对各阶PF分量进一步作CSP分解,获取EEG信号的空间分布特征,从而将时-频域特征扩展为时-频-空域特征。最后采用相关向量机(RVM)[16]对特征选择和分类函数进行联合优化,剔除冗余特征并得到最优的特征组合。基于BCI竞赛数据集的实验结果表明,所提方法可以实现对EEG信号的多角度描述,能够获得更优的分类性能。

1 实验数据

本文实验选用与文献[10]一致BCI2003 Competition Ⅱ中DatasetⅢ标准数据集。实验中对一名25岁健康女性左右手运动想象(MI)EEG信号进行采集和记录。EEG信号采用国际标准的10～20导联系统的C3和C4通道,采用差分电极。采样频率为128 Hz,并对频率为0.5～30 Hz范围内的信号进行滤波。实验中,受试者按要求坐在屏幕前,根据屏幕上出现的光标提示进行想象左右手运动,单次实验流程如图1所示。具体说明如下:

1) 0 s≤T<2 s,受试者前方屏幕显示空白,提示受试者放松。

2)T=2 s,受试者前方屏幕出现“+”光标,提示受试者做好准备。

3) 3 s≤T<9 s,受试者前方屏幕随机出现“←”或“→”,受试者根据提示想象左手或者右手运动。

4)T=9 s,受试者前方屏幕再次显示空白,提示受试者单次实验结束。

图1 单次实验数据采集时序

整个数据集中包含7组,每组40次实验,总共280组实验数据,包含想象左右手运动数据各140组。由于每次实验从第3 s正式开始,因此本文取每组实验数据中第3～9 s对应的768个采样点作为有效数据。

2 时-频-空域特征提取

2.1 基于LMD的时-频域特征提取

EEG信号的幅度和频率均随着时间的变化而变化,是一种典型的调频调幅信号,其瞬时频率为信号中调频成分的瞬时频率,幅度为信号中调幅成分的瞬时幅度,因此如果能将EEG信号中的调频成分和调幅成分进行分离,就能进一步对信号瞬时频率和瞬时幅度的细节信息进行分析,从而可以提取更多对EEG信号进行分类的有用信息。

LMD是在EMD基础上发展起来的一种新的非平稳复杂信号自适应分解算法,LMD可以将信号分解为若干个具有物理意义PF分量之和的形式,而获得PF分量的本质是从原始信号中分离出纯调频信号和包络信号,因此LMD非常适合于分析EEG信号。

对任意给定EEG信号,对其进行LMD分解的步骤如下:

(1) 计算均值函数m11(t)。

Step1对输入信号s(t)进行分析,计算其中包含的全部极值点ni,i=1,2,…,W。

Step2根据式(1)计算s(t)的均值序列mi:

Step3利用移动平均法对mi构成的曲线进行平滑,得到输入信号s(t)的均值函数m11(t)。

(2) 计算包络估计函数a11(t)。

Step1根据式(2)计算包络估计值ai:

Step2利用移动平均法对ai构成的曲线进行平滑,得到包络估计函数a11(t)。

(3) 计算差值信号h11(t)。

从原始信号中减去均值函数m11(t),得到差值信号h11(t)=s(t)-m11(t)。

(4) 计算纯调频信号s1j(t)。

Step1对差值信号h11(t)按式(3)进行幅度解调处理,得到s11(t):

Step2将s11(t)作为步骤(1) 的Step 1中的输入信号s(t),重复进行步骤(1)-步骤(3),直到第j次迭代后,步骤(2)得到的包络函数a1j(t)=1,此时得到的s1j(t)即为纯调频信号。

(5) 计算幅值函数a1(t)。

(6) 计算第一个PF分量f1(t)。

将a1(t)与s1j(t)相乘,得到PF分量f1(t)=a1(t)s1j(t)。

(7) 计算剩余信号u1(t)。

从原始信号s(t)中减去f1(t)得到残留信号u1(t)=s(t)-f1(t)。

(8) 计算所有L个PF分量。

令s(t)=u1(t),重复进行步骤(1)-步骤(7),直到获得L个PF分量,且剩余信号U(t)为单调函数时,迭代终止,此时原始信号可以表示为:

(4)

图2(a)和图2(b)分别给出了利用上述LMD方法对想象左手运动和想象右手运动脑电信号进行分析得到的分解结果。其中第1行子图为原始信号,第2-第5行子图为分解得到的PF1-PF4,最后一行子图为剩余信号。从图2可以看出,随着分解层数的增加,PF信号的频域特性逐渐减少,时域特性逐渐增加,因此PF分量中包含了信号的时域和频域特性。同时对比图2(a)和图2(b)可以看出,两种脑电信号分解结果呈现出了一定的差异性,因此我们可以提取以下时-频域特征对LMD分解结果的差异性进行描述。

(a) 想象左手运动

(b) 想象右手运动图2 LMD分解结果

特征1:波形宽度特征:

式中:find(·)为计算满足括号内要求的样本点个数,i=1,2,…,4。波形宽度特征能够体现信号频率分布差异。

特征2:波形熵特征:

(6)

特征3:波形四阶距特征:

(7)

对PF1-PF4分量分别提取上述三维特征,构成12维特征向量如表1所示。

表1 LMD分解得到的时-频域特征向量

2.2 基于CSP的空域特征提取

CSP是一种针对两类任务的经典空域滤波方法,能够从多通道EEG数据中提取每一类信号的空间分布模式,其基本原理是通过矩阵对角算法构建空间滤波器,将多通道EEG数据投影到低维子空间并使两类样本的协方差矩阵差异最大化。

对于实验中获得的N通道想象左、右手运动脑电信号数据矩阵分别为EL∈RN×T和ER∈RN×T,T为每个通道的采样点数。对其进行CSP分解的步骤如下:

(1) 计算归一化协方差矩阵。根据式(8)可以计算得到想象左、右手运动脑电信号数据矩阵对应的归一化协方差矩阵RL和RR。

式中:trace(·)表示矩阵求迹运算。

(2) 构建复合协方差矩阵并进行特征值分解。

Step1对RL和RR求和得到复合协方差矩阵R=RL+RR。

Step2对R进行特征值分解得到特征值矩阵Σ。

R=UΣU

(9)

式中:U为特征向量矩阵;Σ为特征值构成的对角矩阵。

(3) 对RL和RR进行白化处理。

Step1根据式(10)计算白化矩阵P。

(10)

Step2根据式(11)利用白化矩阵P对RL和RR进行白化处理。

(4) 对SL和SR进行特征分解。

根据式(12)对SL和SR进行特征分解:

(12)

根据式(13)利用空间滤波器W对多通道信号矩阵E进行空域滤波得到特征矩阵:

Z=WE

(13)

计算特征矩阵前m行和后m行对应的方差构成空域特征向量:

特征4:空域特征。

(14)

式中:var(·)为求方差运算。本文取m=9,即总共有18维特征构成空域特征向量。

3 基于RVM的特征选择和分类

得到特征向量后,需要设计合适的分类算法才能获得期望的分类性能,目前常用的分类算法有线性分类器、支撑向量机分类器和贝叶斯分类器等。本文第2节提取的30维时-频-空域特征虽然实现了特征域的扩展,但是高维特征向量中难免会存在一些冗余特征,这些特征的存在会增加分类算法的复杂度,降低实时性,因此本文采用RVM进行特征选择和分类识别的联合优化,降低算法复杂度的同时,提升分类性能。

RVM理论是在SVM基础上发展起来的一种基于贝叶斯框架的监督类学习方法,其基本原理也是利用核函数将低维空间中的线性不可分问题映射到高维空间中的线性可分问题,但不同于SVM的是,RVM具有更高的稀疏性,核函数的选择不受摩西准则[16]限制,并且能够提供概率式的预测。

对于想象左手运动和想象右手运动脑电信号分类问题,完整的RVM模型可以表示为:

式中:xn(n=1,2,…,N,N为样本总数)为第n个训练样本的特征向量;K(x,xn)为核函数;w=[w1,w2,…,wn]T为服从0均值、协方差矩阵为α-1I高斯分布的权重向量;ε为服从0均值、协方差矩阵为γ-1I高斯分布的噪声分量,与信号分量相互独立。为了构建完整的贝叶斯模型,RVM进一步假设α和γ服从超参数为a0、b0、c0和d0的伽马分布。

图3利用RVM算法对上述时-频-空域33维特征向量进行特征选择迭代终止时的归一化权值,可以看出经过特征选择后,权值较大的有6维特征,对应4维时-频域特征Feature11、Feature12、Feature24和Feature33,以及2维空域特征Feature47和Feature416,体现了本文方法信息提取维度的多元化。

图3 RVM特征选择结果

4 实验及结果分析

图4给出了本文所提基于LMD和CSP的多域融合脑电信号分类方法流程,可以看出算法包含训练和测试两个阶段。

图4 本文分类算法流程

在训练阶段,首先对训练样本进行LMD分解,得到PF分量并提取12维时-频域特征,然后将PF分量作为CSP的多通道数据进行空域滤波并提取18维空域特征,然后利用RVM分类器对30维时-频-空域特征进行特征选择和分类。

在测试阶段,对测试样本进行LMD分解,然后提取训练阶段RVM分类器选出的时-频域特征,同样将PF分量作为CSP的多通道数据进行空域滤波,然后提取训练阶段RVM分类器选出的空域特征,并将其与时-频域特征一起构成特征向量送入RVM分类器进行分类。

4.1 单一维度特征和时-频-空三维特征对比

为了验证时-频-空域融合特征的分类性能,将数据集随机分成5个子集,采用5组交叉验证法进行实验,即依次选取1组子集作为训练样本,剩余4组子集作为测试样本,记录每次实验的分类结果如图5所示,表2给出了5次实验的统计结果。为了进行对比,图5和表2给出了在相同条件下分别单独采用时-频域特征、空域特征进行分类的结果。从图5和表2的结果可以看出,所提时-频-空域特征相对于单一的时-频域特征和空域特征分类性能提升了5%,并且对不同数据集的分类稳定性更好。同时只采用RVM分类器选择的6维时-频-空域特征与直接采用30维特征的分类性能接近、方差更小,表明RVM选择的6维特征对分类性能的贡献率要远超其他24维特征。

图5 不同特征单次分类结果

表2 5次实验平均分类结果

4.2 不同算法分类性能对比

表3给出了本文方法与近年来其他文献中介绍方法在相同条件下的分类结果和F1得分指标,其中文献[17]采用的是小波包能量特征以及线性分类器,文献[18]采用的是小波熵特征以及SVM分类器;文献[19]采用的是5阶AR模型系数以及线性分类器。从表3可以看出,本文方法可以获得最优的分类性能。

表3 不同方法的分类性能

除了分类性能外,算法执行时间也是评价算法优劣的一项重要指标。表3给出了不同方法在MATLABr2015平台的算法执行时间,可以看出本文算法的实时性明显优于文献[18]和文献[19]所述方法,仅比文献[17]采样线性分类器慢0.4 s,表明本文方法兼具高分类性能和高算法实时性的优点。

4.3 不同算法低信噪比条件下鲁棒性对比

微弱性是脑电信号的一个显著特点,其电压通常为μV量级,容易受到噪声污染。因此噪声稳健性是评估分类算法优劣的另一项关键指标。为了验证本文方法在不同信噪比条件下的分类性能,采取向实验数据中加入高斯白噪声的方式进行验证。

图6给出了上述方法在不同信噪比条件下的分类性能。可以看出,本文方法在低信噪比条件下的表现要明显优于其他对比方法,当信噪比高于15 dB时,本文方法可以获得优于90%的分类正确率。这是由于本文方法从时-频-空多个维度进行特征提取,降低了噪声分量的影响,同时RVM分类器具有的噪声鲁棒性也提升了低信噪比条件下的鲁棒性。

图6 分类性能随信噪比的变化曲线

5 结 语

针对运动脑电信号非平稳、时变和微弱性的特点,本文提出一种基于LMD和CSP的多域特征融合分类方法,综合考虑了EEG信号的时-频-空域信息。算法首先利用LMD对EEG信号进行时-频域特征提取,并将LMD分解得到的PF分量作为CSP的多通道数据进行空域特征提取,从而得到30维时-频-空域特征向量。然后利用RVM对高维特征向量进行自适应特征选择与分类,从30维特征中自动选取对分类性能贡献最大的6维特征,在降低算法复杂度的同时保证了分类性能。基于实际数据的实验结果表明,相对于单一维度时-频特征和空域特征,本文方法识别率更高、标准差更小、鲁棒性更好,与现有方法的对比结果表明,本文方法可以获得更高的分类性能和噪声稳健性,具有较好的应用前景。