基于双树复小波变换和GBDT的运动想象脑电识别

(东南大学 仪器科学与工程学院,江苏 南京 210096)

脑电信号(Electroencephalogram,EEG)是人类通过大脑在进行想象活动时，大量的脑神经细胞主动参与到脑部信息处理过程中产生的一种生物电信号[1]。脑机接口(Brain-Computer Interface,BCI)就是利用大脑产生的脑电信号不经过外周神经和肌肉组织直接与计算机等外围控制设备建立起的信息交流和控制通道[2]。近年来，随着人口老龄化趋势的上升，以上肢功能运动障碍最为常见的偏瘫患者日益增多，基于脑机接口的康复医疗技术现已成为康复领域中的研究热点[3]。

BCI系统如何快速准确地识别出EEG信号所表达的信息，已经成为衡量一个BCI系统是否优良的重要指标，其中，特征提取和特征分类是BCI研究方向上的重要技术。在特征提取方面，由于EEG信号是一种非平稳随机信号，且具有时变敏感性和个体差异性大等特点，因此优良的特征提取方法对BCI系统至为关键。目前针对EEG信号常用的特征提取方法有：① 针对单一的时频域特性的方法[4]，由于易受噪声的影响而造成模式的识别率较低。② 自回归(Auto Regression,AR)模型[5]和自适应自回归(Adaptive Auto Regression,AAR)模型[6]。AR模型主要针对平稳的随机信号，在处理非平稳的EEG信号时具有其局限性；AAR模型虽然弥补了AR模型在短信号分析时的局限性，但其更适合平稳的随机信号。③ 快速傅里叶变换(Fast Fourier Transformation,FFT)[7]，其在处理平稳信号中起到了重要的作用，但在像EEG信号这种非平稳随机信号的处理中受到了许多局限。④ 小波变换[8]。该方法是一种典型的时频分析算法，具有多分辨率特性，在频域和时域均保持良好的分辨率，非常适用于处理EEG等非平稳信号。

在特征分类方面，目前常采用线性分类LDA[9]、支持向量机SVM[10]等的分类方法。线性分类器LDA在进行一些非线性特征分类任务时会存在分类准确度不高的问题。支持向量机SVM通过将低维的特征映射到高维空间中，通过构建超平面的方式从而达到特征可分，这种分类器的分类性能往往受到多种核函数和相关参数的影响较大，且分类器训练时间一般较长。而梯度提升决策树(Gradient Boosting Decision Tree,GBDT)分类算法是一种基于梯度提升和分类回归决策树相结合的分类算法，通过多分类器融合决策的方式可以显著提高分类器的性能，还具有防止过拟合等优点[11]。

由于双树复小波变换具有较好频率抗混叠和抗噪能力的特点，适用于处理具有非平稳时变的EEG信号[12]。本文利用双树复小波变换作为EEG信号特征的提取方法，而后采用4个典型的时间分段作为特征提取的时间段，分别计算相应时间段内信号的能量均值作为GBDT分类器的特征向量进行分类比较。

1 双树复小波变换

双树复小波变换(Dual-Tree Complex Wavelet Transform,DTCWT)的概念最初是由Kingsbury[13]提出，而后Selesnick[14]等又在2005年进一步完善了双树复小波变换的分解与重构算法。相对于小波包变换(Wavelet Packet Transform,WPT)，双树复小波变换通过引入两路离散小波变换的方式对原始信号进行分解和重构，其一为实部树，另一为虚部树，如图1所示。由于虚部树变换相对于实部树变换的具有延时的特点，使得虚部树变换采样恰好能够采集到实部树变换采样时丢失的原始信息，实现信息的有效互补，较为完善地保留原始信号的时频域信息，从而实现双树复小波变换能够较为完整地分解重构出原始信号。

图1 双树复小波变换的分解和重构

双树复小波变换由两路离散小波变换组成，实部树的离散小波变换系数为

(1)

(2)

式中，j=1,2,3,…,J。

虚部树的离散小波变换系数为

(3)

(4)

式中，j=1,2,3,…,J，ψh和ψg为变换的两个实小波，φh和φg为变换的尺度函数。

由此可得DTCWT在各层的实小波变换系数和虚小波变换系数：

(5)

(6)

式中，j=1,2,…,J。

双树复小波变换的细节系数和近似系数重构为

(7)

(8)

经重构后即可得到双树复小波变换在各层下的细节系数和近似系数，若设本文的原始EEG信号的采样频率为fs，则原始信号经过L层分解和重构后的各个分量AL，DL，DL-1，…，D1所对应的频率段分别为[0,fs/2L+1],[fs/2L+1,fs/2L],[fs/2L,fs/2L-1],…,[fs/22,fs/2]。

2 梯度提升决策树(GBDT)

梯度提升(Gradient Boosting)算法是机器学习中的一种算法，常用于处理分类和回归问题。基于梯度提升的分类算法产生的预估模型建立在一个弱分类器的基础上，在进行训练时，通过不断地对弱分类器的残差损失函数进行迭代和优化，其中每次迭代的结果产生一个新的弱分类器，这个弱分类器将对前一个弱分类器在负梯度方向上的残差损失进行优化，减小上次迭代结果产生的近似残差，最终总分类器的结果就是由每轮迭代的弱分类器加权和得到，形成较为理想的预估模型。在实际应用中，弱分类器的选取较为典型的就是分类回归决策树(Classification and Regression Tree,CART)，本文所使用的梯度提升决策树(Gradient Boosting Decision tree，GBDT)分类算法就是基于梯度提升结合分类回归决策树的一种算法[15]。梯度提升决策树分类器算法具体训练步骤可总结如下：

(2) 初始化模型的初始值F0(x)。

(9)

其中，可以通过最小化残差的方式来求取γ参数值。

(3) 经过对模型赋予初始值，接下来就可以对模型进行M次迭代，设迭代次数的变量为m(m=1,2,3,…,M)。则可由下式计算第m次沿负梯度方向的近似残差：

(10)

式中，i=1,2,3,…,n。

(5) 计算hm(x)分类器的权重系数。

(11)

式中，γm为当前样本残差下决策树分类器hm(x)的权重系数。

(6) 更新训练模型，可得输出函数Fm(x)，迭代M次后，即可输出最终的预估分类模型FM(x)。

Fm(x)=Fm-1(x)+γmhm(x)

(12)

训练完成后，求得各个决策树分类器在不同迭代次数下的权重系数，梯度提升分类决策树分类器FM(x)即可对输入的测试集数据进行测试，所得的结果即为不同决策树分类器在不同权重下的最接近于真实标签的结果。

3 实验研究与结果分析

3.1 实验数据来源

本文采用的实验数据是由Graz大学提供的2003年著名的BCI竞赛DataSetⅢ的标准数据集。数据集采集对象是一位25岁的正常女性，该女性会根据屏幕上随机出现的左右手提示进行想象。实验分为7组进行，每组进行40次实验，共280次实验，所有组的实验都会在同一天内完成。每次实验用时9 s，时序图如图2所示。前2 s为静息时间，2 s后会有一个声音提示并且屏幕上出现持续时间为1 s 的‘+’提示符，提示1 s后将开始想象左右手运动。实验进行到3 s时，屏幕上将随机出现一个指向左或右的箭头提示，在3～9 s时间段内，要求被测试者想象左右手运动。

图2 实验时序图

被测试者的脑电信号由C3、CZ和C4三个脑电导联采集，采样频率为128 Hz并经过0.5～30 Hz的粗略去除噪声信号的滤波处理。数据集中包含140个训练集样本和140个测试集样本，其中70次为想象左手运动，70次为想象右手运动，另测试集样本的真实标签已在赛后公布。

3.2 双树复小波变换应用于脑电信号的特征提取

当人体在进行运动想象时，大脑对侧的运动感知皮层就被激活，该区域中的代谢及血流量均增加，神经元活动加强，脑电信号独立性增强，同步化程度降低，进而造成μ节律(8～12 Hz)和β节律(14～35 Hz)能量减少；而同时大脑同侧的运动感知皮层的神经元活动被抑制，进而导致μ节律和β节律能量值升高，这种现象被称为事件相关去同步(Event-Related Desynchronization,ERD)和事件相关同步(Event-Related Synchronization,ERS)。ERD/ERS规律的发现为脑电信号的运动想象分类提供了一种思路。

双树复小波变换相对于小波包分解在频域上更具有抗混叠能力，为了验证这个结论，本文选取了对非平稳脑电信号适应性良好的db5小波作为基波函数，对训练集中的一个原始脑电信号分别做双树复小波和小波包的4层分解重构，由于8～16 Hz包含了μ节律和部分的β节律信息，因此本文对8～16 Hz对应的近似系数D3进行重构。原始脑电信号波形、小波包分解重构的D3波形和双树复小波分解重构的D3波形在时域上的对比图如图3所示，从图中可以看出重构后的D3波形在时域上相近。

图3 原始波形、WPT和DTCWT重构的D3波形对比

为了进行对比分析，本文对原始信号通过数字带通滤波器进行了8～16 Hz的滤波处理，滤波后的原始信号的频谱分别与DTCWT和WPT的分解重构D3波形的频谱进行对比，如图4所示。

图4 DTCWT和WPT重构的D3波形频谱对比

从图中可以看出DTCWT通过分解重构D3波形能够较为完整地复现原始脑电信号在8～16 Hz频段的信息，具有较好的频率抗混叠能力。相反地，WPT分解重构的D3波形造成了频率混叠，不能完整地复现原始信号在8～16 Hz频段的信息。因此，相比于WPT，DTCWT利用实部树和虚部树的隔点采样和信息互补有效抑制了频率的混叠现象。

为了确定用于脑电信号特征提取的时间段，本文对140个训练集样本中的70个想象左手运动样本和70个想象右手样本分别经DTCWT提取8～16 Hz的脑电信号分量并求取能量均值，想象左右手运动的C3和C4导联能量均值的对比图如图5所示。

图5 想象左右手运动C3和C4导联能量均值对比

从图中可以看出，在想象左手运动时，C3侧脑部区域的活动能量增加，而C3对侧的C4侧脑部区域的活动能量则相应减少。在想象右手运动时，则相反。这一现象有效的证实了ERD/ERS现象。

3.3 特征提取与分类

根据ERD/ERS现象，本文将脑电信号经DTCWT提取8～16 Hz脑电信号分量的能量均值作为C3和C4导联的特征，以此构成二维的特征向量作为样本的输入训练集。相同地，测试集的原始C3和C4导联脑电信号数据也同样经过训练集数据的预处理，构成的二维特征向量作为本文GBDT的输入测试集。其中，能量均值计算表达式为

(13)

观察想象左右手运动的C3和C4导联的能量均值对比图可以发现，C3和C4导联在4～6 s的时间段具有较高的类间区分度，一般地，在3～9 s想象左右手运动的时间段内，4～6 s时间段往往具有较好的类间区分度，所以本文选取了4个时间段4～6 s，4.5～5.5 s，4～5 s和5～6 s作为特征提取的时间段。其中，由于原始脑电信号的采样频率为128 Hz，所以4～6 s时间段包含原始脑电信号D3重构分量的256个采样点，其余3个时间段4.5～5.5 s、4～5 s、5～6 s包含原始脑电信号D3重构分量的128个采样点，分别对这些采样点求取能量均值。最后，通过GBDT对4个时间段内脑电信号进行分类，分类结果如表1所示。

表1 GBDT在各时间段分类准确度

通过比较这4个时间分段的测试集分类准确度可知，在4～6 s时间段训练的GBDT分类器具有较好的分类效果，测试集的分类准确度最高为82.14%，此时训练集的分类准确度为87.86%，测试分类结果达到了与BCI竞赛优胜组第3名成绩相近的水平(第1名89.29%，第2名84.29%，第3名82.86%)。

4 结束语

针对EEG信号的非平稳时变和个体差异性大等特点，本文提出了基于双树复小波变换结合GBDT的脑电信号模式识别的方法。首先经过实验证明了DTCWT相比于WPT在频率域上具有更好的抗混叠能力；接着通过DTCWT提取左右手想象运动的能量均值，验证了ERD/ERS现象；最后通过比较4个时间段的分类准确度，验证了双树复小波变换结合GBDT在想象运动脑电信号识别应用中的有效性，为GBDT分类算法在脑电信号分类的应用中提供了一定的参考价值。