在空间站中,由于失重,航天员不能像在地面一样自如地控制机器完成指定操作,而脑机接口技术可以解码脑电,利用脑电信号(EEG signal)控制空间机械臂运动。脑机接口(brain computer interface,BCI)[1]是一种不依赖于外围神经和肌肉组织的大脑正常输出通路来实现大脑与外部环境之间的直接通信的技术。通过脑机接口可以将操作者的意图传递给机器,达到人脑控制机器运动的效果。
脑电信号的特征提取方法主要有时域分析法、频域分析法、时频分析法、空域分析法。张绍荣等[2]提出基于稀疏贝叶斯逻辑回归运动想象脑电信号分类模型,提取方差、峰度、峭度、6阶自回归模型(auto regression,AR)系数、通带频率以及小波能量作为特征。Mahmood等[3]采用共空间模式方法进行特征提取,用支持向量机(support vector machine,SVM)分类。Park等[4]采用特征向量中心性特征选择方法,通过小波包分解、共空间模式提取特征,并采用内核极限学习机分类。Sun等[5]采用共空间模式提取特征,将概率神经网络(probabilistic neural network,PNN)分成6类脑电信号。
由于运动想象信号主要表现为特定频率的波段信号幅值变化,只采用时域特征无法体现,因此,常用时频分析法和空域分析法进行分析。基于快速傅里叶变换的功率谱密度(power spectral density,PSD)表示信号功率在各频率点的分布情况,该方法简单易行,但事件和频率不能同时满足最佳效果[6]。小波包分解(wavelet packet decomposition,WPD)是进行信号时频分析和处理的常用工具,通过伸缩平移对信号进行多尺度细化,达到高频处时间细分、低频处频率细分的效果,适用于非平稳信号处理[7]。共空间模式(common space pattern,CSP)是一种有效提取运动想象脑电信号特征的方法,于1999年首次用于运动想象信号分类[8],目前在脑机接口系统中已得到广泛应用[9]。通过共空间模式计算空间滤波器,使不同类别的脑电信号特征差异最大化。
以上所提及的相关文献中分别提取了运动想象信号的时域特征、频域特征、空域特征,但没有将3类特征结合起来,且未考虑由于受试者个人条件不同选择不同特征的情况。本文针对小波包分解脑电信号,采用递归特性消除方法(recursive feature elimination,RFE)选择特征,将经过共空间模式滤波得到的特征向量输入长短期记忆网络进行分类,识别4类脑电信号,平均识别准确率达93.69%。将识别的脑电信号作为控制信号,控制虚拟环境中的空间机械臂的中臂和大臂顺时针或逆时针运动。
采用基于负熵最大的FastICA方法对脑电信号进行盲源分离。由于相邻通道采集的信号之间相互干扰,且这些信号相互独立,因此,常用独立成分分析(independent component analysis,ICA)的方法对脑电信号进行分离[10]。
基于神经生理经验知识,选择大脑皮层区域感知区22通道采集的脑电信号。运动想象信号采集通道如图1所示。
图1 运动想象信号采集通道分布
Figure 1 Distribution of electrodes for MI
以C3、C4、Cz通道信号为例,未经独立成分分析方法分离的运动想象信号如图2所示。经过FastICA盲源分离后得到C3、C4、Cz独立成分如图3所示。
图2 C3、C4、Cz通道脑电信号
Figure 2 EEG signals of C3,C4,Cz channels
图3 经过FastICA盲源分离后C3、C4、Cz通道脑电信号
Figure 3 EEG signals of C3,C4,Cz channels after blind source separation by FastICA
小波包中提供多种类型的小波基。对跳变较多的信号,Haar-Wsalsh基比较适用;局部三角函数基适用于含正弦分量或高频振荡的信号;对于运动想象脑电信号这类由分段多项式结构组成的离散型的信号,常采用Daubechies小波基[11]。本文采用Db3小波(滤波器长度为3)对4类运动想象信号进行3层分解。
小波包分解原理如下[12]:
(1)
式中:s为尺度因子;b为平移因子;f(t)为源信号;ψ(·)为小波基函数。离散表达式为
(2)
式中:cl为逼近系数;gj,l为细节系数;φ(·)为尺度函数。
cl=x(t)φ(t-l);
(3)
gj,l=x(t)2-j/2ψ(2jt-l)dt;
(4)
(5)
(6)
式中:h0为低通滤波器;h1为高通滤波器。
小波包分解如图4所示。本文中数据集采集运动想象信号经过Butterworth滤波后的频率为0.5~100 Hz,因此第3层小波包节点对应的频率如表1所示。
图4 小波包3层分解图
Figure 4 Three-level sub-band tree of WPD
表1 第3层小波包节点对应的频率
Table 1 Frequency domain corresponding to nodes of level 3
节点频率/Hzaaa0.5^12.5aad12.5^25.0add25.0^37.5ada37.5^50.0dda50.0^62.5ddd62.5^75.0dad75.0^87.5daa87.5^100.0
当进行运动想象时,不同的受试者子频带变化不同。针对不同的受试者选择不同的通道和频带信号可以去除每个受试者不重要的子频带,使剩余子频带有更好的可分性,可提高识别准确率[13]。特征选择具体步骤如下。
步骤1 计算第3层中每个节点的功率值:
(7)
式中:l表示第l个通道;j表示第j层;n表示第n个节点;r表示采样点个数。
步骤2 对于每个样本,形成一个新的向量其中,N=n×l=22×8=176。
步骤3 通过递归特性消除(RFE)的方法选择特征。每轮模型训练时,根据特征重要性排序,删除一个最不重要的特征,共训练10轮,去掉10个最不重要的特征,得到166个子频带
步骤4 将被选入V′中子频带的系数重新构成一个新的矩阵E165×q,q为每个节点系数个数。
二分类共空间模式的原理为求取空间滤波器使得一类信号的投影矩阵方差最大,另一类信号的投影矩阵方差最小[14],算法如下。
第k类脑电信号Xk可表示为
(8)
式中:N为脑电仪的通道数;T为总采样点数。
求取脑电信号的规范化协方差矩阵Rk:
(9)
式中:表示对矩阵
取迹。设2类任务的规范化协方差矩阵分别为R1和R2。
复合的规范化协方差矩阵RC:
RC=R1+R2。
(10)
对规范化协方差矩阵RC进行特征值分解:
(11)
式中:λC为特征值矩阵;UC为特征向量法矩阵。
构造白化变换矩阵Q:
(12)
利用白化变换矩阵Q对协方差矩阵R1和R2进行处理:
S1=QR1QTS2=QR2QT。
(13)
式中:S1、S2具有相同的特征向量。构造对角矩阵λ1、λ2和特征向量矩阵B,对S1、S2进行主分量分解:
S1=Bλ1BTS2=Bλ2BT,λ1+λ2=I。
(14)
若λ1中特征值按降序排列,则λ2中特征值按升序排列,即当S1特征值最大时,S2特征值最小,保证了2类信号差异最大化。
求取投影矩阵,即空间滤波器M:
M=BTQ。
(15)
通过空间滤波器M得到第k类脑电信号特征矩阵Zm×T:
Zm×T=MT×XN×T。
(16)
式中:N为脑电仪的通道数;T为总采样点数;m为生成空间滤波器时选取特征个数。本文中m=30,计算平均频带功率:
(17)
将4类运动想象信号两两组合,作为CSP滤波器的输入,得到6个空间滤波器。将2类信号组合,通过空间滤波器的方法增大信号与另3类信号的差异性,使4类信号具有可分性,为后续分类提供良好的基础。CSP提取4类脑电特征如图5所示。
图5 CSP提取脑电特征示意图
Figure 5 CSP extraction of EEG detail coefficients
本文选择长短期记忆(long short term memory,LSTM)网络进行分类。长短期记忆网络是一种改进的循环神经网络算法,它使用一种被称为 LSTM的记忆单元来判别哪些信息应该被保留,控制信息从前一时刻到下一时刻进行传输,是目前应用最为广泛的具有记忆功能的网络[15]。
每个空间滤波器选取30维特征,信号与共空间模式滤波器相乘,计算平均频带功率,6个空间滤波器可组成30×6的特征矩阵。选择长短期记忆网络搭建分类模型。LSTM的记忆单元共100个神经元,输出层选择Softmax作为激活函数,得到4类分类结果,共迭代30次。训练模型损失函数选择交叉熵函数,优化器选择Adam优化器。
本文介绍了用于验证的公开数据集,并将本文方法在数据集上测试得到的识别准确率与其他特征提取和分类方法比较。
本文选用的数据集为奥地利大学提供的第4届脑-机接口大赛的data set A。该数据集中包含9位受试者4类运动想象的数据。每个受试者进行288次单次实验,其中每类运动想象各72次。每个实验持续时间约为8 s。在第2 s时,屏幕出现上、下、左、右箭头,受试者分别进行舌头运动想象、双脚运动想象、左手运动想象、右手运动想象。3~6 s时,受试者保持该动作想象。信号采样频率为250 Hz,并且通过0.5~100 Hz的带通滤波器和50 Hz的陷波滤波。实验共采集25个通道信号,包括22个脑电通道和3个眼电(electro-oculogram,EOG)通道[16]。
选择前3位受试者的3~6 s的22个脑电采集通道的训练数据随机划分,其中75%作为训练集,25%作为测试集。与其他方法对比结果如表2所示。
表2 不同的特征提取方法和分类方法的识别准确率对比
Table 2 Recognition accuracy comparison of different feature extraction and classification methods %
受试者WPD+CSP+RFE+LSTMWPD+CSP+RFE+LDAWPD+CSP+LSTMCSP+LSTMWPD+LSTM190.2878.1688.8943.0627.78294.2577.0194.4462.5025.00396.5572.2295.8345.8326.35平均93.6975.8093.0550.4626.38
表2中LDA为线性判别分析。由表2可知,在采用WPD和CSP进行特征提取并用LSTM进行分类时,采用RFE选择特征的识别准确率略高于未选择特征的识别准确率。在采用WPD和CSP进行特征提取并用RFE进行特征选择时,LSTM分类效果好于LDA。若未采用CSP提取特征,则几乎无法正确分类。由于采用RFE选择特征时,需要遍历所有特征,且对准确率影响较小,因此当特征数量较大时,不建议采用。本文将不同类别信号两两组合进行CSP空域滤波器计算,即需要计算k(k-1)/2个空域滤波器(其中k为类别数),若类别过多,则需要计算的空域滤波器会大量增加。
将上述分类结果作为控制信号,通过串口在Unity 3D中控制空间机械臂运动,仿真实验流程如图6所示。
图6 仿真实验流程
Figure 6 Simulation experiment process
空间机械臂是由底座(基座)、大臂、中臂、小臂、末端手爪组成,将其安装在空间站上,主要完成空间在轨抓取目标物体任务。航天员利用脑电解码指令控制空间站外的空间机械臂。当脑电信号为左手运动想象时,控制空间机械臂大臂逆时针转动。当脑电信号为右手运动想象时,控制械臂大臂顺时针转动。当脑电信号为舌头运动想象时,控制机械臂中臂顺时针转动。当脑电信号为双脚运动想象时,控制机械臂中臂逆时针转动。各类脑电信号对应控制指令及作用如表3所示,空间机械臂初始位置如图7所示,空间机械臂按指令转动效果如图8所示。
表3 4类脑电信号对应控制指令及机械臂动作
Table 3 Four types of EEG corresponding control instructions and space manipulator actions
脑电信号控制指令机械臂动作左手运动想象0001大臂逆时针转动右手运动想象0010大臂顺时针转动双脚运动想象0100中臂逆时针转动舌头运动想象1000中臂顺时针转动
图7 空间机械臂初始位置
Figure 7 Initial position of the space manipulator
图8 空间机械臂按指令转动后位置
Figure 8 Position of the space manipulator as instructions
本文提出一种基于小波包分解的共空间模式方法提取脑电信号特征,采用递归特性消除方法选择特征,通过长短期记忆网络进行分类。首先,用快速独立成分分析(FastICA)方法对22通道的脑电信号进行盲源分离。然后,用Db3小波分解脑电信号,计算每个节点的功率值,采用RFE选择特征,将选择的子频带和通道组成特征矩阵,通过共空间模式对脑电信号进行空间滤波,并采用长短期记忆神经网络对提取的脑电特征进行分类。将本文方法与其他特征提取方法在data set A公开数据集上进行测试。测试结果表明,本文方法识别率均优于常用的几种方法。最后,采用解码的脑电信号作为控制信号,控制空间机械臂中臂和大臂顺时针、逆时针转动,验证了本文方法可以用于解码控制。
该方法可用于帮助航天员控制空间站外的空间机械臂在轨操作。脑电信号易受周围环境干扰,而且可实现任务较少,人与人之间具有较大差异性,如何将训练好的模型通过迁移学习的方法解决脑机接口系统的局限性[17]是下一步研究的重点。
[1] VAUGHAN T M.Guest editorial brain-computer interface technology:a review of the second international meeting[J].IEEE transactions on neural systems and rehabilitation engineering,2003,11(2):94-109.
[2] 张绍荣,朱志斌,冯宝,等.基于组稀疏贝叶斯逻辑回归运动想象脑电信号分类模型的通道选择与分类新算法[J].仪器仪表学报,2019,40(10):179-191.
[3] MAHMOOD A,ZAINAB R,AHMAD R B,et al.Classification of multi-class motor imagery EEG using four band common spatial pattern[C]//2017 39th Annual International Conference of the IEEE Engineering in Medicine and Biology Society(EMBC).Piscataway:IEEE,2017:1034-1037.
[4] PARK H J,KIM J,MIN B,et al.Motor imagery EEG classification with optimal subset of wavelet based common spatial pattern and kernel extreme learning machine[C]//2017 39th Annual International Conference of the IEEE Engineering in Medicine and Biology Society(EMBC).Piscataway:IEEE,2017:2863-2866.
[5] SUN A Q,FAN B H,JIA C C.Motor imagery EEG-based online control system for upper artificial limb[C]//Proceedings of 2011 International Conference on Transportation,Mechanical,and Electrical Engineering(TMEE).Piscataway:IEEE,2011:1646-1649.
[6] 庞存锁,刘磊,单涛.基于短时分数阶傅里叶变换的时频分析方法[J].电子学报,2014,42(2):347-352.
[7] 王艳景,乔晓艳,李鹏,等.基于小波包熵和支持向量机的运动想象任务分类研究[J].仪器仪表学报,2010,31(12):2729-2735.
[8] RAMOSER H,MULLER-GERKING J,PFURTSCHELLER G.Optimal spatial filtering of single trial EEG during imagined hand movement[J].IEEE transactions on rehabilitation engineering,2000,8(4): 441-446.
[9] 杨帮华,陆文宇,何美燕,等.脑机接口中基于WPD和CSP的特征提取[J].仪器仪表学报,2012,33(11):2560-2565.
[10] BAI X P,WANG X Z,ZHENG S H,et al.The offline feature extraction of four-class motor imagery EEG based on ICA and Wavelet-CSP[C]//Proceedings of the 33rd Chinese Control Conference.Piscataway:IEEE,2014:7189-7194.
[11] 乔晓艳,王春晖,任兆麟.小波统计方法提取想象运动诱发脑电特征[J].测试技术学报,2013,27(3):201-207.
[12] 张焕,乔晓艳,王飞.基于小波-共空间模式的多类想象运动诱发脑电特征提取[J].测试技术学报,2015,29(6):523-528.
[13] WANG H X.Optimizing spatial filters for single-trial EEG classification via a discriminant extension to CSP:the fisher criterion[J].Medical & biological engineering & computing,2011,49(9):997-1001.
[14] YANG B H,LI H R,WANG Q.Subject-based feature extraction by using fisher WPD-CSP in brain computer interfaces[J].Computer methods and programs in biomedicine,2016,129:21-28.
[15] 王丙琛,司怀伟,谭国真.基于深度强化学习的自动驾驶车控制算法研究[J].郑州大学学报(工学版),2020,41(4):41-45,80.
[16] 曹胜海,彭尧,张建平,等.基于CSP和RR的多类运动想象脑电信号的识别分类研究[J].软件,2017,38(12):223-228.
[17] 曾焕生.基于迁移学习的运动想象脑电分类方法研究[D].广州:华南理工大学,2019.