在左右手运动想象的脑电(EEG)分析方法中,目前大多针对多通道采集的EEG数据,难以应用到单通道脑机接口(BCI)中。本文采用一种改进的独立成分分析(ICA)方法,实现了对EEG数据有效的预处理。首先,对数据进行线性漂移校正,去除数据漂移;然后采用延时窗口数据增加虚拟通道数,利用ICA除去EEG数据中的伪迹,即眼电和心电;最后利用希尔伯特-黄变换(HHT)后的瞬时幅值,求平均瞬时能量特征并分类。实验证明,该方法测试性完成了EEG数据的预处理工作,提高了单通道EEG信号的分类率,可为单通道的便携式BCI研究打下基础。
引用本文: 李松, 伏云发, 杨秋红, 刘传伟, 孙会文. 基于左右手运动想象单通道脑电信号的预处理研究. 生物医学工程学杂志, 2016, 33(5): 862-866. doi: 10.7507/1001-5515.20160139 复制
引言
脑机接口(brain-computer interface, BCI)是一种新型的人机交互技术,是指不依赖肌肉的参与,把大脑活动的信号转化为电信号,通过对电信号的解码分析,实现与外部设备的通信与控制[1-3]。近年来BCI实现了由离线到在线的飞速发展,尤其是在第一届和第二届中国BCI竞赛中,基于运动想象的BCI在控制机器人比赛中取得了优异成绩,但采用多通道采集的脑电(electroencephalograpm, EEG)信号并不适宜应用到便携式设备中。
EEG信号是电极记录下来的脑细胞群活动,表现微弱,只有微伏数量级,极其容易受到数据漂移、眼电(electrooculogram, EOG)和心电(electrocardiograpm, ECG)的干扰[4-6]。因此,在基于单通道EEG信号采集中,有效去除伪迹又不使原EEG信号失真是一项重要的研究内容。目前,主要的伪迹去除方法有独立成分分析法(independent component analysis, ICA)和平均伪迹回归分析方法[7-9]。前者主要是一种盲源分离方法,假设采集的数据为伪迹信号和EEG信号混合而成,则该方法可通过机器学习的方法分离并去除伪迹信号,再重构EEG信号。但这种方法有一定的限制条件:①混合信息必须是线性的;②信息源是瞬时相互独立的;③通道数大于源数目。实际情况中,一般前两个条件都能满足,但随着BCI设备的便携式化发展,将难以满足第③个条件。而平均伪迹回归分析方法虽应用较早且简单,但需要增加通道以实时测出EOG和ECG成分,也难于应用到便携式EEG设备中。
通过以上比较说明,当前具有代表性的两种方法均难以适用于少通道或单通道的BCI设备,因此针对ICA通道数目限制问题,本文提出采用一种改进的ICA方法,即延时窗口数据以增加虚拟通道数目来满足ICA的条件。本文通过对改进ICA方法进行EEG数据预处理研究,一定程度上验证了该方法的有效性,为研究者提供了一些新的研究思路,可在便携式BCI应用领域研究中深入和推广。
1 资料和方法
1.1 受试者信息
本文共召集5名健康受试者(S1~S5)参与基于左右手运动想象的EEG数据采集试验,受试者信息如下:男,23~30岁,本科以上学历,工科背景,所有受试者皆身体健康,均没有感觉运动疾病和心理病史,且无EEG数据采集的试验经验。所有受试者在试验前均签署了知情同意书。
1.2 试验范式和试验过程
本次试验设计了左右手想象运动,试验时要求受试者根据屏幕随机提示的左右箭头进行想象左手或右手运动,并实时采集EEG数据,每个试验环节受试者进行10个测试。在正式采集数据前,要求去除影响试验的一些因素,并且熟悉试验环境,理解试验步骤。试验时,试验室环境保持安静,受试者身体放直,端坐在舒适的椅子上,双手平放在桌子上,面对呈现刺激的电脑屏幕(距离受试者0.6~1 m),按屏幕提示完成试验。
每个受试者一次试验完成10个测试,每个测试分为练习和正式。正式的一个采集时序如图 1所示,包括:开始时,1.5 s的固定十字(含0.5 s的蜂鸣声),为提示准备试验;第1.5秒开始有1.5 s的左右箭头任务提示,第3秒开始有5 s的*表示任务执行时间和第8秒开始随机6~8 s的空屏表示休息时间。该刺激范式由常用于心理学、认知神经科学和生物医学工程领域的刺激软件E-prime完成。

1.3 EEG数据采集
本次试验的EEG采集装置采用具有16导联的EEG放大器(Mipower-UC, EEG Collection_V2, 清华大学神经工程实验室)进行试验,该放大器常规参数如下:采样率为1 000 Hz;24位模数(analog-digital,AD)转换器;信号频带:0~250 Hz,输入信号范围:±200 mV;同步事件输入:8位输入。脑电帽采用国际标准的10~20系统的16导电极帽(Ag-AgCl粉末电极, 武汉格林泰克科技有限公司)。此次试验只用覆盖运动功能区的C3和C4电极,Fpz电极接地,M1电极接左侧乳突,试验要求所有电极的阻抗均小于5 kΩ。
1.4 数据处理
1.4.1 ICA基本原理
ICA算法主要是从多个通道数据中分离出相互独立的源,即盲源分离(blind source separation, BSS)[10-13]。设x(t)为N个通道采集的数据,该数据中包含了伪迹成分,s(t)为真实的源信号,x(t)=As(t),即观测信号x(t)是经过s(t)线性混合后而成,其中A为混合矩阵。即该方法在混合矩阵A和s(t)未知情况下,利用源瞬时相互独立,找到一个线性变换分离矩阵W,使得输出信号y(t)=Wx(t)尽可能地逼近真实信号s(t),可以表示为如式(1)、(2)所示:
$ x\left( t \right) = \boldsymbol{A}s\left( t \right) $ |
$ x\left( t \right) = \boldsymbol{W}x\left( t \right) $ |
式中x(t)=[x1(t), x2(t), …, xn(t)]T为N个通道的观测数据,s(t)=[s1(t), s2(t), …, sn(t)]T为N个通道相互独立的源信号。由以上可知找到分离变换矩阵W是关键,本文采用基于负熵的FastICA算法寻找最佳W。
1.4.2 通道变换
用h(t)表示单次单通道提取的左右手运动想象的EEG数据,其包含了M s的真实EEG信号和伪迹信号,本文采用延时N s的窗口数据作为虚拟通道数,提取L s的数据作为处理对象,由单通道变为i个通道。则通道模型可以表示为如式(3)~(5)所示[14-15]:
$ {V_i}\left( t \right) = {h_{0 \sim L}}\left( t \right)\;\;\;\;\;\;i = 1 $ |
$ {V_i}\left( t \right) = {h_{\left( {i - 1} \right) \cdot N \sim \left( {i - 1} \right) \cdot N + L}}\left( t \right)\;\;\;\;\;\;i \ge 1 $ |
$ h\left( t \right) = \left[ {{v_1}\left( t \right),{v_2}\left( t \right), \cdots ,{v_i}\left( t \right)} \right] $ |
则ICA算法可以表示如式(6)所示:
$ \left[ \begin{array}{l} {v_1}\left( t \right)\\ {v_2}\left( t \right)\\ \vdots \\ {v_i}\left( t \right) \end{array} \right] = A\left[ \begin{array}{l} {s_1}\left( t \right)\\ {s_2}\left( t \right)\\ \vdots \\ {s_i}\left( t \right) \end{array} \right] $ |
式中vi(t)表示虚拟的通道,si(t)表示相互独立的源信号。在本次实验中M=5, N=1, L=3, 即包含5 s的数据,延时1 s,取3 s的数据作为处理对象时效果较好。
1.4.3 功率比计算
本文首先采用在时域和频域同时具有很高分辨率的希尔伯特-黄变换(Hilbert-Huang transform, HHT)提取瞬时幅值,并求瞬时能量,然后选用功率比来观察运动想象所特有的事件相关去同步(event-related desynchronization, ERD)现象。步骤如下[16]:
(1)根据本征模函数(intrinsic mode function, IMF)进行经验模态分解(empirical mode decomposition, EMD),直到分解满足IMF条件定义,即局部极值点和过零点数目必须相等或者相差一个,局部最大值包络和最小值包络平均为零。分解后的信号可看成由式(7)组成
$ x\left( t \right) = \sum\limits_{i = 1}^n {{h_i}\left( t \right) + {r_n}\left( t \right)} $ |
式(7)中x(t)为原始信号,hi(t)表示第i个IMF分量,rn(t)为剩余分量。
(2)选取分类效果较好的前三阶,按式(8)进行Hilbert谱分析,并按其解析信号式(9)、(10)求其瞬时幅值:
$ {y_i}\left( t \right) = \frac{1}{{\rm{\pi }}}\int\limits_{ - \infty }^{ + \infty } {\frac{{{h_i}\left( \tau \right)}}{{t - \tau }}{{\rm{d}}_\tau }} $ |
$ {x_i}\left( t \right) = {h_i}\left( t \right) + j{y_i}\left( t \right) $ |
$ {A_i}\left( t \right) = \sqrt {{y_i}{{\left( t \right)}^2} + {h_i}{{\left( t \right)}^2}} $ |
(3)根据式(11)计算平均瞬时能量值
$ AE{C_N} = \frac{1}{N}\sum\limits_n^N {C_i^2} $ |
式中N表示采样点数目,Ci2表示第i个采样点瞬时能量值。
(4)功率比计算公式如下[17]:
$ A{E_{\max }} = \max \left( {{c_3}{e_{\max }},{c_4}{e_{\max }}} \right) $ |
$ A{E_{\min }} = \max \left( {{c_3}{e_{\min }},{c_4}{e_{\min }}} \right) $ |
$ {I_n} = \frac{{e_i^n - A{E_{\min }}}}{{A{E_{\min }} - A{E_{\min }}}} $ |
式(12)中c3emax, c4emin分别代表一个测试时,c3和c4通道瞬时能量的最大值和最小值;AEmax是求c3、c4通道瞬时能量的最大值;同理,式(13)中可求得AEmin。ein为c3或者c4通道的瞬时能量值。式(14)即为功率比求解。
1.4.4 EEG信号模式分类
为了进一步验证改进的ICA预处理方法对EEG信号分类效果的影响,本文采集了5名受试者,共计150个测试(每名受试者30个测试)的离线EEG数据,首先分别对数据进行了改进的ICA预处理和未改进的ICA预处理,然后提取瞬时能量特征,最后用支持向量机(support vector machine, SVM)和基于Fisher核的线性判别分析(linear discriminant analysis, LDA)进行分类。
2 结果与分析
此次试验数据首先经过降采样处理,由1 000 Hz降为250 Hz,采集的5名受试者的部分原始数据如图 2所示,表示c3通道2个测试的EEG数据,从图中可以看出数据受到漂移、EOG、ECG的严重污染。

针对原始数据存在的严重污染,本文先通过阈值法漂移校正,去除线性漂移,校正后如图 3所示,表示c3通道4个测试的EEG数据,可以看出经过校正后,数据去除了漂移现象,但还是存在ECG和EOG的干扰。

然后把原始数据中一个测试的数据进行窗口延时1 s,取3 s的数据作为虚拟通道数据,使用ICA去除EOG和ECG的干扰。把c3通道进行通道变换后虚拟成3个通道,同理c4通道进行变换后也虚拟成3个通道,即可看成6导联的数据进行ICA处理, 两个通道处理后的效果分别如图 4所示,第一排表示c3通道一个测试的数据变换成3个通道经ICA处理后的效果,第二排表示c4通道一个测试的数据变换成3个通道经ICA处理后的效果。

如图 4所示,为某受试者的某个测试数据进行通道变换后的伪迹去除效果图,最后选择3个虚拟通道中效果较好的结果进行功率比计算。计算时先对数据进行8~30 Hz的带通滤波,对每N个采样点进行平均瞬时能量的计算,然后求取功率比。结果如图 5所示,为受试者在进行左右手运动想象时的某个测试的能量图,验证了运动想象ERD现象。

5名受试者基于SVM和LDA的分类准确率如图 6所示,表明经过改进的ICA预处理与未改进的ICA相对比,提高了EEG信号的分类率。但总体分类准确率可以通过融合多种特征,进一步提高。

3 结论
基于运动想象的单通道或少通道的便携式EEG信号采集设备正在成为BCI技术的一个发展趋势。近年来,ICA和共空间模式(common spatial pattern, CSP)是分析基于运动想象EEG预处理和分类的主流方法,但这两种方法都是基于多通道模式,不适合少通道EEG数据的分析。本文基于改进的ICA方法,针对EEG信号易受到漂移、EOG信号、ECG信号的干扰,采用了漂移校正,用改进ICA方法去除EOG信号和ECG信号,通过带通滤波后提取HHT特征,求取功率比和分类率,试验数据得到了较好的预处理;同时比较功率比和分类率,验证了左右手运动想象所产生的ERD现象,提高了EEG分类率。本研究可望为促进便携式BCI应用领域的研究提供一些思路。
引言
脑机接口(brain-computer interface, BCI)是一种新型的人机交互技术,是指不依赖肌肉的参与,把大脑活动的信号转化为电信号,通过对电信号的解码分析,实现与外部设备的通信与控制[1-3]。近年来BCI实现了由离线到在线的飞速发展,尤其是在第一届和第二届中国BCI竞赛中,基于运动想象的BCI在控制机器人比赛中取得了优异成绩,但采用多通道采集的脑电(electroencephalograpm, EEG)信号并不适宜应用到便携式设备中。
EEG信号是电极记录下来的脑细胞群活动,表现微弱,只有微伏数量级,极其容易受到数据漂移、眼电(electrooculogram, EOG)和心电(electrocardiograpm, ECG)的干扰[4-6]。因此,在基于单通道EEG信号采集中,有效去除伪迹又不使原EEG信号失真是一项重要的研究内容。目前,主要的伪迹去除方法有独立成分分析法(independent component analysis, ICA)和平均伪迹回归分析方法[7-9]。前者主要是一种盲源分离方法,假设采集的数据为伪迹信号和EEG信号混合而成,则该方法可通过机器学习的方法分离并去除伪迹信号,再重构EEG信号。但这种方法有一定的限制条件:①混合信息必须是线性的;②信息源是瞬时相互独立的;③通道数大于源数目。实际情况中,一般前两个条件都能满足,但随着BCI设备的便携式化发展,将难以满足第③个条件。而平均伪迹回归分析方法虽应用较早且简单,但需要增加通道以实时测出EOG和ECG成分,也难于应用到便携式EEG设备中。
通过以上比较说明,当前具有代表性的两种方法均难以适用于少通道或单通道的BCI设备,因此针对ICA通道数目限制问题,本文提出采用一种改进的ICA方法,即延时窗口数据以增加虚拟通道数目来满足ICA的条件。本文通过对改进ICA方法进行EEG数据预处理研究,一定程度上验证了该方法的有效性,为研究者提供了一些新的研究思路,可在便携式BCI应用领域研究中深入和推广。
1 资料和方法
1.1 受试者信息
本文共召集5名健康受试者(S1~S5)参与基于左右手运动想象的EEG数据采集试验,受试者信息如下:男,23~30岁,本科以上学历,工科背景,所有受试者皆身体健康,均没有感觉运动疾病和心理病史,且无EEG数据采集的试验经验。所有受试者在试验前均签署了知情同意书。
1.2 试验范式和试验过程
本次试验设计了左右手想象运动,试验时要求受试者根据屏幕随机提示的左右箭头进行想象左手或右手运动,并实时采集EEG数据,每个试验环节受试者进行10个测试。在正式采集数据前,要求去除影响试验的一些因素,并且熟悉试验环境,理解试验步骤。试验时,试验室环境保持安静,受试者身体放直,端坐在舒适的椅子上,双手平放在桌子上,面对呈现刺激的电脑屏幕(距离受试者0.6~1 m),按屏幕提示完成试验。
每个受试者一次试验完成10个测试,每个测试分为练习和正式。正式的一个采集时序如图 1所示,包括:开始时,1.5 s的固定十字(含0.5 s的蜂鸣声),为提示准备试验;第1.5秒开始有1.5 s的左右箭头任务提示,第3秒开始有5 s的*表示任务执行时间和第8秒开始随机6~8 s的空屏表示休息时间。该刺激范式由常用于心理学、认知神经科学和生物医学工程领域的刺激软件E-prime完成。

1.3 EEG数据采集
本次试验的EEG采集装置采用具有16导联的EEG放大器(Mipower-UC, EEG Collection_V2, 清华大学神经工程实验室)进行试验,该放大器常规参数如下:采样率为1 000 Hz;24位模数(analog-digital,AD)转换器;信号频带:0~250 Hz,输入信号范围:±200 mV;同步事件输入:8位输入。脑电帽采用国际标准的10~20系统的16导电极帽(Ag-AgCl粉末电极, 武汉格林泰克科技有限公司)。此次试验只用覆盖运动功能区的C3和C4电极,Fpz电极接地,M1电极接左侧乳突,试验要求所有电极的阻抗均小于5 kΩ。
1.4 数据处理
1.4.1 ICA基本原理
ICA算法主要是从多个通道数据中分离出相互独立的源,即盲源分离(blind source separation, BSS)[10-13]。设x(t)为N个通道采集的数据,该数据中包含了伪迹成分,s(t)为真实的源信号,x(t)=As(t),即观测信号x(t)是经过s(t)线性混合后而成,其中A为混合矩阵。即该方法在混合矩阵A和s(t)未知情况下,利用源瞬时相互独立,找到一个线性变换分离矩阵W,使得输出信号y(t)=Wx(t)尽可能地逼近真实信号s(t),可以表示为如式(1)、(2)所示:
$ x\left( t \right) = \boldsymbol{A}s\left( t \right) $ |
$ x\left( t \right) = \boldsymbol{W}x\left( t \right) $ |
式中x(t)=[x1(t), x2(t), …, xn(t)]T为N个通道的观测数据,s(t)=[s1(t), s2(t), …, sn(t)]T为N个通道相互独立的源信号。由以上可知找到分离变换矩阵W是关键,本文采用基于负熵的FastICA算法寻找最佳W。
1.4.2 通道变换
用h(t)表示单次单通道提取的左右手运动想象的EEG数据,其包含了M s的真实EEG信号和伪迹信号,本文采用延时N s的窗口数据作为虚拟通道数,提取L s的数据作为处理对象,由单通道变为i个通道。则通道模型可以表示为如式(3)~(5)所示[14-15]:
$ {V_i}\left( t \right) = {h_{0 \sim L}}\left( t \right)\;\;\;\;\;\;i = 1 $ |
$ {V_i}\left( t \right) = {h_{\left( {i - 1} \right) \cdot N \sim \left( {i - 1} \right) \cdot N + L}}\left( t \right)\;\;\;\;\;\;i \ge 1 $ |
$ h\left( t \right) = \left[ {{v_1}\left( t \right),{v_2}\left( t \right), \cdots ,{v_i}\left( t \right)} \right] $ |
则ICA算法可以表示如式(6)所示:
$ \left[ \begin{array}{l} {v_1}\left( t \right)\\ {v_2}\left( t \right)\\ \vdots \\ {v_i}\left( t \right) \end{array} \right] = A\left[ \begin{array}{l} {s_1}\left( t \right)\\ {s_2}\left( t \right)\\ \vdots \\ {s_i}\left( t \right) \end{array} \right] $ |
式中vi(t)表示虚拟的通道,si(t)表示相互独立的源信号。在本次实验中M=5, N=1, L=3, 即包含5 s的数据,延时1 s,取3 s的数据作为处理对象时效果较好。
1.4.3 功率比计算
本文首先采用在时域和频域同时具有很高分辨率的希尔伯特-黄变换(Hilbert-Huang transform, HHT)提取瞬时幅值,并求瞬时能量,然后选用功率比来观察运动想象所特有的事件相关去同步(event-related desynchronization, ERD)现象。步骤如下[16]:
(1)根据本征模函数(intrinsic mode function, IMF)进行经验模态分解(empirical mode decomposition, EMD),直到分解满足IMF条件定义,即局部极值点和过零点数目必须相等或者相差一个,局部最大值包络和最小值包络平均为零。分解后的信号可看成由式(7)组成
$ x\left( t \right) = \sum\limits_{i = 1}^n {{h_i}\left( t \right) + {r_n}\left( t \right)} $ |
式(7)中x(t)为原始信号,hi(t)表示第i个IMF分量,rn(t)为剩余分量。
(2)选取分类效果较好的前三阶,按式(8)进行Hilbert谱分析,并按其解析信号式(9)、(10)求其瞬时幅值:
$ {y_i}\left( t \right) = \frac{1}{{\rm{\pi }}}\int\limits_{ - \infty }^{ + \infty } {\frac{{{h_i}\left( \tau \right)}}{{t - \tau }}{{\rm{d}}_\tau }} $ |
$ {x_i}\left( t \right) = {h_i}\left( t \right) + j{y_i}\left( t \right) $ |
$ {A_i}\left( t \right) = \sqrt {{y_i}{{\left( t \right)}^2} + {h_i}{{\left( t \right)}^2}} $ |
(3)根据式(11)计算平均瞬时能量值
$ AE{C_N} = \frac{1}{N}\sum\limits_n^N {C_i^2} $ |
式中N表示采样点数目,Ci2表示第i个采样点瞬时能量值。
(4)功率比计算公式如下[17]:
$ A{E_{\max }} = \max \left( {{c_3}{e_{\max }},{c_4}{e_{\max }}} \right) $ |
$ A{E_{\min }} = \max \left( {{c_3}{e_{\min }},{c_4}{e_{\min }}} \right) $ |
$ {I_n} = \frac{{e_i^n - A{E_{\min }}}}{{A{E_{\min }} - A{E_{\min }}}} $ |
式(12)中c3emax, c4emin分别代表一个测试时,c3和c4通道瞬时能量的最大值和最小值;AEmax是求c3、c4通道瞬时能量的最大值;同理,式(13)中可求得AEmin。ein为c3或者c4通道的瞬时能量值。式(14)即为功率比求解。
1.4.4 EEG信号模式分类
为了进一步验证改进的ICA预处理方法对EEG信号分类效果的影响,本文采集了5名受试者,共计150个测试(每名受试者30个测试)的离线EEG数据,首先分别对数据进行了改进的ICA预处理和未改进的ICA预处理,然后提取瞬时能量特征,最后用支持向量机(support vector machine, SVM)和基于Fisher核的线性判别分析(linear discriminant analysis, LDA)进行分类。
2 结果与分析
此次试验数据首先经过降采样处理,由1 000 Hz降为250 Hz,采集的5名受试者的部分原始数据如图 2所示,表示c3通道2个测试的EEG数据,从图中可以看出数据受到漂移、EOG、ECG的严重污染。

针对原始数据存在的严重污染,本文先通过阈值法漂移校正,去除线性漂移,校正后如图 3所示,表示c3通道4个测试的EEG数据,可以看出经过校正后,数据去除了漂移现象,但还是存在ECG和EOG的干扰。

然后把原始数据中一个测试的数据进行窗口延时1 s,取3 s的数据作为虚拟通道数据,使用ICA去除EOG和ECG的干扰。把c3通道进行通道变换后虚拟成3个通道,同理c4通道进行变换后也虚拟成3个通道,即可看成6导联的数据进行ICA处理, 两个通道处理后的效果分别如图 4所示,第一排表示c3通道一个测试的数据变换成3个通道经ICA处理后的效果,第二排表示c4通道一个测试的数据变换成3个通道经ICA处理后的效果。

如图 4所示,为某受试者的某个测试数据进行通道变换后的伪迹去除效果图,最后选择3个虚拟通道中效果较好的结果进行功率比计算。计算时先对数据进行8~30 Hz的带通滤波,对每N个采样点进行平均瞬时能量的计算,然后求取功率比。结果如图 5所示,为受试者在进行左右手运动想象时的某个测试的能量图,验证了运动想象ERD现象。

5名受试者基于SVM和LDA的分类准确率如图 6所示,表明经过改进的ICA预处理与未改进的ICA相对比,提高了EEG信号的分类率。但总体分类准确率可以通过融合多种特征,进一步提高。

3 结论
基于运动想象的单通道或少通道的便携式EEG信号采集设备正在成为BCI技术的一个发展趋势。近年来,ICA和共空间模式(common spatial pattern, CSP)是分析基于运动想象EEG预处理和分类的主流方法,但这两种方法都是基于多通道模式,不适合少通道EEG数据的分析。本文基于改进的ICA方法,针对EEG信号易受到漂移、EOG信号、ECG信号的干扰,采用了漂移校正,用改进ICA方法去除EOG信号和ECG信号,通过带通滤波后提取HHT特征,求取功率比和分类率,试验数据得到了较好的预处理;同时比较功率比和分类率,验证了左右手运动想象所产生的ERD现象,提高了EEG分类率。本研究可望为促进便携式BCI应用领域的研究提供一些思路。