王金甲 胡 备

（燕山大学信息科学与工程学院，秦皇岛 0 66004）

引言

脑机接口（brain-computer interface，BCI）主要目的就是希望建立人脑和计算机（或仪器设备）的直接通信，进而实现人脑对外部设备的控制［1］。从复杂脑电EEG信号中提取可以有效分类的特征就成为了脑机接口中一个严峻挑战。

对于脑机接口左右手想象运动的脑电EEG信号特征选择和提取，很多文献中都有研究报道。研究者采取很多脑电特征对想象运动进行分类，如共空间模式（common spatial pattern，CSP）、核共空间模式（kernel CSP，KCSP）、核线性判别分析（kernel discriminant analysis，KDA）、自回归模型参数等［2－5］。CSP方法针对脑机接口中想象运动的分类效果已被广泛认可［2］，但其缺点在于对两类样本严格的线性假设，针对这一问题文献［3］提出了KCSP，引入核函数很好解决了线性要求限定，新问题是计算复杂性；文献［6］采用聚类方法解决计算复杂性问题，这在一定程度上丢失了部分特征。对此本研究提出了一种基于核方法和广义奇异值分解［7］（generalized singular value decomposition，GSVD）的广义核线性判别分析方法（GKLDA）。首先通过一个非线性映射，将训练数据映射到一个高维特征空间，解决了CSP中对原始数据线性的要求，使得不同类别EEG信号在该空间线性可分（或近似线性可分），采用线性判别分析，寻找一个最佳投影方向，即为最优空域滤波器。在解该滤波器时，为了最大限度的释放不同类间的非线性模式，同时解决小样本采样问题，采用GSVD方法解除了奇异性限定。最后将数据通过该最优空域滤波器，提取最佳非线性特征进行分类。分类器采用fisher线性判别分析分类器。文献［8］提出了一种广义判别分析（generalized discriminant analysis，GDA），本研究提出将该方法应用在脑电分类中。在实验部分，采用了线性、多项式及高斯等三种不同的核函数，将GKLDA与KCSP、KDA、GDA及普通 LDA分别用于相同的脑机接口公开数据进行特征提取对比分类结果，结果表明 GKLDA能取得很好的效果，值得推广。

1 CSP方法及其核推广KCSP

CSP方法所基于的事实是，假定大脑中每一类运动都是由一种或几种源所触发的。假设有两类运动记录下来的EEG信号为X1和X2，对其各自空间协方差矩阵R1，R2取均值，得到平均协方差矩阵和。对复合空间协方差矩阵 R =+做特征值分解，U和λ分别为其特征向量和特征值。定义白化矩阵 P，P=λ－1/2UT。分别对 PPT和PPT做特征值分解，可以证明两者有共同的特征向量，且特征值对角阵之和为单位阵。两者最大m个特征值所对应的特征向量记为V1m和V2m，即可得空域滤波器 F，F=［V1mV2m］P。

KCSP方法中，首先通过一个非线性映射将EEG数据映射到高维空间，即X→Φ。类比于CSP，在这个新空间中有R1=Φ1W1Φ1，从而R=R1+R2=ΦtWΦtT。W是所有样本数n的倒数组成的对角阵，其目的在于求协方差矩阵的平均值。为了在高维空间完成CSP方法，首先要解决R的特征值分解问题。但是直接分解是不可能的，因为通常情况下R维数特别高。文献［6］中给出一个思路，即通过构造特征值分解因式，用已知或可求变量代替未知或不易求变量，通过构造核矩阵的形式代替高维空间的运算。最终滤波器表达式为：

式中，V是KW特征向量，K=ΦtTΦt是核矩阵；P=Λ－1Ψ－1/2VTWΦtT为白化矩阵，Λ是KW特征值；Ψ是一个对角矩阵，它的第j列元素为vjTWvj，vj是V的第j列。细节内容详见文献［3，6］。

2 LDA方法及其核推广KDA

以两类数据为例。设有一集合X，包含n个d维样本x1，x2，…xn，其中n1个属于w1类，n2个属于w2类，n1+n2=n，LDA目标是找到线性投影方向，使得投影后类内散度最小，类间散度最大。Sw为总类内离散度矩阵，Sb为类间离散度矩阵，最佳投影方向w通过下式目标函数最大化：

式（2）为广义Rayleigh商形式，可以用Lagrange乘子法求解，转化为方程Sbw*=λSww*求解问题。解该方程即为式（2）的极值。细节内容见文献［9］。

3 基于核方法和GSVD的广义核线性判别分析

矩阵与其转置的乘积一样，离散度矩阵也分解成其平方根因子乘积的形式。

此时，在映射空间中，目标函数可写为如下Rayleigh商形式

式（6）可转化为求解下式广义特征值问题

据再生核理论，可知式（7）的解空间可限定在由{φ（x）}张成的空间中。令

式中 Y =［φ（x1）…φ（xn）］，α=［α1…αn］T。将φ分别代入式（7），左乘 YT，得

此时，GTY，HTY可以通过下式计算：

本研究提出一个广义核线性判别分析方法。首先，对于EEG数据X，大小为M×T×n，M为通道数，T为每次试验的采样点数，n为试验次数。两类试验次数分别为n1、n2，即n1+n2=n，n为EEG信号所有试验的次数。将数据重新排列为二维数据，每一列为一次实验的所有采样点，即（M×T）×n。假定通过一个非线性的映射函数将输入空间的数据映射到一个希尔伯特空间φ：R→F即X→φ（X），据式（3），在映射后的空间中投影方向向量目标函数变为

散度矩阵与协方差矩阵相似，都是半正定的。正像在许多信号处理问题中协方差矩阵分解为一个

式（9）等价于：

为了释放KG，KH的非奇异限定，对，)采取广义奇异值分解 G SVD，即UTα=［Γ0］，VTα = ［Γ0］，由此可得，

式（11）和式（12）中的α即为式（5）的解，也是所求最优空域滤波器。

对于新输入EEG信号Z，可获得该数据特征，

式中，j=1，2，…m，m为α的列数，一般情况下等于n，但适当的选取m有利于正确率的提高，笔者对m的选择做了进一步实验。式（13）中滤波后的EEG信号其能量可作为特征向量用于分类。

4 实验方法

实验分为两个部分，实验一采用BCI竞赛一、二和三的部分数据分别进行了实验，实验二采用自己的数据做了实验。分别采用了GKLDA、KCSP、KDA、GDA等4种方法进行对比；并且给出相应数据已知的最优正确率作为对比；此外采用了线性核函数，多项式核函数和高斯核函数等作了对比。其中多项式核函数k（x，y）=（xTy）d，d=1.5，d的选择是经过大量实验后的最优结果。高斯核函数k（x，y）=exp（－x－y2/2σ2），σ由十倍交叉验证获得。实验中的分类投影图为GKLDA方法提取特征后的一维投影，其中横轴为样本数，纵轴为提取的特征降到一维后的值。

实验步骤如下：首先是数据预处理，即将三维数据（通道×采样点×试验次数）变为两维（［通道×采样点］×试验次数），符合GKLDA的处理格式；其次对数据实施前文所述的GKLDA以及需要比较的几种方法，提取各类别数据的特征；最后将各类的特征送入分类器，得出该数据集的正确率。

5 实验结果

5.1 实验一

1.数据采用BCI竞赛一的数据集，本数据集来自一名受试者，其任务是根据提示用手指敲击指定键盘，共有413个训练样本，其中两类样本数分别为194和219，100个测试样本。27导信号，采样率为100Hz。GKLDA降维后的投影见图1，其中横轴为试验次数n，纵轴为提取特征后将特征维数降为一维时的投影幅值，不同类别分别用星号和圆点表示，而最后分类并没有降维。

4种方法下的训练和测试结果见表1。高斯核中KCSP的参数2σ2=1×1010，KDA的参数2σ2=4×106，GDA的参数2σ2=1×106，GKLDA的参数2σ2=1×106。GKLDA正确率较高。本数据集用普通LDA特征和CSP提取特征的分类错误率结果分别为为11.0％和23.0％。竞赛第一名的错误率为4.0％。

图1 数据一的投影结果。（a）Linear核训练数据；（b）Linear核测试数据；（c）Polynomial核训练数据；（d）Polynomial核测试数据；（e）Gaussian核训练数据；（f）Gaussian核测试数据Fig.1 The projection of data I.（a）Training data with linear kernel；（b）Testing data with linear kernel；（c）Training data with polynomial kernel；（d）Testing data with polynomial kernel；（e）Training data with Gaussian kernel；（f）Testing data with Gaussian kernel

表1 数据一在4种方法下的分类错误率（％）Tab.1 The error rate（％）on datasetⅡ with four different methods

2.数据采用BCI竞赛二的数据集Ⅳ，该数据集为左右手运动，受试者坐在椅子上，手臂放松放在桌子上，其任务是敲击指定按键。共有样本416个，其中316个为训练集，左右手分别为159和157，其他100个作为测试集。电极数为28，采样率100Hz，时间长度为0.5 s。图2所示GKLDA降维后的投影结果，其中横轴为试验次数n，纵轴为提取特征后将特征维数降为一维时的投影幅值，不同类别分别用星号和圆点表示。而最后分类并没有降维。

4种方法下的训练和测试结果见表2。高斯核中 K CSP的参数 2 σ2=6.4×107，KDA的参数 2 σ2=2.5×107，GDA的参数 2 σ2=6.4×107，GKLDA的参数2σ2=6.4×107。GKLDA正确率较高。本数据集用普通LDA特征提取直接分类，结果错误率为28.0％；文献［11］中采用 C SSD和 F DA结合的方法，错误率为16.0％。

3.数据采用 BCI竞赛三的数据集ⅢB中的S4b，该数据集为在线反馈的想象左右手运动，反馈时间为第4 s到第7 s。共有样本1 080个，其中540个为训练集，左右手各为270，其他540个作为测试集。电极数为2，采样率125Hz，时间长度为3 s。图3所示为GKLDA降维后的投影结果，其中横轴为试验次数n，纵轴为提取特征后将特征维数降为一维时的投影幅值，其中（c）和（d）为多项式核提取特征后的二维投影，不同类别分别用星号和圆点表示。而最后分类并没有降维。

表2 数据二在4种方法下的分类错误率（％）Tab.2 The error rate（％）on dataset II with four different methods

4种方法下的训练和测试结果见表3，高斯核中KCSP的参数2σ2=2 500，KDA的参数2σ2=2 500，GDA的参数2σ2=10 000，GKLDA的参数2σ2=10 000。GKLDA正确率较高。实验中试着改变（13）式中α的列数，在线性核下有明显效果，减少α的列数能明显提高正确率，当α全选时错误率为40.0％，而当α为210列时，错误率降为21.3％，所以在此α的列数取为210；多项式核情况下没有明显效果，所以在此α全选；高斯核函数α列数取为450，比全部的540正确率提高近15.0％，所以在此α保留450列。本数据集用普通LDA特征提取的分类错误率为39.1％；文献［12］中采用Chiu＇s FIS（模糊推理系统）方法，错误率为25.3％；文献［13］中采取自回归的方法，错误率为17.4％。

图3 数据三的投影结果。（a）Linear核训练数据；（b）Linear核测试数据；（c）Polynomial核训练数据；（d）Polynomial核测试数据；（e）Gaussian核训练数据；（f）Gaussian核测试数据Fig.3 The Projection of data III.（a）Training data with linear kernel；（b）Testing data with linear kernel；（c）Training data with polynomial kernel；（d）Testing data with polynomial kernel；（e）Training data with Gaussian kernel；（f）Testing data with Gaussian kernel

5.2 实验二

数据格式说明：受试者为一健康男性，采用gtec信号放大器，C3和C4导联，时间长为3 s，采样率为125Hz，任务为根据提示想象左右手运动。共采集样本300，其中训练样本200，左右手各100，测试样本100。4种方法下的训练和测试结果见表4，多项式核中设置d=1.5，分类效果较为明显；高斯核中KCSP的参数2σ2=10 000，KDA的参数2σ2=2 500，GDA的参数2σ2=10 000，GKLDA的参数2σ2=10 000。GKLDA正确率较高。本数据集用普通LDA特征的分类错误率为13.0％。

表3 数据三在4种方法下的分类错误率（％）Tab.3 The error rate（％）on dataset Ⅲ with four different methods

表4 数据四在3种方法下的分类错误率（％）Tab.4 The error rate（％）on dataset Ⅳ with four different methods

6 讨论和结论

上述3组公开数据集中，所提出的GKLDA方法均能在适当选择核函数和其适当参数后取得较好结果。相较于其他核方法（KCSP、KDA、GDA）和普通线性判别分析，竞赛二中本方法略优，而竞赛一竞赛三中本方法正确率明显超过其他核方法。每组实验结果中，本方法并不是最优的，最优的方法往往只适合一组数据，不具有通用性。而GKLDA往往能在适当调节参数后适应各类数据（也即个体差异），这是其他无参数方法所不具备的。所以从通用性角度来讲，GKLDA更为突出。

核函数的参数选择方面，多项式核d=1.5；高斯核函数的σ参数选择过程如下，在交叉验证时首先以10的幂次为阶数，大致找到一个范围。不失一般性以数据三为例，首先将训练数据按9∶1的比例分为训练数据和测试数据，然后参数σ的选择范围在10－4～105，据此，可以判断该数据的参数范围大致应该在101～104之间，且最有可能在102附近，故下一步可进行交叉验证选取最优参数。α的列数选择方面，根据多组实验结果，多数情况下是选取全部列数为优，数据三中α的列数选取也是经过大量测试选取的，因为训练数据中某些样本的识别率较低（这从其降维后的投影中可以看出来），故在选取滤波器时稍有差异。从3种核函数的实验结果可以看出，多项式核则较为不稳定，结果稍差；线性核没有可调参数，故而有较强的普遍性，实验结果较好；而高斯核有参数可调（σ），对不同的样本（即不同的人）有较强的针对性，故实验结果最好。

脑电信号分类中，由于信号与信号源之间的复杂性，其大多数有用特征跟源之间都是非线性的关系，所以这就限制了CSP和线性判别分析的应用。本研究提出的广义核线性判别分析方法，核心在于结合了核方法和广义奇异值分解，核方法的加入解决了一般方法对线性要求的局限性，在核函数的选择上，由于高斯核函数的参数可变性，针对个体差异性找出相应参数值，取得了较好的结果，又由于GSVD的应用而减弱了对非奇异的严格要求，从而解决了小样本采样问题，使分类结果更加有效。实验表明该方法在脑电信号处理中值得推广。未来进一步的工作可以放在将该方法推广至多类别分类中。

［1］Wolpaw.JR，Birbaumer N，McFarlandDJ，etal.Braincomputer interfaces for communication and control［J］.Clinical Neurophysiology，2002，113（6）：767－791.

［2］Wang Yijun，Gao Shangkai，Gao Xiaorong.Common spatial pattern method for channel selection in motor imagery based brain-computer interface［C］//IEEE EMBS 27th Annual Conference 2005.Shanghai：IEEE，2005：5392－5395.

［3］Sun Shiliang，Zhang Changshui.An optimal kernel feature extractor and its application to EEG signal classification［J］.Neurocomputing 2005，69（13－15）：1743－1748.

［4］高湘萍，许丹，吴小培.基于核Fisher判别分析的意识任务识别新方法［J］.计算机技术与发展，2006，16（9）：82－84.

［5］王金甲，周丽娜，赵玉超.基于MEG的脑机接口特征提取方法研究［J］.仪器仪表学报，2010，31（7）：1460－1465.

［6］Nasihatkon B，Boostani R，Jahromi M.An efficient hybrid linear and kernel CSP approach for EEG feature extraction［J］.Neurocomputing.2009，73（1－3）：432－437.

［7］赵海龙，穆志纯，张霞.基于LDA/GSVD和支持向量机的人耳识别［J］.上海理工大学学报，2009，31（6）：601－608.

［8］Baudat G，Anouar F.Generalized discriminant analysis using a kernel approach［J］.Neural Comput，2000，12（10）：2385－2404.

［9］边肇祺，张学工.模式识别［M］.北京：清华大学出版社，2000：84－90.

［10］Wang Yanxia，Ruan Qiuqi.Kernel fisher discriminant analysis forpalmprintrecognition［C］// IEEE 18th International Conference on Pattern Recognition.Washington DC：IEEE，2006：457－460.

［11］Wang Yijun；Zhang Zhiguang，Li Yong，et al.BCI competition 2003—data set IV：an algorithm based on CSSD and FDA for classifying single-trial EEG［J］.IEEE transactions on biomedical engineering 2004，51（6）：1081－1086.

［12］Lotte F，Lécuyer A，Lamarche F，et al.Studying the use of fuzzy inference systems for motor imagery classification［J］.IEEE Transactions on Neural System and Rehabilitation Engineering，2007，15（2）：322－324.

［13］Brodu N，Lotte F，Lécuyer A.Comparative study of bandpower extraction techniques for motor imagery［C］，IEEE Symposium Serieson Computational Intelligence，Cognitive Algorithms，Mind，and Brain，2011：95 －100.

［14］刘广权，黄淦，朱向阳.共空域模式方法在多类别分类中的应用［J］.中国生物医学工程学报，2009，28（6）：935－938.