基于噪声辅助快速多维经验模式分解的运动想象脑电信号分类方法
郑潜 1 , 乔丹 1 , 郎恂 2 , 谢磊 1 , 李东流 3 , 王琪冰 3 , 苏宏业 1
1 浙江大学智能系统与控制研究所,浙江 杭州 310027
2 云南大学信息学院电子工程系,云南 昆明 650091
3 森赫电梯股份有限公司,浙江 湖州 313000
【摘 要】 脑机接口是一项新兴的技术,它可以处理分析采集到的运动想象脑电信号,从而实现对外部辅助设备的控制。针对目前运动想象脑电信号处理方法计算效率低、分类准确率不高等问题,提出了一种新的基于噪声辅助快速多维经验模式分解(NA-FMEMD)的运动想象脑电信号分类方法。该方法首先利用 NA-FMEMD 得到全部的多维本征模式函数和趋势项;接着,根据平均频率选取特定的信号层,构建出新的多维信号;然后,通过共空间模式提取出脑电信号的特征向量;最后,将特征向量输入支持向量机分类器中进行分类。分别采用仿真数据和BCI Competition IV数据进行测试,并与基于噪声辅助多维经验模式分解(NA-MEMD)的方法进行比较,验证了所提方法的有效性和优势。 【关键词】 脑电信号 ; 运动想象 ; 噪声辅助快速多维经验模式分解 ; 共空间模式
郑潜,乔丹,郎恂, 等. 基于噪声辅助快速多维经验模式分解的运动想象脑电信号分类方法[J]. 智能科学与技术学报, 2020, 2(3): 240-250.
ZHENG Q, QIAO D, LANG X, et al. Classification of motor imagery signals using noise-assisted fast multivariate empirical mode decomposition[J]. CHINESE JOURNAL OF INTELLIGENT SCIENCE AND TECHNOLOGIE, 2020, 2(3): 240-250.
脑机接口(brain-computer interface,BCI)是一项新兴的技术,可以对大脑皮层神经元活动产生的电信号进行采集,分析信号特性并产生控制信号,进而实现对外部设备的直接控制。其中,脑电图(electroencephalogram,EEG)由于具有非侵入性、经济实惠、便于实时操作等特征,应用最为广泛。利用EEG信号控制外部辅助设备,不仅可以为严重残疾的患者提供与外界交流的渠道,改善他们的生活环境,还可以促进生物工程领域、神经修复学等多学科交叉技术的发展。 Silva G P H L D等人的研究表明,当人对肢体动作进行想象时,大脑的主感觉运动皮层处于4~8 Hz频段的θ节律、8~13 Hz频段的μ节律以及13~30 Hz频段的β节律会发生与该动作真实实施时相似的电生理响应,诱发运动想象(motor imagery,MI)电位。人类对肢体动作的准备或者想象将使皮质运动中枢中大量神经细胞的活动状态发生变化,导致脑电信号中某些频率增强或减弱。Nikouline V V等人的研究证明,当人做单侧肢体实际运动或者单侧肢体想象运动时(如左手实际运动或运动想象),大脑对侧(即右侧)主感觉运动区的 μ节律幅值会明显降低,即事件相关去同步(event-related desynchronization,ERD),大脑同侧(即左侧)主感觉运动区的 μ节律幅值会明显增高,即事件相关同步(event-related synchronization,ERS)。 由于脑电信号中噪声较多,原始脑电信号的信噪比低,可靠性差,如何提取脑电信号的特征成为脑电信号分析、分类非常核心的难点。目前常用的方法有傅里叶变换(fourier transform,FT)、功率谱密度(power spectral density,PSD)、自回归模型(autoregressive model,AR)、小波变换(wavelet transform,WT)、相关分析、共空间模式(common spatial pattern,CSP)等。CSP方法是一种利用主成分分析和时空域变换的特征提取方法。该方法没有考虑 EEG 信号中的频率特性,而运动想象脑电信号的ERD和ERS发生在特定区域的特定频段。因此,结合共空间模式和信号时频分解的方法,成为运动想象脑电信号处理研究的趋势。 本文提出了一种新的基于噪声辅助快速多维经验模式分解(noise-assisted fast multivariate empirical mode decomposition,NA-FMEMD)的运动想象脑电信号分类方法。该方法的主要贡献在于解决了噪声辅助的多维经验模式分解(noise-assisted empirical mode decomposition,NA-MEMD)的模式混叠、计算效率低等问题,分解结果的可解释性更强,能够更加精准地定位相应频段,从而提高脑电信号的分类准确率。本文分别采用仿真数据和 BCI Competition IV 数据进行测试,验证了该方法的有效性和优势。 2.1 噪声辅助快速多维经验模式分解
1998年,Huang N E等人提出了经验模式分解(empirical mode decomposition,EMD)算法,该方法完全由数据驱动,适用于非线性、非平稳数据。2009年,Wu Z H等人提出了一种噪声辅助的改进算法——集合经验模式分解(ensemble empirical mode decomposition,EEMD)算法,该算法利用白噪声具有频率均匀分布的统计特性,避免了EMD算法的模式混叠。2009 年,Rehman N等人提出了多维经验模式分解(multivariate empirical mode decomposition,MEMD)算法,可同时对多通道数据进行分解,能够避免本征模式函数(intrinsic mode function,IMF)在各个数据通道之间频率不一致的问题。2011 年,Ur Rehman N 等人提出了NA-MEMD 算法。但随着信号维度的增加, NA-MEMD 算法的计算效率明显下降,高计算负荷和过度分解在一定程度上限制了该算法的应用。为了解决这个问题,Lang X等人首先提出了直接和间接的多维本质时间尺度分解(multivariate intrinsic time-scale decomposition,MITD)算法。2018年,Lang X等人又提出了快速多维经验模式分解(fast multivariate empirical mode decomposition, FMEMD)算法以及其噪声辅助的版本NA-FMEMD,计算效率得到明显提升。 NA-FMEMD 算法通过优化多维运算空间至一维运算空间的映射关系,很大程度上缩短了计算时间,尤其在信号维度增加时效果愈发显著。对于由p 维输入与 m 维白噪声组成的多变量信号 ,NA-FMEMD 算法的具体过程如下。 步骤一 设定 。同时设置分解的层数为i=1。利用低差异序列生成K维空间内均匀分布的方向向量集合 。 步骤二 依次计算信号 沿着第 k 个方向矢量 的投影函数 。 步骤三 使用标准EMD算法提取K个投影函数 的第一层本征模函数 。 步骤四 结合全部的第一层本征模函数 及其对应的方向向量 ,通过求解式(1)的超定方程得到多维本征模式函数(multivariate intrinsic mode function,MIMF) :
步骤五 将该多维本征模态函数 从原信号中分离,得到序列 。
步骤六 如果 单调,或者极值少于3个,则停止迭代过程,趋势项 ,IMF层数i =M;否则,更新 ,返回步骤二。
步骤七 得到全部的多维本征模态函数信号 和趋势项 。
2.2 平均频率计算 如前所述,脑电运动想象信号的变化主要体现在3个不同的频段上,分别被称为θ波(4~8 Hz)、μ波(8~13 Hz)、β波(13~30 Hz)。本文拟利用NA-FMEMD 算法分解脑电信号,根据一定指标自动提取对应于这些频段的MIMF。由NA-FMEMD/FMEMD的滤波器特性可以看出,每个分量的频率响应呈现窄带带通滤波器特性,围绕一个中心频率波动,且阻带衰减速度快。因此,可以根据每个MIMF的中心频率定位其所处的频带。根据统计学知识,MIMF的瞬时频率的平均值可以被作为其中心频率,表征其所处频带的数值信息。 因此,通过求取MIMF的平均频率,在时频分解的全部输出结果中筛选有效成分,以进行后续的特征提取等操作。本文利用希尔伯特变换构造的复解析信号的相位的导数来定义信号的瞬时频率。通过各层本征模态函数信号的瞬时频率可计算出平均频率。 计算某一个本征模态函数c(t)的平均瞬时频率时,首先对c(t)做希尔伯特变换:
其中,*表示卷积。进而得到c(t)的解析信号z(t)=c(t)+ jH[c(t)]。则c(t)可被表示为:c(t)=a(t)cosφ(t)。其中, 为幅值函数, 为相位函数。可计算出平均瞬时频率:
2.3 特征提取
共空间模式是运动想象脑机接口系统中提取运动想象信号特征非常有效的方法之一。其原理为寻找一组空间滤波器,使得在这组滤波器的作用下,一类信号的方差达到极大值,另一类信号的方差达到极小值,从而实现特征提取。该方法主要包含以下步骤。
步骤一 对于某类运动想象任务,计算其规范的平均空间协方差: 。其中, ϕ c 对应该类运动想象任务的集合,包含属于该类别的所有试验; 表示该类运动想象任务的第i次实验采集到的数据。 步骤二 计算混合空间协方差矩阵 ,对其进行特征值分解: 。 表示矩阵 R分解后的特征向量矩阵, 表示矩阵R 分解后与 对应的特征值的对角矩阵。 步骤三 求出白化特征值矩阵 ,使用白化特征值矩阵P对两类运动想象任务的平均规范化空间方差 分别进行白化处理: 。 步骤四 分解白化后的矩阵,求出共同特征向量 及对角矩阵 和 ,, ,其中I表示单位矩阵。 步骤五 得到空间投影矩阵W,即最优空间滤波器为: 。
步骤六 计算两类信号投影: 。
步骤七 选取投影后信号的前m行和后m行,得到两类任务的特征向量,分别为:
其中, 。
值得注意的是,CSP算法中m值对应特征向量的维度。由于信号投影矩阵(特征矩阵) 或者 的信息是不等效的,特征信号主要集中于特征矩阵的头部和尾部 ,因此需要引入变量m来确保特征向量中信息的完整性。 2.4 支持向量机分类器
本文采用径向基函数(radial basis function, RBF)的非线性支持向量机(nonlinear support vector machine)算法对特征进行分类。以固定比例将数据集随机分成训练集和测试集,利用 SVM 对训练集数据进行训练以建立分类模型,最后利用测试集数据计算预测准确率。 本文提出基于NA-FMEMD的运动想象脑电信号分类方法,其流程如图1所示。首先,对原始的EEG信号进行8~30 Hz的带通滤波,消除来自非接触式采集系统的高频噪声。其次利用NA-FMEMD算法对经过预处理的信号进行分解,得到全部的本征模态函数信号和趋势项。然后,根据本征模态函数信号的各层平均频率选取信号层,构建新多维信号,通过CSP算法进行特征提取。最后,将特征向量输入SVM分类器中进行分类,根据分类准确率选取CSP中的最优模型参数m,利用此最优参数下的脑电特征对脑电信号进行分类。
图1 基于NA-FMEMD的运动想象脑电信号分类方法流程
由图1可知,本文没有通过运动想象脑电信号本身的特性确定模型参数(即 CPS 算法中的 m值,对应于特征向量的维度),而是通过遍历m的值(m=1,2,3,...) ,根据最终的分类准确率寻求最优m值。
3.1 仿真数据实验
使用一组仿真数据来验证本文提出的方法。仿真信号由正弦信号和随机高斯白噪声合成而来,4个通道合成信号的构造方法如下。
每个通道的 3 个正弦波信号的频率f i,j 是随机生成的,用来模拟不同频段的脑电信号,分别对应 θ波(4~8 Hz)、μ波(8~13 Hz)、β波(13~30 Hz);初始相位角ϕi,j 在区间 [0,2π)随机生成;幅值Ai,j 在区间 [0,1]随机生成。4个通道的频率、初始相位角、幅值的生成相互独立。ωx i 为信噪比为 20 dB 的随机高 斯白噪声。 每次生成长度为2 s的4个通道信号,将重复4次的结果拼接为长度为8 s的脑电信号仿真数据,如 图2 所示。
利用NA-MEMD算法、NA-FMEMD算法将此4 个通道的仿真数据与两个通道的噪声信号进行时频分解,分解结果分别如图3、图4所示。其中第一行为4个通道的原始信号,d 1 ~ d8 分别为各通道数据经过时频分解得到的IMF分量,r为分解的残留趋势项。 两种分解的结果图中,不同层的IMF分量包含不同频段的时间序列信号,随着层数的提高,信号频率降低。
图3 利用NA-MEMD算法对仿真信号进行分解的结果
图4 利用NA-FMEMD算法对仿真信号进行分解的结果 计算两种时频分解方法得到的 IMF 的平均频率,并与仿真信号的真实正弦波信号的频率进行比较,结果见表1、表2、表3,“/”分隔的4个数值分别为0~2 s、2~4 s、4~6 s、6~8 s的频率。
由分析结果可知,NA-MEMD算法成功将β波、α波、θ波分别分离至第2层、第3层、第5层;NA-FMEMD算法将它们分离至第3层。而NA-MEMD算法存在过度分解的情况,第4层信号包含部分μ波频段、部分θ波频段,该算法对各频段信号的提取效果不如NA-FMEMD算法。同时,NA-FMEMD算法的分解运行时间为1.05 s,与NA-MEMD算法的6.89 s相比,计算效率提升了5.56倍。 3.2 脑电信号数据实验
实验数据来自德国柏林研究中心 2008 年的BCI Competition IV运动想象脑电数据集,共包含4个受试者a、b、f、g的实验数据。每个受试者选择两个运动想象任务(a 选择了左手和脚,b 选择了左手和右手,f选择了左手和脚,g选择了左手和右手),进行了200次实验,每个任务各100次。本文使用受试者g的实验数据,共200次实验(左手/右手运动想象任务各100次)。每组数据含59个通道的EEG信号,采样频率为100 Hz。关于这59个通道的具体描述请见参考文献[24]。由于人脑运动想象的反应与大脑中心区域有关,而 59 个通道中有11个通道位于大脑的中心区域,因此本文选择这11个通道进行分析,分别为FC3、FC4、Cz、C3、C4、C5、C6、T7、T8、CCP3、CCP4。 为了增强特定频段的脑电特征并且降低外界因素的干扰,本文首先采用8~30 Hz巴特沃斯带通滤波器对所有受试者(a、b、f、g)的每组运动想象脑电数据进行滤波预处理,提高信噪比。然后,对每组 11 个通道的数据进行时频分解,由于篇幅限制,本文仅选取4个通道进行展示。 将每次实验得到的多维信号和两个通道的噪声信号利用NA-FMEMD算法进行时频分解,其中某次右手运动想象的 4 个通道(C3、C4、CCP3、CCP4)实验结果的分解结果如图5所示。其中第一行为4个通道的原始信号,d 1 ~ d7 分别为各通道数据经过时频分解得到的频率从高至低排序的 IMFs分量,r 为分解的残留趋势项。另外,由于NA-FMEMD算法拥有模式对齐属性,某一层MIMF中各个通道之间的频带基本一致,因此笔者同样选择4个通道的平均频率进行展示,结果见表4。可见,利用NA-FMEMD算法对运动想象任务采集的脑电信号进行分解,结果中第1层、第2层信号分别对应β波、μ波,选择这两层的信号进行特征提取,可以更准确地提取出特征信息。
图5 利用NA-FMEMD算法对实际信号进行分解的结果
C3(大脑对侧感觉运动区)、C4(大脑同侧感觉运动区)采集到的脑电信号数据,经分解后得到的μ节律(d 2 信号)如图6所示,其中横坐标为采样时间,纵坐标为脑电信号幅值。从2 s处开始,受试者开始右手运动想象任务,C3的μ节律幅值明显小于C4的 μ节律幅值,体现了ERD现象与ERS现象。
此外,笔者统计了 NA-FMEMD 算法和NA-MEMD算法对受试者g的200组数据分解所需的平均运行时间。其中,前者所需的分解运行时间为 1.26 s,后者所需的分解运行时间为 7.86 s, NA-FMEMD算法的计算效率比NA-MEMD算法提升了5.24倍。 将选取出的新的信号矩阵通过 CSP 空间滤波器,根据给定 m 的取值提取出脑电信号的特征向量。从200组实验数据中随机选取140组作为训练集,其余60组作为测试集。将训练集通过CSP空间滤波器得到的特征向量分别输入 SVM 分类器中进行分类,计算不同m值(m=1,2,3,4,5)得到的模型在测试集上的分类准确率,通过比较分类准确率来反向确定最优m值。需要说明的是,本次实验仅对数据集进行一次划分,重复实验的本质是对m的取值进行遍历,以最大化分类准确率。除此之外,本实验的重点在于说明本文提出的方法在计算效率和提取有用信息方面的优势,因此弱化了不同分类器对现有实验框架造成的影响,关于不同分类器对分类准确率的影响会在未来的工作中进行分析。本次分类实验结果见表5,每种方法对于每个受试者的最优分类准确率已加粗显示。
从每个受试者最优的分类准确率来看,本文提出的方法对于每个受试者的分类性能都优于NE-MEMD算法。本文提出的方法的平均分类准确率(4 个受试者最高准确率的均值)为85.4%,比基于NA-MEMD算法的方法提高了2.1%。另外,可以观察到每种方法中 CSP 的最优特征维度 m 都是不固定的。这是由脑电信号数据本身的性质决定的,有些数据的特征信息在数据矩阵中相对靠前(比如受试者b的实验数据),有些则分布较散(比如受试者 g 的实验数据)。值得注意的是,在受试者 f 的实验中, NA-FMEMD算法与NA-MEMD算法都收获了一致的最高准确率,但是在不同的特征维度下,NA-FMEMD算法性能普遍比NA-MEMD算法更加优越。 为了进一步验证基于NA-FMEMD的方法在特征提取上的优势,笔者将 NA-FMEMD 算法和NA-MEMD算法分别应用于受试者a和b,并且利用不同的分类器对所提取的特征进行分类。这些分类器分别是SVM、分类回归树(classification and regression tree,CART)、K-最近邻(K-Nearest Neighbor,KNN)和快速线性判别分析(fast linear discriminant,FLD)。通过对特征维度m的遍历,选择最优的分类准确率作为分类结果,不同特征维度下的脑电信号数据集分类准确率见表6。
由表6 可以看出,相较于 NA-MEMD 算法, NA-FMEMD 算法仍收获了最高的分类准确率。对于受试者a来说,NA-FMEMD+FLD获得了93.3%的分类准确率,比NA-MEMD+SVM高出6.6%;对于受试者b来说,NA-FMEMD+FLD获得了86.7%的分类准确率,比NA-MEMD+FLD高出3.4%。除了受试者a中的CART和受试者b中的KNN,在同一个分类器下, NA-FMEMD算法的分类性能基本上高于NA-MEMD算法。另外,从分类器的角度来看,FLD分类器比原先使用的SVM分类器更加适用于基于FMEMD的脑电分类方法。在基于FMEMD的方法框架下,对于受试者a的数据集来说,FLD获得的准确率比SVM高出3.3%,而对于受试者b来说,FLD比SVM高出6.7%。 本节针对受试者g的数据集,根据前述的分类框架进行了另一个实验。笔者打乱选取顺序,随机将受试者g的数据集中的140组数据作为训练集, 60组作为测试集,根据第3.2.2节实验得出的最优m 值(m=5),分别针对基于 NA-FMEMD 与基于NA-MEMD的运动想象脑电信号分类方法重复50次分类实验,利用箱形图统计分析 50 次实验的分类准确率,如图7所示。由图7可知,本文所提方法的分类准确率总体位于基于 NA-MEMD 的方法之上,且异常值少,算法鲁棒性更强。 本文提出了新的基于 NA-FMEMD 的运动想象脑电信号分类方法。相较于 NA-MEMD 算法, NA-FMEMD 算法能够更加准确地估计分解子信号的频率信息(见表2、表3、表4),保留更精确的信号时域特征(见图3、图4、图5),同时拥有更低的运算复杂度。借助这些优势,基于 NA-FMEMD 的方法在处理运动想象脑电信号时显示出更加优越的时频分析性能(见图5、表4),更准确地反应了特定节律(β波、μ波)中的ERD和ERS现象(见图6)。因此,本文对不同的受试者数据进行分类,结果表明所提方法能够获得一致的高分类准确率(平均准确率为85.4%),相较于NA-MEMD算法提高了2.1%(见表5)。另外,本文针对受试者g进行了50次分类实验,结果显示基于NA-FMEMD的方法的性能更加稳定,异常值少,更适用于临床分析(见图7)。 值得注意的是,本文选择了 11 个通道参与所有分类实验。实际上,NA-FMEMD 算法同样适用于密集通道的EEG信号。由于NA-FMEMD算法将多通道运算映射到单通道中,因此,在处理密集通道的信号时,相较于NA-MEMD算法,NA-FMEMD算法在运算效率更有优势。其次,单通道运行能够减少算法本身的过分解,减少信号提取畸变,且增加算法鲁棒性,因此NA-FMEMD处理密集通道时同样能够保持有用信息提取的高效性与精准性,从而获得一致的分类性能。未来的工作将集中于NA-FMEMD在硬件上的实现,以及不同分类器在该框架下的运行效果等。
本文基于NA-FMEMD算法提出了一种新的脑电信号特征提取方法,并将其应用于运动想象任务的分类问题。NA-FMEMD 提高了以往分解算法的运算速度,能够进行更加准确的成分评估。本文以分类准确率为优化目标,指导特征维度的选择,一定程度上最大化了当前算法的分类性能。采用仿真数据对算法进行测试,结果表明,相较于NA-MEMD算法,NA-FMEMD算法在脑电运动想象有用信息提取方面运算效率更高,且能够更加准确地刻画信号瞬时变化信息。最后,结合CSP方法,对BCI Competition IV数据集进行实际测试,并与基于 NA-MEMD 的方法进行比较,结果显示本文所提方法可以获得更好的分类效果。本文提出的方法分类结果更加准确,计算效率高,更符合便携式设备对实时性的要求。下一步工作是提高方法的性能,并且结合嵌入式设备完成在线检测与分类。 作者简介 About authors
郑潜(1995-),男,浙江大学智能系统与控制研究所博士生,主要研究方向为生物医疗信号处理、时频分析、机器学习等 。 乔丹(1995-),女,浙江大学智能系统与控制研究所硕士生,主要研究方向为生物医疗信号处理、时频分析、因果分析等 。 郎恂(1994-),男,云南大学信息学院电子工程系博士后,主要研究方向为工业过程性能评估、信号处理、时频分析等 。 谢磊(1979-),男,浙江大学智能系统与控制研究所教授、博士生导师,主要研究方向为先进控制理论与应用、工业数据挖掘与人工智能、工业控制系统性能评估等。 李东流(1966-),男,高级工程师,森赫电梯股份有限公司董事长,浙江省电梯行业协会会长 。 王琪冰(1978-),男,博士,正高级工程师,IETFellow,森赫电梯股份有限公司研究院院长,英国英中智能机电装备与机器人联合研究中心首席科学家 。 苏宏业(1969-),男,浙江大学智能系统与控制研究所教授、博士生导师,主要研究方向为工业智能与优化控制、机器人与智能无人系统等 。
转载本文请联系原作者获取授权,同时请注明本文来自王晓科学网博客。 链接地址: https://wap.sciencenet.cn/blog-951291-1266421.html
上一篇:
[转载]基于无人机的边缘智能计算研究综述 下一篇:
车载双目摄像头,为什么特斯拉还在迟疑?