基于HHT高频组合编码稳态视觉诱发电位脑机接口方法

摘要

基于HHT高频组合编码稳态视觉诱发电位脑机接口方法，先进行硬件连接，再通过n个高频刺激频率经时序排列组合形成的nn个刺激序列表示不同目标，并对刺激序列编码形成编码库，然后采用基于希尔伯特-黄变换(Hilbert-Huang Transform，HHT)的变频脑电信号特征提取方法提取脑电数据特征，并得到脑电数据的Hilbert时-频图，最后使用局部频谱极值目标辨识方法，通过提取的脑电数据特征与编码库比较来实现量化目标识别正确率。本发明具有操作简单、电极数目少、目标数目多、减少受试者疲劳和降低了诱发受试者癫痫可能的优点。

说明书

技术领域

本发明涉及脑-机接口(Brain-Computer Interface，BCI)技术领域，具体涉及基于HHT高频组合编码稳态视觉诱发电位脑机接口方法。

背景技术

脑-机接口是基于脑电信号实现大脑和计算机或其他电子设备直接交流通讯和控制的系统。脑-机接口(BCI)作为人-机接口(Human-computer interface，简称HCI)中的一种，由于不依赖常规的大脑输出通路，为人的大脑开辟了一条全新的与外界进行信息交流和控制的途径，使得利用人脑信号直接控制外部设备的想法成为可能。近年来，脑-机接口(BCI)技术发展非常迅速，已经在生物医学、虚拟现实、游戏娱乐、康复工程以及航天、军事等领域体现出重要的价值。

稳态视觉诱发电(Steady State Visually Evoked Potential，SSVEP)是大脑视觉系统对外部持续周期性视觉刺激的响应，是脑-机接口BCI系统中常用的输入信号，相比P300、事件相关同步、自发脑电等信号，稳态视觉诱发电位(SSVEP)由于具有操作简单、记录电极少、训练时间短、信息传输率高及抗干扰能力强等优点，是BCI系统中具有广泛应用前景和应用价值的重要研究范式。

但目前一般的SSVEP-BCI系统常规采用6～25Hz的低频段，基本都采用一个任务用一个频率表示的简单范式；而由于频率分辨率及屏幕刷新率的限制，在有限的频带范围内使得可呈现的任务目标数有限；另外，在实际应用中，低频段SSVEP的低频闪烁刺激容易使受试者产生视觉疲劳，甚至存在诱发受试者癫痫发作的可能，因此无法保证系统长时间使用的舒适性。

高频SSVEP已经接近人眼的光驱动临界频率，一般频率范围可达30～60Hz或更高，与低频不同的是高频刺激会产生闪烁熔合(flicker fusion)效应，使用者主观上已经感觉不到闪烁，但在脑电中仍可检测到SSVEP的高频响应，这样就大大降低了视觉疲劳和诱发癫痫的可能性，实现了无损；同时，虽然高频SSVEP幅度较低，该频段内背景脑电的幅度也很低，因此高频SSVEP仍然具有一定的信噪比，保障了刺激响应信号的基本辨识准确率。

希尔伯特-黄变换(Hilbert-Huang Transform，HHT)是一种自适应的信号处理方法，非常适合处理非线性、非平稳信号。这种方法包括两个过程：经验模式分解(empirical mode decomposition，EMD)和Hilbert变换(Hilbert Transform，HT)，其中最关键的部分是EMD方法，EMD方法基于信号的局部特征时间尺度，能把复杂信号分解为本征模式函数(Intrinsic Mode Functions，IMF)之和，每一个IMF所包含的频率成分不仅与分析频率有关，而且最重要的是随信号本身变化而变化，因此，HHT方法非常适合处理非平稳、非线性的脑电信号。

发明内容

为了克服上述现有技术的缺点，本发明提供了基于HHT高频组合编码稳态视觉诱发电位脑机接口方法，在传统SSVEP简单刺激模式基础上提出了基于高频时序组合编码SSVEP范式，并针对变频组合编码信号的时序特征提出基于HHT特征提取方法，具有操作简单、电极数目少和目标数多的优点。

为了达到上述目的，本发明采取的技术方案为：

基于HHT高频组合编码稳态视觉诱发电位脑机接口方法，包括以下步骤：

步骤1，在受试者头部D的枕区Oz位置安放测量电极A，在受试者头部D的单侧耳垂位置安放参考电极B，在受试者头部D的前额部Fpz位置安放地电极C，测量电极A的输出端接入脑电放大器E的第一输入端E1，参考电极B的输出端接入脑电放大器E的第二输入端E2，地电极C的输出端接入脑电放大器E的第三输入端E3，脑电放大器E的输出端和计算机F的输入端相连，计算机F的输出端和计算机屏幕G连接；

步骤2，在SSVEP响应高频段，即闪烁频率高于25Hz，选择n个刺激频率，n为大于1的整数，将n个刺激频率经时序排列组合形成的nⁿ个刺激序列通过计算机屏幕G呈现在受试者面前，并对刺激序列编码形成编码库；

步骤3，受试者注视步骤2的高频时序组合编码刺激序列时，通过步骤1连接的脑电信号采集系统获得脑电信号，并采用基于希尔伯特-黄变换(HHT)的变频脑电信号特征提取方法实现脑电信号的特征提取，具体步骤如下：

◆对采集到的脑电数据进行叠加平均，提高信噪比，得到预处理的脑电数据；

◆利用带通滤波对预处理的脑电数据进行带通滤波，得到滤波后的脑电数据；

◆采用EMD方法对滤波后的脑电数据进行分解；

在使用EMD方法时，针对端点问题，采用边界预测法，边界预测法确定的极值计算公式如下：

极值＝MAX{|端点值|，|一阶近似点|}

针对停止准则问题，采用固定筛选次数，EMD算法是一个筛选的过程，筛选的数学公式如下：

$x\left(t\right)=\underset{1}{\overset{n}{\Sigma }}{c}_j=x\left(t\right)-\underset{1}{\overset{k}{\Sigma }}{m}_j+(\underset{1}{\overset{k}{\Sigma }}{m}_j-\underset{1}{\overset{p}{\Sigma }}{m}_{2j})+....$

其中：x(t)为待分解信号；

c_j为EMD分解后的IMF信号；

n为EMD分解后IMF个数；

m_j，m_2j为筛选过程通过上下包络计算得到的平均值曲线；

k，p为筛选次数；

为了保证筛选出的IMF在幅值和频率上都具有足够的物理意义，即保证相邻IMF尺度之比接近2，对筛选过程的次数k，p必须加以限制，具体根据脑电信号特点在2³-2⁴间优化选择，即选择k＝p＝l，其中

l＝相邻IMF尺度之比为2时的筛选次数；

◆根据EMD分解结果，选择目标频率所在频段的IMF重构数据；

◆使用GZC(Generalized Zero-Crossing)方法对IMF重构数据计算其瞬时频率f_gzc；

瞬时频率f_gzc的计算公式如下：

$f_{\mathrm{gzc}}=\frac{4f_1+2(f_{21}+f_{22})+f_{41}+f_{42}+f_{43}+f_{44}}{12};$

其中：$f_1=\frac{1}{{4T}_1};$

$f_{2i}=\frac{1}{{2T}_{2i}},i=\mathrm{1,2};$

$f_{4i}=\frac{1}{T_{4i}},i=\mathrm{1,2,3,4},$

其中T为零点与峰值之间的时间段；

◆根据瞬时频率f_gzc，进而得到脑电数据的Hilbert时-频图；

步骤4，将步骤3获得的相关IMF重构数据根据总时长平均分3段，分别计算出其FFT频谱，进而计算出频谱极值，将获得的频谱极值与预先设定的阈值比较并且编码，将获得的编码与步骤2中形成的编码库中的相对应的编码比较，如果在编码库中，则标记非零，否则标记为零，其中非零数为目标辨识正确次数，最终实现局部频谱极值目标辨识方法量化识别率。

由于本发明针对现有SSVEP-BCI系统的范式简单和SSVEP低频段等不足所带来的容易使受试者疲劳及诱发癫痫、目标数目少等缺点，综合系统地引入编码的思想，利用高频SSVEP的无损特点以及潜在的优势机制，提出高频时序组合编码稳态视觉诱发电位脑-机接口新范式及其变频脑电特征提取方法，解决单频率呈现目标数量受限的问题，在可用的n个基本频率单元确定后，使得可表示的目标数目由n个增加到nⁿ个，故而增加了SSVEP的可选目标数；利用高频SSVEP的刺激闪烁融合效应，提高了使用者的舒适度，降低了诱发疲劳和癫痫的可能，故而保障了高频SSVEP脑-机接口系统的高效无损功能的实现。鉴于高频SSVEP的以上诸多优点，瞄准低频SSVEP在实际应用中的不足。本发明提出高频时序组合编码SSVEP的BCI新范式，通过多个稳态高频时序排列组合实现刺激编码的思路，n个频率可以形成nⁿ个目标，解决单频率呈现目标数量受限的问题，相比n个频率可以形成nⁿ个目标的简单范式，大大扩展了目标数和选择范围；同时，变频组合编码信号具有较强的抗干扰能力，其时序特征也增加了除频率外的目标信号信息，为进一步提高辨识正确率和缩短辨识时间奠定了基础。显然，显著增加的可呈现目标数和进一步提升的辨识效率和准确率，会大大提高人机交互的通讯传输率；再加上高频SSVEP自身无损的优点，基于高频时序组合编码SSVEP的脑-机接口技术会大大提高现有BCI技术的交互能力、可靠性及使用的舒适性和可行性，可应用于宇航、军事等众多领域，具有重要的学术理论研究和实际应用价值。

附图说明

图1是本发明硬件连接示意图。

图2是本发明的刺激序列示意图。

图3是本发明实施例的目标刺激时序图及其频率图，其中图3(a)是实施例中27个目标的视觉刺激序列，图3(b)是图3(a)相应刺激序列的频率。

图4是实施例中27个目标脑电信号叠加平均后的时域波形及FFT频谱图，其中图4(a)是27个目标的SSVEP响应，图4(b)是SSVEP响应的FFT频谱。

图5是基于HHT的高频时序组合编码SSVEP脑电信号特征提取方法的流程图。

图6是边界预测法原理图。

图7是基于GZC方法计算瞬时频率的示意图。

图8是实施例中27个目标刺激序列及其高频时序组合编码SSVEP响应的Hilbert时-频图，其中：图8(a)是刺激序列Hilbert时-频图，图8(b)是SSVEP响应Hilbert时-频图。

具体实施方式

下面结合附图对本发明作进一步的详细说明。

基于HHT高频组合编码稳态视觉诱发电位脑机接口方法，包括以下步骤：

第一步，参照图1，在受试者头部D的枕区Oz位置安放测量电极A，在受试者头部D的单侧耳垂位置安放参考电极B，在受试者头部D的前额部Fpz位置安放地电极C，测量电极A的输出端接入脑电放大器E的第一输入端E1，参考电极B的输出端接入脑电放大器E的第二输入端E2，地电极C的输出端接入脑电放大器E的第三输入端E3，脑电放大器E的输出端和计算机F的输入端相连，计算机F的输出端和计算机屏幕G连接；

第二步，参照图2，在SSVEP响应高频段，即闪烁频率高于25Hz，选择n个刺激频率，n为大于1的整数，将n个刺激频率经时序排列组合形成的nⁿ个刺激序列通过计算机屏幕G呈现在受试者面前，并对刺激序列编码形成编码库；

第三步，参照图5，受试者注视第二步的高频时序组合编码刺激序列时，通过第一步连接的脑电信号采集系统获得脑电信号，并采用基于希尔伯特-黄变换(HHT)的变频脑电信号特征提取方法实现脑电信号的特征提取，具体步骤如下：

◆对采集到的脑电数据进行叠加平均，提高信噪比，得到预处理的脑电数据；

◆利用带通滤波对预处理的脑电数据进行带通滤波，得到滤波后的脑电数据；

◆采用EMD方法对滤波后的脑电数据进行分解；

在EMD方法应用时，针对端点问题，采用边界预测法.如图6所示，边界预测法确定的极值计算公式如下：

极值＝MAX{|端点值|，|一阶近似点|}；

针对停止准则问题，采用固定筛选次数，EMD算法是一个筛选的过程，描述筛选过程的数学公式如下：

$x\left(t\right)=\underset{1}{\overset{n}{\Sigma }}{c}_j=x\left(t\right)-\underset{1}{\overset{k}{\Sigma }}{m}_j+(\underset{1}{\overset{k}{\Sigma }}{m}_j-\underset{1}{\overset{p}{\Sigma }}{m}_{2j})+....$

其中：x(t)为待分解信号；

c_j为EMD分解后的IMF信号；

n为EMD分解后IMF个数；

m_j，m_2j为筛选过程通过上下包络计算得到的平均值曲线；

k，p为筛选次数；

为了保证筛选出的IMF在幅值和频率上都具有足够的物理意义，即保证相邻IMF尺度之比接近2，对筛选过程的次数k，p必须加以限制，具体根据脑电信号特点在2³-2⁴间优化选择，即选择k＝p＝l，其中

l＝相邻IMF尺度之比为2时的筛选次数；

◆根据EMD分解结果，选择目标频率所在频段的IMF重构数据；

◆使用GZC(Generalized Zero-Crossing)方法对IMF重构数据计算其瞬时频率f_gzc，如图7所示；

瞬时频率f_gzc的计算公式如下：

$f_{\mathrm{gzc}}=\frac{{4f}_1+2(f_{21}+f_{22})+f_{41}+f_{42}+f_{43}+f_{44}}{12};$

其中：$f_1=\frac{1}{{4T}_1};$

$f_{2i}=\frac{1}{{2T}_{2i}},i=\mathrm{1,2};$

$f_{4i}=\frac{1}{T_{4i}},i=\mathrm{1,2,3,4};$

其中T为零点与峰值之间的时间段；

◆根据瞬时频率f_gzc，进而得到脑电数据的Hilbert时-频图，如图8所示；

第四步，将第三步获得的相关IMF重构数据根据总时长平均分3段，分别计算出其FFT频谱，进而计算出频谱极值，将获得的频谱极值与预先设定的阈值比较并且编码，将获得的编码与第二步中形成的编码库中的相对应的编码比较，如果在编码库中，则标记非零，否则标记为零，其中非零数为目标辨识正确次数，最终实现局部频谱极值目标辨识方法量化识别率，辨识正确率计算公式如下：

下面结合具体实施例对本发明作进一步的详细说明。

基于HHT高频组合编码稳态视觉诱发电位脑机接口方法，包括以下步骤：

第一步，在受试者头部D的枕区Oz位置安放测量电极A，在受试者头部D的单侧耳垂位置安放参考电极B，在受试者头部D的前额部Fpz位置安放地电极C，测量电极A的输出端接入脑电放大器E的第一输入端E1，参考电极B的输出端接入脑电放大器E的第二输入端E2，地电极C的输出端接入脑电放大器E的第三输入端E3，脑电放大器E的输出端和计算机F的输入端相连，计算机F的输出端和计算机屏幕G连接；

第二步，参照图2和图3，在SSVEP响应高频段，即闪烁频率高于25Hz，选择25、33.33和40Hz的n＝3个频率作为基本的刺激频率单元，分别编号1、2、3，通过时序排列组合，可形成nⁿ＝3³＝27个代表不同目标的时间序列，分别为

111、112、113、121、122、123131、132、133

211、212、213、221、222、223、231、232、233

311、312、313、321、322、323、331、332、333

27个目标的视觉刺激序列如图3所示，例如，132代表25Hz、40Hz和33.33Hz的时间序列组成的刺激序列，如图2所示，每个频率采用20个脉冲触发，每个时序刺激序列时间长度在1.5s和2.4s不等，通过计算机屏幕G按照以上所设计的代表不同目标时序刺激分别重复呈现；

第三步，受试者坐在离计算机CRT屏幕G正前方约70厘米的位置，注视第二步所述的不同高频时序组合编码刺激序列目标，通过第一步连接的脑电信号采集系统获得相应的典型变频变幅非线性脑电信号，传统的频谱方法难于有效提取其特征，如图4(b)所示，为此，采用基于希尔伯特-黄变换(HHT)来实现变频脑电信号的特征提取，如图5所示，具体包括以下步骤：

◆对采集到的脑电数据进行叠加平均，提高信噪比，得到预处理的脑电数据；

◆利用带通滤波对预处理的脑电数据进行带通滤波，得到滤波后的脑电数据；

◆采用EMD方法对滤波后的脑电数据进行分解；

在EMD方法应用时，针对端点问题，采用边界预测法.如图6所示，边界预测法确定的极值计算公式如下：

极值＝MAX{|端点值|，|一阶近似点|}；

针对停止准则问题，采用固定筛选次数，EMD算法是一个筛选的过程，筛选的数学公式如下：

$x\left(t\right)=\underset{1}{\overset{n}{\Sigma }}{c}_j=x\left(t\right)-\underset{1}{\overset{k}{\Sigma }}{m}_j+(\underset{1}{\overset{k}{\Sigma }}{m}_j-\underset{1}{\overset{p}{\Sigma }}{m}_{2j})+....$

其中：x(t)为待分解信号；

c_j为EMD分解后的IMF信号；

n为EMD分解后IMF个数；

m_j，m_2j为筛选过程通过上下包络计算得到的平均值曲线；

k，p为筛选次数；

为了保证筛选出的IMF在幅值和频率上都具有足够的物理意义，即保证相邻IMF尺度之比接近2，对筛选过程的次数k，p必须加以限制，具体根据脑电信号特点在2³-2⁴间优化选择，即选择k＝p＝l，其中

l＝相邻IMF尺度之比为2时的筛选次数；

◆根据EMD分解结果，选择目标频率所在频段的IMF重构数据；

◆使用GZC(Generalized Zero-Crossing)方法对IMF重构数据计算其瞬时频率f_gzc，如图7所示；

瞬时频率f_gzc的计算公式如下：

$f_{\mathrm{gzc}}=\frac{{4f}_1+2(f_{21}+f_{22})+f_{41}+f_{42}+f_{43}+f_{44}}{12}$

其中：$f_1=\frac{1}{{4T}_1};$

$f_{2i}=\frac{1}{{2T}_{2i}},i=\mathrm{1,2};$

$f_{4i}=\frac{1}{T_{4i}},i=\mathrm{1,2,3,4},$

T为零点与峰值之间的时间段；

◆根据瞬时频率，进而得到脑电数据的Hilbert时-频图，如图8所示；

第四步，将第三步获得的相关IMF重构数据根据总时长平均分3段，分别求其FFT频谱，并计算频谱极值，假设x为频谱极值，

将获得的频谱极值按照如下阈值设定分别编码：

a)20＜x＜＝29，x编码为1；

b)29＜x＜＝37，x编码为2；

c)37＜x＜＝37，x编码为3；

最终将获得的编码与第二步中形成的编码库中的元素比较，如果在编码库中，标记非零，否则标记为零，其中非零数为目标辨识正确次数，最终实现局部频谱极值目标辨识方法量化识别率。辨识正确率计算公式如下：

最后的辨识率为77.78％。

附图中：

A为测量电极；B为参考电极；C为地电极；D为受试者头部；E为脑电放大器；F为计算机；G为计算机屏幕；E1为第一输入端；E2为第二输入端；E3为第三输入端；

■代表端点；

●代表一阶近似点：靠近端点的两个极值点直线延拓和通过端点纵轴的交点；

▲代表边界预测法确定的最后的极值点；

T₁；

T_2i，i＝1，2；

T_4i，i＝1，2，4，5；以上三个T均为零点与峰值之间的时间段。