基于变分模态分解与深度信念网络的运动想象分类识别研究

何 群，杜 硕，王煜文，陈晓玲，谢 平

(燕山大学 河北省测试计量技术与仪器重点实验室, 河北 秦皇岛 066004)

1 引 言

脑机接口(brain-computer interface, BCI)能够在人脑和电脑或者是外部设备之间建立直接通路；借此，人无需身体上的动作就可以将命令和控制指令直接发送给外部设备[1]。运动想象(motor imagery, MI)分类研究是BCI技术的重要组成部分; 但由于MI脑电信号的非平稳、非线性以及低信噪比的特点[2,3]，其特征提取与分类识别的效果均不稳定，传统人工确定最优时段及最优频段的方法又难免造成信息遗漏进而导致识别率的降低；因此基于脑电信号的MI分类研究成为难点问题。

迄今为止，研究者已经基于受试者进行MI时在mu频带以及beta频带所发生的事件相关去同步/事件相关同步(event-related desynchronization/event-related synchronization, ERD/ERS)[4]现象, 提出了利用时-频特征、共空间模式(common spatial pattern, CSP)以及非线性动力学特征在内的多种特征提取方法[5～7]用于MI脑电信号的特征提取； 同时也将支持向量机(support vector machine, SVM)、线性判别分析(linear discriminant analysis, LDA)和logistic分类器[8～10]等多种机器学习方法用于MI分类中。 Bashar S K等[11]利用多元经验模态分解(multivariate empirical decomposition, MEMD)将多通道脑电信号分解后，选取最优模态进行短时傅里叶变换(short-time Fourier transform, STFT)，将频谱峰值作为MI特征输入k-最近邻(knearest neighbor, KNN)分类器对MI信号进行分类。李明爱[12]等利用主成分分析(principal component analysis, PCA)进行特征提取，以IGHSOM神经网络为分类器识别想象模式。此外，有学者将深度网络运用于MI分类问题。Ren Y F[13]等利用CSP算法得到特征，并以此作为输入来训练一个4层的深度信念网络(deep belief network, DBN)，得到MI的分类结果。 Lu N等[14]获得脑电信号的频域表示, 再以此作为输入来训练一个4层的频率深度信念网络(frequency deep belief network, FDBN)得到MI的分类结果。Tabar Y R等人[15]使用卷积神经网络(convolution neural network, CNN)提取特征后使用堆叠自编码器(stacked autoencoders, SAE)进行MI分类。针对MI分类问题，现有的方法已经取得了一定的效果，但脑电信号的识别率依然是限制MI脑机接口发展的瓶颈；其根源在于现有特征提取及模式识别方法难以兼顾受试者个体差异性，同时手动提取最优想象时段以及最优频段的方法难以避免信息遗漏导致的识别率降低问题。

变分模态分解(variational mode decomposition, VMD)[16]可将脑电信号分解为多个窄带分量，各分量包含了信号不同时-频尺度上的特征，同时解决了模态混叠的问题。深度信念网络(deep brief network, DBN)具有自动学习有效特征的能力[17]，从而达到特征降维的目的，已有研究者基于这一特性将其应用到睡眠分期、情绪识别和癫痫诊断等脑电信号处理领域[18]。

DBN能够对高维特征降维得到MI的最优特征，从而避免手动确定最优想象时段及最优频段导致的识别率降低问题。因此，本文提出使用VMD对脑电信号分解；利用希尔伯特变换提取边际谱(marginal spectrum, MS)、特征频带下的瞬时能谱(instantaneous energy spectrum, IES)以及时-频联合特征(joint time-frequency features, JT-F)来涵盖完整MI时程下的时频域特征，并融合3种特征；引入DBN对融合后的高维特征降维得到MI的最优特征并实现MI脑电模式的识别，避免了手动确定最优时段及最优频段导致的识别率降低问题。

2 基于VMD的高维时频特征提取

由于脑电信号具有非线性、非平稳和低信噪比的特点，经验模态分解(empirical mode decomposition, EMD)及局部均值分解(local mean decomposition, LMD)等自适应分解方法被广泛应用于脑电信号特征的提取[19, 20]。

但现有的自适应分解框架存在多种缺陷[16]：(1)缺乏数学理论支撑；(2)多数分解方法基于递归筛选，误差因多次递归被放大，进而出现模态混叠现象；(3)缺乏噪声鲁棒性。

VMD是一种新的自适应分解方法，与EMD及LMD等递归式分解方式不同；VMD假设每个模态为带宽有限，中心频率不同的窄带分量，将模态的估计问题转化为求解变分问题，改善了EMD及LMD等方法的模态混叠问题。

VMD算法中，定义本征模态分解函数(intrinsic mode function，IMF)为调幅-调频信号，记为：

uk(t)=Ak(t)cos(φk(t))

(1)

(2)

式中：{u1,…,uK}及{ω1,…,ωK}分别为分解得到的K个模态分量及每个模态对应的中心频率。

式(2)引入增广拉格朗日函数得到式(3)：

(3)

式中：α为惩罚参数；λ为拉格朗日乘子。

(4)

采用Parseval/Plancherel傅里叶等距变换将式(4)转化到频域中，得到各模态分量及中心频率的的频域更新分别为：

(5)

(6)

同时利用式(7)更新：

(7)

本文首先对MI脑电信号进行8～30 Hz带通滤波，然后对滤波后的脑电信号进行VMD得到K个IMF分量。

为提取运动想象脑电信号的时频特征，对分解得到的各IMF分量进行希尔伯特变换：

(8)

Zk(t)=uk(t)+jYk(t)=ak(t)ejθk(t)

(9)

式中： PV是柯西主值；t′为时域中的某一时刻；ak(t)是瞬时幅值;θk(t)是瞬时相位。则瞬时频率定义为：

(10)

则原始脑电信号可以表示为各IMF分量的和值形式：

(11)

定义希尔伯特谱为：

(12)

本文按照式(13)～式(15)分别求取希尔伯特边际谱(MS)，特征频带下的瞬时能谱(IES)及时-频联合特征(JT-F)。

(13)

(14)

(15)

对于式(13)，为避免人工选择造成的时、频域信息遗漏，将式中的T定为运动想象结束时刻，由于C3与C4通道分别分布在大脑的左侧和右侧，CZ通道则分布在大脑的正中央，所以受试者在进行左右手MI实验时，仅C3与C4通道脑电信号能够体现出ERS/ERD现象。因此，本文在进行特征提取时，不对CZ通道脑电信号进行处理, 可以得到完整MI时程下C3、C4通道的边际谱特征F1：

F1={MSC3,MSC4}

(16)

式(14)中的ω1及ω2分别对应MI的mu及beta特征频带，分别取两组对应值，其中mu频带下的ωmu1与ωmu2分别为8 Hz和13 Hz，beta频带下的ωbeta1与ωbeta2分别为16 Hz和28 Hz，则得到两个特征频带下的瞬时能量特征IESmu及IESbeta，至此得到特征频带下C3、C4通道的瞬时能量特征F2为：

F2={IESmuC3,IESmuC4,IESbetaC3,IESbetaC4}

(17)

计算JT-F特征，则首先对时间及频率分别加窗，定义twindow1及twindow2分别为时间窗的起止时间，ωwindow1及ωwindow2分别为频率窗的起止频率。本文定义0.25 s时间窗将完整MI时程分割为n个时段，以2 Hz频率窗将8～30 Hz的特征频带分割为11个子频带，然后以式(15)计算第i(≤n)个时段第j(j≤11)个频段下的JT-Fi,j值，则JT-F特征可以表示为：

J={JT-F1,1,JT-F2,1,…JT-Fn,1,…JT-F1,11,

JT-F2,11,…JT-Fn,11}

(18)

则C3、C4通道的时频联合特征F3为：

F3={JC3,JC4}

(19)

融合3种特征形成高维特征向量F={F1,F2,F3}作为MI特征。

3 基于DBN的高维特征分类

DBN由底层的若干个限制玻尔兹曼机(restricted boltzmann machines，RBMS)和顶层的分类器组成，相邻2层即为一个RBM,各层内部相互独立，但数据可通过特定的激活函数按照RBM的学习规则在层间相互转换，图1为5层DBN网络结构图。

图1 DBN结构图Fig.1 The structure of DBN

图1中，每个RBM由可视层(v)和隐含层(h)组成，层与层之间通过权重w彼此互联，但层内无连接，假设v层包括n个可见单元，h层包括m个隐含单元。则1组确定状态的RBM系统的能量为：

(20)

式中：vi为第i个可见单元状态；hj为第j个隐含单元状态；对∀i,j,vi∈{0,1},hj∈{0,1}，θ={wij,aj,bj}为RBM模型的参数；wij为可见单元i与隐含单元j间的连接权重；ai为可见单元i的偏置；bj为隐含单元j的偏置。定义可视层和隐含层的联合概率密度为：

P(v,h;θ)=e-E(v,h;θ)/Z(θ)

(21)

(22)

由式(22)，隐含节点与可视节点被激活的概率为：

(23)

式中:σ为激活函数。

为使RBM良好拟合训练样本，将训练目标定为最大化似然函数：

(24)

Hinton利用对比散度(contrastive divergence，CD)方法，只需少数次状态转移就可以得到对可视层的较好估计。该算法首先计算所有隐含层单元条件概率，并使用Gibbs抽样确定隐含单元状态；再计算可视层状态，完成可视层的重构，并求解式(24)。

最后得出权重、偏置的更新过程为：

(25)

式中:〈·〉data为原始数据的期望；〈·〉recon为重构数据的期望；ε为学习率。针对参数ε，Hinton提出学习动量项参数Momentum，使收敛更为精准，Momentum定义为：

(26)

式中:ρ为动量学习率。

则权重偏置的更新法则变为：

(27)

在图1所示的DBN中，上一个RBM的隐含层即为下一个RBM的可视层，上一个RBM的输出即为下一个RBM的输入，输入数据经过低层RBM学习后的输出结果作为高一层RBM的输入，直至最后一层；因此高层RBM能够学习到更抽象，更稀疏的特征表示。

以上过程为前向的无监督学习过程，利用数据的标签，采用反向传播方法从网络的softmax分类器开始对网络参数进行反向微调，可以进一步优化DBN网络参数，避免网络陷入局部最优。

针对运动想象高维细节特征F，本文构建4层DBN模型对其进行特征降维得到运动想象的最优特征，进而实现对受试者左右手MI意图的识别。至此，本文完成了如图2所示的基于VMD与DBN的脑电信号特征提取与识别框架。

图2中，利用DBN中的4层网络结构对高维特征F降维得到其稀疏表达，降低MI高维特征中的噪声、非线性和非相关信息，达到特征降维的目的，从而达到更好的分类效果。

图2 基于VMD与DBN的特征提取与识别框架Fig.2 The framework of feature extraction and recognition based on VMD and DBN

4 实验研究与结果分析

4.1 受试对象与实验流程

为验证本文提出的方法的普适性及有效性，实验采用BCI 2003竞赛Dataset Ⅲ数据集[21](S1)、BCI 2005竞赛Dataset Ⅲb数据集[22](S2～S4)以及BCI 2008竞赛Dataset IIb数据集[23](S5～S7)。

各数据集信息如表1所示。

表1 数据集信息Tab.1 The information of datasets

信号采集过程中，由于采集设备原因，个别时刻信号未被成功采集，这些数据点在信号中表现为“NAN”点，为提高信号信噪比，实验中对数据进行预处理工作，对数据中的“NAN”数据点置零，同时为减小数据处理的工作量，对受试者S5～S7对应数据进行了降采样处理，将采样频率降为125 Hz。由于运动想象的ERD/ERS现象主要出现在mu节律(8～13 Hz)与beta节律(16～28 Hz)，因此对原始信号进行8～30 Hz带通滤波，所用滤波器为6阶巴特沃斯滤波器，设置阻带截止频率分别为6 Hz和32 Hz。

实验中将VMD分解模态数K定为7，选择sigmoid函数作为DBN网络的激活函数，隐藏层神经元个数定为[400，200]，动量学习率定为0.1，批处理容量为全体训练集，各受试者学习率、预训练次数及微调次数等信息如表2所示。

4.2 实验结果与分析

4.2.1 特征组合及对比

为比较不同特征提取方法对识别率的影响，实验设计3个参照组与所提方法进行对比，分别为：

(1)采用4层DBN网络对边际谱特征F1降维并识别MI脑电信号；(2)采用4层DBN网络对瞬时能量谱特征F2降维并识别MI脑电信号；(3)采用4层DBN网络对时-频联合特征F3降维并识别MI脑电信号。表3呈现了各受试者在3个参照组与本文所提方法下的最大识别率。

表2 各受试者参数信息Tab.2 The parameters information of subjects

表3 不同特征下的最大识别率Tab.3 The maximum recognition rates of different features (%)

由表3中发现，采用边际谱特征所得平均最大识别率为60.51%，采用瞬时能量谱特征所得平均最大识别率为73.03%，采用时-频联合特征所得平均最大识别率为77.59%，采用融合特征所得平均最大识别率为79.00%。

其中，边际谱特征取得的分类效果相对较差，这是因为实验范式的设定导致MI信号强弱随时间变化，完整MI时程下的边际谱特征的时域分辨率较低，不能准确表征MI，但边际谱特征的频域分辨率高，对最后的识别结果有辅助作用；与边际谱特征相反，特征频带下的瞬时能谱则具有优秀的时域分辨率和较差的频域分辨率，取得了更加良好的识别效果；由于兼顾了时、频分辨率，3种单一特征中，时-频联合特征取得了最优的识别效果；而将3种特征融合后，由于时、频分辨率均较高，MI的识别率得到提升，同时不同受试者识别率间的标准差降低，分类框架的稳定性得到了进一步提升。

4.2.2 不同分解方法对比

实验对预处理后的脑电信号进行VMD分解，S1的C3通道脑电信号VMD分解结果及其分量频谱图如图3。

实验中作为对比，将脑电信号进行EMD分解，S1的C3通道脑电信号EMD分解结果及分量频谱图如图4。

在图3及图4中，由于对脑电信号预先进行了滤波处理，自适应分解的IMF分量频率集中在8～30 Hz范围内。由图3及图4可以看出，EMD分解得到的分量频带混叠现象严重，而VMD的分解结果均为窄带分量，不会出现相同成分信息出现在不同分量的情况。实验中发现EMD分解所得结果，特征频带主要分布在前3阶分量中。

图3 VMD分解结果及各分量频谱图Fig.3 The decomposition results of VMDand spectral diagrams of components

图4 EMD分解结果及各分量频谱图Fig.4 The decomposition results of EMD and spectral diagrams of components

本文选取EMD分解的前3阶分量，提取相同高维特征，利用DBN对2种方法所得特征进行特征降维与分类，分别对2种分解方法所得识别率进行比较，识别结果如表4所示。表4中可以看出，由于更优的分解效果，相比于用EMD提取的特征，VMD所提取的特征具有更佳的识别率。

表4 EMD与VMD所得识别率比较Tab.4 Comparison of recognition rates between EMD and VMD (%)

4.2.3 不同方法对比

基于相同的脑电数据，将本文所述方法与其它方法进行比对，表5及表6给出了比较结果。

表5中，空白处表示文献[13]未采用此位受试者的MI脑电数据验证所提方法。表5中，文献[24]选用功率谱密度、Hjorth参数、AR模型系数和小波系数等多种特征对运动想象分类。文献[25]利用总体经验模态分解(ensemble empirical mode decomposition, EEMD)提取边际谱、瞬时能谱及近似熵作为特征输入线性判别分类器识别运动想象脑电信号。

表6中，文献[26]采用滤波器组共空间模式(filter bank common spatial pattem, FBCSP)提取MI特征，利用多分类器降低训练集与测试集之间的分布差异，提升了MI识别率；文献[27]则利用新型相空间重构算法计算非线性动力学特征，选用集成分类器识别想象类别。此类人工确定MI最优时段及最优频段的方法难免遗漏细节信息造成识别率的降低，本文利用DBN网络对高维时、频特征降维得到最优特征并识别MI模式，明显提高了分类准确率，且S2、S4及S7的识别率提升明显。

文献[13]同样将深度框架应用于MI问题中，利用CSP算法提取特征，并构建2层DBN网络对特征降维并识别脑电模式，从表5的比较可以看出本文所构建的方法在特征提取与识别框架方面优于文献[13]所提的方法。

表5 S1～S4所得识别率比较Tab.5 The comparison of recognition rates of S1～S4 (%)

表6 S5～S7所得识别率比较Tab.6 The comparison of recognition rates of S5～S7 (%)

BCI 2003竞赛Dataset Ⅲ数据集中，组织者将互信息与识别率同时定为判断指标；BCI 2005竞赛Dataset Ⅲb数据集中，组织者将互信息作为判断指标；BCI 2008竞赛Dataset IIb数据集中，组织者将Kappa值与识别率作为判断指标；实验中将本文方法所得分类结果与BCI竞赛冠军指标进行对应比较，比较结果见表7。表7中空白处表示此届比赛组织者未采取此项指标作为评判标准。

表7 BCI竞赛冠军结果与所提方法结果比较Tab.7 Comparison of results between the BCI competition’s winner and the proposed method

从表7中可以看出：对于S2、S4、S5～S7，本文所述方法所得指标均优于BCI冠军所得结果；对于S1，本文方法所得互信息低于BCI冠军结果，但所得识别率优于BCI冠军结果。

5 结 论

本文应用VMD方法对脑电信号分解，利用希尔伯特变换提取边际谱、特征频带下的瞬时能谱和时-频联合特征表征MI，引入4层DBN网络对融合后的高维特征降维得到最优特征并分类，避免了人工确定最优时段及最优频段导致的信息遗漏，提升了MI脑电信号的识别率。运用所述方法采用脑机接口竞赛数据集进行分类实验，平均准确率达到79%。本文方法有效提升了MI脑电信号的识别率，为MI脑机接口的应用奠定了基础。