CN103690160A - 一种基于非高斯时序模型的脑电特征提取方法 - Google Patents

Abstract

本发明公开了一种基于非高斯时序模型的脑电特征提取方法，包括以下步骤：获取待处理脑电数据和两组训练脑电数据，去除伪迹，将得到的有效频段分为若干数据段；提取每个数据段的时频特征值、形态特征值和复杂度特征值，每个数据段的特征值构成一个特征向量；为第一组训练脑电数据中的每个特征向量标记状态值，并利用标记结果训练支持向量机；将第二组训练脑电数据的特征向量输入支持向量机中，得到第二组训练脑电数据的状态值序列；建立观察方程和状态转移方程，利用第二组训练脑电数据的特征向量和状态值序列，确定方程中的参数；利用待处理脑电数据的特征向量以及两个方程，得到待处理脑电数据的状态值。本方法能够准确辨别不同的大脑状态。

Description

一种基于非高斯时序模型的脑电特征提取方法

技术领域

本发明涉及脑电数据分析领域，具体涉及一种基于非高斯时序模型的脑电特征提取方法。

背景技术

头皮脑电信号数据包含了大量与大脑特征以及状态相关的信息，是进行大脑状态判断的重要工具。由于大脑的时变性和非线性性，其产生的信号也是时变的、非线性的，同时由于头皮脑电信号容易受到眼电、肌电以及心电信号的干扰，有效而鲁棒的脑电信号特征提取方法成为难题。

申请号为200910196746.3的发明公开了一种脑电波分析方法，应用了经典的时频域分析和主成份分析方法解决了脑电信号特征提取的难题，成功地提取了与人体紧张，疲劳与放松息息相关的时频域参数，并且将其映射到主成分空间中去，还运用了支持向量机在主成份空间中高效地分析非线性关系，提高了判读的精确性和有效性。但是，该发明能够处理的频率范围有限，不能满足更宽频段的需求，而且对于脑电信号中的噪声信号没有进行有效的处理。

申请号为200710163302.0的发明公开了一种脑电超慢波分析方法及其数据处理装置，它通过对0.5-50Hz范围内的脑电波进行系列分析，提取出其中超慢波，并分析超慢谱系功率，用以反映所对应的脑神经递质的功能状况，为建立一种无创伤脑神经递质功能检测方法提供了条件。但是，该发明对于脑电数据特征具有针对性，当特征值不一致时效果较差。

通过对脑电数据信号的监测和分析，可以得到大脑的状态，如深度睡眠，大脑疲劳，大脑警觉等状态，同时，也为人类与计算机的交互控制提供了新的手段。

对脑电数据信号进行特征提取，得到可识别的，可靠程度高的特征值是脑电数据分析的基础环节，具有重要价值，因此，需要提供一种脑电特征的提取方法，提取得到能够反应大脑状态变化信息的特征值。

发明内容

本发明提供了一种基于非高斯时序模型的脑电特征提取方法，结合脑电信号的时频特征、形态特征以及复杂度特征等多种特征，对大脑状态进行有效刻画，同时通过状态空间模型对大脑状态进行实时跟踪，有效去除伪迹，能够有效反应大脑状态的变化，敏感性高，抗噪性好。

一种基于非高斯时序模型的脑电特征提取方法，包括以下步骤：

（1）获取待处理脑电数据和两组训练脑电数据，去除待处理脑电数据和两组训练脑电数据中的伪迹，分别获得待处理脑电数据的有效频段和两组训练脑电数据的有效频段，再分别将待处理脑电数据的有效频段和每组训练脑电数据的有效频段分为若干数据段；每组训练脑电数据中均包括两种大脑状态下的脑电数据。

不同的大脑状态对应具有不同特征的脑电数据，例如深度睡眠时，脑电数据波动较为平缓；剧烈运动时，脑电数据波动很大。

作为优选，每组训练脑电数据中的两种大脑状态所对应的时长相同，且每种大脑状态所持续的时间不少于1分钟。

例如，需要区分的两种大脑状态分别为深度睡眠和剧烈运动时的大脑状态，深度睡眠和剧烈运动的时长分别为2分钟。

作为优选，所述步骤（1）采用带通滤波去除伪迹，去除伪迹后获得的有效频段频率为1.6～70Hz。带通滤波可以采用2阶巴特沃兹滤波器（Butterworth），有效频段的频率可以根据需要通过滤波参数进行选择。

步骤（1）中数据分段采用滑动时间窗的方法，时间窗的长度为1s，步长为0.5s。数据分段的目的在于确定数据处理的最小单元。滑动时间窗的长度以及滑动步长可以根据需要选取，时间窗的长度越短，滑动步长越短，划分的数据段越多，则提取得到的时频特征值和形态特征值越精确，但是，相应的计算量也大。

（2）提取步骤（1）中每个数据段的时频特征、形态特征和复杂度特征，得到相应的时频特征值、形态特征值和复杂度特征值，每个数据段的时频特征值、形态特征值和复杂度特征值构成一个特征向量。

步骤(2)中的时频特征提取采用经验模态分解，取前三个固有模态函数，利用下式计算得到每个数据段的时频特征值VoIMF：

${\mathrm{VoIMF}}_M\left(t\right)=\frac{1}{n-1}{\mathrm{\Σ}}_{i=1}^n{(y_i-\stackrel{\‾}{y})}^2,M=\mathrm{1,2,3};$

其中，y_i为每个数据段中每个数据点的固有模态函数值；

为每个数据段中所有数据点的固有模态函数平均值；

i为每个数据段数据点的序数；

n为每个数据段数据点的个数；

M为固有模态函数x的序数；

t为数据段的序数。

脑电数据经过经验模态分解(Empirical Mode Decomposition，简称EMD)，被转化为包含了原脑电数据信号的不同时频尺度特征的固有模态函数(Intrinsic Mode Function，简称IMF)，能够反映脑电数据在不同时频尺度上的特征。

一般来说，序数较小的固有模态函数代表频率较高的信号分量，M可以根据需要选择所需的固有模态函数的数目。当M为3时，即选取前三个固有模态函数，分别得到VoIMF₁、VoIMF₂和VoIMF₃。

对于每一个数据段分别计算VoIMF₁、VoIMF₂和VoIMF₃。例如计算VoIMF₁时，选取第一个固有模态函数，计算该数据段中每个数据点的第一个固有模态函数值y_i，计算该数据段中所有数据点的第一个固有模态函数值y_i的平均值得到

利用式

计算得到VoIMF₁，式中，i为每个数据段的数据点的序数；n为每个数据段的数据点的个数；t为数据段的序数。

步骤(2)中形态特征提取的步骤如下：

a、对每个数据段进行均值滤波，得到平滑后的数据；

b、求取平滑后的数据段的极大值，连接极大值得到上包络线E_upper；求取平滑后的数据段的极小值，得到下包络线E_lower；利用下式求取上包络线和下包络线的之间的包络范围：

Envelope_Range(t)＝E_upper(t)-E_lower(t)

其中，t为数据段序数；

Envelope_Range(t)为第t个数据段的包络范围；

利用下式计算每个数据段的形态特征值VoE，

$\mathrm{VoE}\left(t\right)=\frac{1}{n-1}{\mathrm{\Σ}}_{i=1}^n{({\mathrm{Envelope}\_\mathrm{Range}\left(t\right)}_i-\stackrel{\‾}{\mathrm{Envelope}\_\mathrm{Range}\left(t\right)})}^2;$

其中，t为数据段的序数；

n为每个数据段数据点的个数；

i为每个数据段中数据点的序数；

Envelope_Range(t)_i为第t个数据段第i个数据点的包络范围值;

为第t个数据段所有数据点包络范围的平均值。

利用均值滤波（即线性滤波），对每个数据段进行平滑，去除频率较高的极值，得到平滑过的数据段，均值滤波可以采用滑动时间窗平均法，时间窗的长度和滑动步长可以根据需要进行选择。

优选地，所述复杂度特征值为样本熵。步骤（2）中复杂度特征提取采用样本熵的方法，样本熵（Sample Entropy，简称SampEn）是时间序列复杂度的一种度量，相比于近似熵而言（Approximate Entropy，简称ApEn），样本熵的计算简单，且与数据长度无关，已在生物信号序列分析中被广泛应用。

（3）两组训练脑电数据分别记为第一组训练脑电数据和第二组训练脑电数据：

3-1、为第一组训练脑电数据中的每一个特征向量标记状态值，利用标记完状态值的第一组训练脑电数据训练支持向量机；所述支持向量机采用线性核，训练结束后得到最优分类超平面（Optimal SeperatingHyperplane）；

3-2、将第二组训练脑电数据的所有特征向量输入步骤3-1中所得的训练好的支持向量机中，得到第二组训练脑电数据的状态值序列。

状态值序列中的每一个状态值对应一个特征向量到最优分类超平面的距离（Decision Value），由于支持向量机在训练过程中能够得到两类大脑状态数据的最优分类超平面，因此，可以对两类大脑状态进行区分。

(4)建立表达特征向量与大脑状态之间关系的观察方程，并利用自回归模型建立状态转移方程，利用第二组训练脑电数据的特征向量和状态值序列，确定观察方程和状态转移方程中的所有参数。

优选地，所述步骤(4)中采用最小二乘法，确定观察方程和状态转移方程中的所有参数。

所述状态转移方程描述当前状态和之前状态的关系，一般认为大脑状态变化是一个相对平缓的过程，因此，基于自回归模型(AutoregressiveModel，简称AR模型)建立状态转移方程。

作为优选，所述状态转移方程的表达式如下： $x_t=\underset{i=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\α}}_i{x}_{t-i}+\mathrm{\β}+v_t,v_t~N(0,{\mathrm{\σ}}_v^2)$

其中；

t为数据段的序数；

x_t为第t个数据段的状态值；

p为自回归模型的阶数；

x_t-i为第t-i个数据段的状态值；

α_i为状态值x_t-i对应的系数；

β为常数；

v_t表示第t个数据段的状态转移噪声，服从均值为0，方差为

的正态分布。

p取1就基本能够反映大脑状态的变化，也可以根据不同的状态转移情况选择更高的阶数。

步骤(3)中的观察方程描述特征向量与大脑状态值之间的关系，由于大脑状态的复杂性、非线性性，因此，建立基于多项式模型的观察方程。多项式模型能够表达大脑状态值与特征向量之间的复杂关系，同时，特征向量中也包含外部噪声，这些噪声由其他生理信号(如肌电，脑电，眼电等)或者电极脱落等原因导致，一般具有突发性，使得脑电信号中出现具有高幅度或高频率的伪迹，称为“脉冲噪声”，该类噪声不符合高斯分布(即正态分布)，一般需要使用重尾噪声进行描述，因此，观察方程中的观察噪声采用重尾的柯西分布来刻画。

作为优选，所述观察方程的表达式如下：

f_t＝Wx′_t+u_t，u_t～C(0，Scale)

其中，

t为数据段的序数；

f_t表示第t个数据段的特征向量；

表示第t个数据段的状态值x_t的i阶乘方，q为多项式模型的阶数；

W为系数矩阵；

u_t为第t个数据段的观察噪声，服从中心为0，尺度为Scale的柯西分布。

(5)利用待处理脑电数据的特征向量以及步骤(4)获得的观察方程和状态转移方程，采用粒子滤波的方法得到待处理脑电数据的状态值。

所述状态转移方程和观察方程共同构成状态空间模型，由于状态空间模型非线性、非高斯，通过粒子滤波的方法能够对状态值进行有效的估计，粒子滤波方法基于时序蒙特卡洛的思想，通过粒子集合以及每个粒子的权重，对状态值的后验概率进行估计，即得到待处理脑电数据的状态值。

本发明将脑电信号的时频特征、形态特征以及复杂度特征相结合，能够准确快速地识别脑电数据的特征变化，检测大脑的状态。

附图说明

图1为待处理脑电数据图；

图2为利用本发明方法处理后的状态值；

图3为本发明方法中提取得到的时频特征值、形态特征值以及复杂度特征值。

具体实施方式

下面结合附图，对本发明基于非高斯时序模型的脑电特征提取方法做详细描述。

一种基于非高斯时序模型的脑电特征提取方法，包括以下步骤：

(1)获取待处理脑电数据和两组训练脑电数据，去除待处理脑电数据和两组训练脑电数据中的伪迹，分别获得待处理脑电数据的有效频段和两组训练脑电数据的有效频段，再分别将待处理脑电数据的有效频段和每组训练脑电数据的有效频段分为若干数据段；每组训练脑电数据中均包括两种大脑状态下的脑电数据。

每组训练脑电数据中包含两种大脑状态下的脑电信号数据，每种大脑状态的信号连续且时长为2分钟。

获取待处理脑电数据和两组训练脑电数据后，对脑电数据进行预处理，将原始的脑电数据通过带通滤波去除伪迹，选取有效频段进行分析和处理，带通滤波采用2阶巴特沃兹滤波器，滤波参数为1.6～70Hz，即得到的有效频段频率为1.6～70Hz。

将得到的有效频段均匀分为若干数据段，确定数据处理的最小单元，采用滑动时间窗的方法对得到的有效频段进行分解，滑动时间窗的长度为1s，滑动步长为0.5s。

例如，采样频率为1000Hz，滑动时间窗的长度为1s，则每个时间窗内包含1000个数据点，滑动步长为0.5s，即下一个时间窗的起点为当前时间窗起点后的0.5s，即当前时间窗的第501个数据点为下一个时间窗的第1个数据点。

两组训练脑电数据分别记为第一组训练脑电数据和第二组训练脑电数据，将第一组训练脑电数据中的第一种大脑状态下的训练脑电数据记为segment1，第二种大脑状态下的训练脑电数据记为segment2。

将segment1划分为m个数据段，得到数据段集合

其中表示segment1被划分后的第i个数据段；同理，将segment2划分为b个数据段，得到数据段集合其中表示segment2被划分后的第i个数据段。

本申请中涉及数据段的序数时，不同组的脑电数据独立编号，例如第二组训练脑电数据划分为T段，则第二组脑电数据划分得到的数据段依次编号为1～T；待处理脑电数据划分为Q段，则待处理脑电数据划分得到的数据段依次编号为1～Q。

（2）提取步骤（1）中每个数据段的时频特征、形态特征和复杂度特征，得到相应的时频特征值、形态特征值和复杂度特征值，每个数据段的时频特征值、形态特征值和复杂度特征值构成一个特征向量。

提取时频特征采用现有的经验模态分解方法，首先，将数据段分解成为具有不同时频特征的固有模态函数，再通过计算固有模态函数的方差得到时频特征值。

经验模态分解步骤如下：

首先，求取第t个数据段中所有极大值点，并利用三次样条插值函数拟合形成上包络线；求取第t个数据段中所有的极小值点，并利用三次样条插值函数拟合形成下包络线；计算上包络线和下包络线的平均值，得到第t个数据段的平均数据序列；

其次，第t个数据段的原始数据序列减去平均数据序列，得到第t个数据段的最终数据序列；

最后，若得到的第t个数据段的最终数据序列仍然存在负的极大值和正的极小值，说明这还不是一个本征模函数，重复前述步骤，直至最终数据序列满足固有模态函数要求。

一般来说，序数较小的固有模态函数代表频率较高的信号分量，根据数据处理经验，选取前三个固有模态函数作为有效时频特征，利用下式计算得到时频特征值； ${\mathrm{VoIMF}}_M\left(t\right)=\frac{1}{n-1}{\mathrm{\Σ}}_{i=1}^n{(y_i-\stackrel{\‾}{y})}^2,M=\mathrm{1,2,3};$

其中，y_i为每个数据段中每个数据点的固有模态函数值；

为每个数据段中所有数据点的固有模态函数平均值；

i为每个数据段的数据点的序数；

n为每个数据段的数据点的个数；

M为固有模态函数y的序数；

t为数据段的序数。

M为固有模态函数y的序数，这里M取3，即可以得到3个时频特征值，分别为VoIMF₁、VoIMF₂、VoIMF₃。

对每个数据段均做时频特征的提取，得到相应的时频特征值VoIMF₁、VoIMF₂、VoIMF₃。

对每个数据段提取形态特征的步骤如下：

A.采用时间窗平均法，时间窗长度为1s，步长为0.5s；将数据窗中点的数据值赋值为时间窗均值，再通过三次样条插值方法得到平滑过的数据段。

B.求取平滑后的数据段的极大值，连接极大值得到上包络线E_upper；求取平滑后的数据段的极小值，得到下包络线E_lower；利用下式求取上包络线和下包络线的之间的包络范围：

Envelope_Range(t)＝E_upper(t)-E_lower(t)

其中，t为数据段序数；

Envelope_Range(t)为第t个数据段的包络范围；

利用下式计算每个数据段的形态特征值VoE，

$\mathrm{VoE}\left(t\right)=\frac{1}{n-1}{\mathrm{\Σ}}_{i=1}^n{({\mathrm{Envelope}\_\mathrm{Range}\left(t\right)}_i-\stackrel{\‾}{\mathrm{Envelope}\_\mathrm{Range}\left(t\right)})}^2;$

其中，t为数据段的序数；

n为每个数据段数据点的个数；

i为每个数据段中数据点的序数；

Envelope_Range(t)_i为第t个数据段第i个数据点的包络范围值。复杂度特征采用现有的样本熵作为度量，对于每个数据段，样本熵计算步骤如下(样本熵的计算参考文献J.Richman and R.Moorman，“Physiologicaltime-series analysis using approximate entropy and sample entropy,”AmericanJournal of Physiology-Heart and Circulatory Physiology,vol.278，no.6，pp.H2039-H2049，2000.)：

A、对于包含了n个数据点的数据段{a₁，a₂,…，a_n}，其中a_i表示第i个数据点，取长度为h的模板向量x_h(i)＝(x_i，x_i+1,…x_i+h-1)，其中i表示模板向量序数，并定义两个向量之间的距离d[x_h(i)，x_h(j)]为切比雪夫距离；

B、定义相似容限r，计算每两个向量之间的距离，将向量间距离小于r的向量的对数记为A；

C、设步骤A中模板向量的长度为h+1，重复步骤A和B，将模板长度h+1时，向量间距离小于r的向量的对数记为B；

D、利用下式计算样本熵SampEn的值：

$\mathrm{SampEn}\left(t\right)=-\log \frac{A}{B}$

其中，t为数据段的序数；

SampEn(t)为第t个数据段的样本熵的值。

样本熵计算中参数设置为：模板向量的长度h＝2，相似容限r＝0.2。

每个数据段提取得到三个时频特征值，一个形态特征值和一个复杂度特征值，五个特征值构成一个五维的特征向量，对segment1进行特征向量的提取，得到segment1的特征向量集

其中，

对应segment1中的第i个数据段，同理，对segment2进行特征向量的提取，得到segment2的特征向量集

对应segment2中的第i个数据段。

(3)3-1、为第一组训练脑电数据中的每一个特征向量标记状态值，利用标记完状态值的第一组训练脑电数据训练支持向量机；

为segment1和segment2中的每个特征向量标记状态，第一种状态标记为1，第二种状态标记为0，标记完毕后，segment1的特征向量集为

segment2的特征向量集为

将这两个特征向量集输入线性支持向量机中，进行支持向量机的训练，训练好的支持向量机可以得到最优分类超平面。

3-2、将第二组训练脑电数据的所有特征向量输入步骤3-1中所得的训练好的支持向量机中，得到第二组训练脑电数据的状态值序列。

对第二组训练脑电数据进行特征提取，得到特征向量集

将特征向量集

输入训练好的支持向量机中，得到状态值序列

状态值序列中的每一个状态值对应一个特征向量与最优分类超平面的距离(Decision Value)。

(4)建立表达特征向量与大脑状态之间关系的观察方程，并利用白回归模型建立状态转移方程，利用第二组训练脑电数据的特征向量和状态值序列，确定观察方程和状态转移方程中的所有参数。

状态转移方程的表达式如下： $x_t=\underset{i=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\α}}_i{x}_{t-i}+\mathrm{\β}+v_t,v_t~N(0,{\mathrm{\σ}}_v^2)$

其中，

t为数据段的序数；

x_t为第t个数据段的状态值；

p为自回归模型的阶数；

x_t-i为第t-i个数据段的状态值；

α_i为状态值x_t-i对应的系数；

β为常数；

v_t表示第t个数据段的状态转移噪声，服从均值为0，方差为的正态分布。

观察方程的表达式如下：

f_t＝Wx′_t+u_t，u_t～C(0，Scale)

其中，

t为数据段的序数；

f_t表示第t个数据段的特征向量，该特征向量即为步骤(2)中提取得到的五维的特征向量；

x′_t为一个q维向量，其中第i维(即)代表第t个数据段的状态值x_t的i-1阶乘方，q为多项式模型的阶数；

W为一个5*q的系数矩阵；

u_t为第t个数据段的观察噪声，服从中心为0，尺度为Scale的柯西分布，第i维u_ti服从中心为0，尺度为Scale i的柯西分布，其中Scale i为特征向量的第i维。

状态值序列

中的每一个状态值即状态转移方程中的x_t，采用最小二乘法，计算状态转移方程中的系数α_i和β，然后利用下式计算状态转移方程的拟合残差

$\widehat{e_i}=x_i-(\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\α}}_j{x}_{i-j}+\mathrm{\β})$

其中，i是数据段的序号；

x_i是第i个数据段的状态值，即

p为状态转移方程的阶数。

利用下式，计算得到状态转移噪声v_t的标准差σ_v：

${\mathrm{\σ}}_v=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(\widehat{e_i}-\stackrel{\‾}{e})}^2}$

其中，

i是数据段的序号；

是第i个数据段的拟合残差；

是拟合残差的平均值；

N为数据段的个数。

利用状态值序列

和相应的特征向量集

通过最小二乘法拟合得到系数矩阵W；

对每一个特征值(特征值包括三个时频特征值、一个形态特征值以及一个复杂度特征值)，分别利用下式计算拟合残差

${\hat{e}}_i^k=f_i^k-\underset{j=0}{\overset{q-1}{\mathrm{\Σ}}}{w}_j{x}_i^j$

其中，

为第i个数据段的特征向量

的第k维；

x_i是第i个数据段的状态值，即

为x_i的j次方；

w_i为系数矩阵W的第k行；

q为多项式模型的阶数；

利用拟合残差集合训练观察噪声u_t，得到柯西噪声u_t的尺度Scale(计算方法参考文献F.Nagy,“Parameter Estimation of the Cauchy DistributioninInformation Theory Approach,”Journal of Universal Computer Science.vol.12，no.9，1332-1344，2006.)

(5)利用待处理脑电数据的特征向量以及步骤(4)获得的观察方程和状态转移方程，采用粒子滤波的方法(计算方法参考文献M.SanjeevArulampalam，S.Maskell，N.Gordon，and T.Clapp，“Atutorial onparticle filters for online nonlinear/non-Gaussian Bayesiantracking,”SignalProcessing,IEEE Transactions on,vol.50,no.2,pp.174-188,2002.）得到待处理脑电数据的状态值。本实施例中的粒子滤波方法中的粒子总数量为400。

本实施例测试一组头皮脑电数据（如图1所示），以检测脑电信号中的高频节律为例，对基于非高斯时序模型的脑电特征提取方法进行测试。

在特征提取阶段采用了3个时频特征、形态特征以及复杂度特征3种共5个特征值对脑电信号进行分析。

图2中共有25秒头皮脑电数据，其中前十秒为静息状态，从第10秒开始，脑电中出现高频节律信号，即要被检测出的状态，由图3可以看到，五个特征值在10秒以后都呈现上升趋势，表明所有特征都能有效表征大脑的两种状态。

图2中通过本发明提供的基于非高斯时序模型的脑电特征提取方法，对大脑状态进行估计，依据获得的状态值能够较好的区分两种状态。

Claims (7)

1.一种基于非高斯时序模型的脑电特征提取方法，其特征在于，包括以下步骤：

（1）获取待处理脑电数据和两组训练脑电数据，去除待处理脑电数据和两组训练脑电数据中的伪迹，分别获得待处理脑电数据的有效频段和两组训练脑电数据的有效频段，再分别将待处理脑电数据的有效频段和每组训练脑电数据的有效频段分为若干数据段；每组训练脑电数据中均包括两种大脑状态下的脑电数据；

（2）提取步骤（1）中每个数据段的时频特征、形态特征和复杂度特征，得到相应的时频特征值、形态特征值和复杂度特征值，每个数据段的时频特征值、形态特征值和复杂度特征值构成一个特征向量；

（3）两组训练脑电数据分别记为第一组训练脑电数据和第二组训练脑电数据：

3-1、为第一组训练脑电数据中的每一个特征向量标记状态值，利用标记完状态值的第一组训练脑电数据训练支持向量机；

3-2、将第二组训练脑电数据的所有特征向量输入步骤3-1中所得的训练好的支持向量机中，得到第二组训练脑电数据的状态值序列；

（4）建立表达特征向量与大脑状态之间关系的观察方程，并利用自回归模型建立状态转移方程，利用第二组训练脑电数据的特征向量和状态值序列，确定观察方程和状态转移方程中的所有参数；

（5）利用待处理脑电数据的特征向量以及步骤（4）获得的观察方程和状态转移方程，采用粒子滤波的方法得到待处理脑电数据的状态值。

2.如权利要求1所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，所述状态转移方程的表达式如下：

$x_t=\underset{i=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\α}}_i{x}_{t-i}+\mathrm{\β}+v_t,v_t~N(0,{\mathrm{\σ}}_v^2)$

其中，

t为数据段的序数；

x_t为第t个数据段的状态值；

p为自回归模型的阶数；

x_t-i为第t-i个数据段的状态值；

α_i为状态值x_t-i对应的系数；

β为常数；

v_t表示第t个数据段的状态转移噪声，服从均值为0，方差为

的正态分布。

3.如权利要求2所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，所述观察方程的表达式如下：

f_t＝Wx′_t+u_t，u_t～C(0，Scale)

其中，

t为数据段的序数；

f_t表示第t个数据段的特征向量；

表示第t个数据段的状态值x_t的i阶乘方，q为多项式模型的阶数；

W为系数矩阵；

u_t为第t个数据段的观察噪声，服从中心为0，尺度为Scale的柯西分布。

4.如权利要求3所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，每组训练脑电数据中的两种大脑状态所对应的时长相同，且每种大脑状态所持续的时间不少于1分钟。

5.如权利要求4所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，所述步骤(4)中采用最小二乘法，确定观察方程和状态转移方程中的所有参数。

6.如权利要求5所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，所述步骤(1)采用带通滤波去除伪迹，去除伪迹后获得的有效频段频率为1.6～70Hz。

7.如权利要求6所述的基于非高斯时序模型的脑电特征提取方法，其特征在于，所述复杂度特征值为样本熵。

Images