基于二维时频图像深度卷积神经网络的癫痫脑电识别方法与流程

本发明涉及脑电信号分类领域，具体来说是一种通过对癫痫脑电信号进行特征提取与分类来对癫痫发作进行识别的方法。

背景技术：

近年来，脑电信号(EEG)已经发展成为癫痫研究的一种主要手段。传统的通过分析EEG信号的癫痫识别方法已有了较为固定的模式，首先是提取EEG信号中的有效频段，再对提取出的有效频段信号进行特征提取，这些特征通过是人工选取的，如样本熵、小波熵、多尺度熵、多尺度LZ复杂度等用于反映信号非线性程度的一些指标，之后再将提取的特征通过支持向量机或神经网络等分类器进行分类来完成对癫痫的识别。

目前，已有许多学者对上述的癫痫识别方法进行了研究，如汪春梅等采用近似熵+NEYMAN-PEARSON准则^[1]进行分类；黄丽亚等采用多尺度样本熵+支持向量机^[2]的形式。这些分类方法固然取得了一定的效果，但这些研究方法中用于分类的特征均是通过人工来选取，存在着一定的随机性，无法完全表现出癫痫发作状态下与未发作状态下的差异，可能无法适用于所有的EEG样本的癫痫识别。所以，如何提取出癫痫发作时与未发作时EEG信号具有最大差异性的特征，是目前研究的一个重要方向。

参考文献：

[1]汪春梅,张崇明,王丽慧.一种癫痫脑电信号分类检测装置及方法:CN102429657A[P].2012.

[2]黄丽亚,郭迪,沈洋洋.基于多尺度样本熵的癫痫脑皮层电图信号分类方法:CN105046273A[P].2015.

技术实现要素：

本发明旨在提供一种新型的通过对脑电信号进行特征提取与分类的癫痫发作状态识别方法。本发明先利用小波变换提取出脑电信号0～32Hz的有效频段信号，再通过时频分析方法中的短时傅里叶变换，将提取出的0～32Hz频段的一维脑电信号转化为二维时频图像，再利用LeNet-5网络结构的深度卷积神经网络对时频图像进行特征提取与分类，接着根据部分样本数据测试出具有患者特异性的最优的五个识别通道并计算权重，最后采用加权求和的方法结合最优的五个通道进行癫痫状态的识别。技术方案如下：

一种基于二维时频图像深度卷积神经网络的癫痫脑电识别方法，包括下列步骤：

步骤1：对原始脑电信号的预处理

分成多段数据，为每段数据设置相应的标签，标签一共分为两类，第一类为发作状态，设置为1；第二类为非发作状态，设置为0；再对各段数据分别进行低通滤波处理，去除高于32Hz的信号。

步骤2：脑电信号的有效频段提取

通过小波变换进行5层分解，提取出经过预处理的脑电信号中多个有效频段的信号，再将提取出的多个信号进行叠加，合成得到0～32Hz频段的一维脑电信号，此频段的脑电信号能反映出患者癫痫发作时的有效信息。

步骤3：脑电信号的时频分析

对合成后的0～32Hz频段的一维脑电信号通过短时傅里叶变换对有效频段的脑电信号进行变换，使一维时域信号，转化为包含时间和频率信息的二维时频图像，时频图像纵轴的频率范围为0～32Hz，横轴的时间范围为窗口长度，所得到的时频图像在癫痫发作状态与非癫痫发作状态具有不同的图像特征，反映不同的时频信息；

步骤4：利用时频图训练深度卷积神经网络

将对某患者的经过步骤1-3处理得到的二维时频图像分为训练数据和测试数据，建立LeNet-5结构的深度卷积神经网络，输入到深度卷积神经网络中对其进行训练，对图像进行特征提取，通过全连接的网络进行数据降维，最终输出用于表示分类结果的二维向量，二维向量为(0,1)则表示癫痫发作状态，(1,0)则表示非癫痫发作状态，将患者的脑电数据的各个通道，均通过这种方式来训练对应的网络参数。

步骤5：选出患者特定的最优的五个通道，并计算权重

根据患者所有通道测试数据的识别结果，选出准确率排在前五的通道，并用其准确率构建该五个通道的权重，这五个准确率相对较高的通道用于最终的癫痫发作识别通道。

步骤6：采用加权求和的方法结合最优的五个通道进行癫痫的识别

进行实时癫痫发作状态识别时，对所采集的原始脑电信号进行步骤1-3的处理，再根据上述步骤选出的五个通道及其权重，将五个通道的识别结果加权求和，所得结果与癫痫发作标签值进行比较，如果所得结果的标签值与癫痫发作标签值更接近，便可判断患者处于癫痫发作状态，反之如果所得结果的标签值与非癫痫发作标签值更接近，则判断患者处于非癫痫发作状态。

本发明提出的用于癫痫脑电实时识别的基于二维时频图像特征的深度卷积神经网络方法，采用时频分析将一维脑电信号转化为二维的时频图像，再通过深度卷积神经网络对时频图像进行分类，从而判别出患者是否处于癫痫发作状态。本发明打破了传统的基于一维脑电信号的人工特征提取并用SVM等进行分类的方法。而采用通过时频分析将脑电信号的有效信息反映在时频图像上，相应的特征通过深度卷积神经网络进行自动提取与分类，而非人工地指定需要提取的特征。在最大程度保留脑电信号有效信息的情况下，充分发挥深度卷积神经网络的优势，提高癫痫发作状态识别的准确率。同时对多个通道通过准确率进行删选，对选出的五个相对最优通道，采用加权求和的方式来综合五个通道的分类结果，使得多个通道的信息得以充分利用。克服了以单个通道进行识别造成其他通道信息浪费的缺点，并使结果的准确率更高。

附图说明

图1：某位患者某次发作前后的某通道脑电信号4分钟数据时域图

图2：分段后的脑电信号数据的标签设置图

图3：脑电信号有效频段提取的小波分解示意图

图4a、图4b、图4c：一位患者某次发作前的单通道delta波时域图

图5a、图5b、图5c：一位患者某次发作时的单通道delta波时域图

图6a、图6b、图6c：一位患者某次发作后的单通道delta波时域图

图7a、图7b、图7c：一位患者某次发作前的单通道时频图

图8a、图8b、图8c：一位患者某次发作时的单通道时频图

图9a、图9b、图9c：一位患者某次发作后的单通道时频图

图10：本发明采用的LeNet-5的深度卷积网络模型

图11：本发明采用的LeNet-5的深度卷积网络实际结构

具体实施方式

现结合实施例、附图对本发明作进一步描述：

步骤1：脑电信号的预处理

脑电信号的原始数据以PhysioNet网站的CHB-MIT数据库为例来进行说明。该数据库共有23个患者，患者的年龄从3岁到22岁不等，每个患者均在停药的情况下，同时监测23个通道，共监控了46个小时。同时，每位患者每次癫痫发作的开始时间和结束时间均已知。

该数据库每个通道脑电信号的采样频率均为256Hz，其中一位患者的某次发作前后的某个通道一维脑电信号4分钟数据时域图如图1所示，其中，两条红线之间为患者癫痫发作的时间段。为了模拟实时采集的情况，以4秒长数据的窗口，每次移动1秒长的数据，将原始脑电数据分割成多段的脑电信号。同时为了更好地进行分类，每一段脑电信号数据设置相应的标签。由于事先已知患者每次发作的开始时间和结束时间，因此对于每段数据，标签设置方法如图2所示。其中，标签一共分为两类，第一类为发作状态，设置为1；第二类为非发作状态，设置为0。只要4秒长的数据的窗口部分或完全处于癫痫发作时间段，则将该段脑电信号的数据设为发作状态，否则设为非发作状态。

由于原始的脑电信号中包含着各类噪声，如50Hz工频干扰、肌电干扰等，因此在进行特征提取之前，需要对脑电信号进行滤波。根据前人的研究，脑电信号的有效频段约为32Hz之内，因此，针对每一段4秒长的脑电信号的数据，通过参数如表1的低通滤波器进行低通滤波，滤除32Hz之后的信号，以减少噪声的干扰。

表1

步骤2：脑电信号的有效频段提取

为了更好地提取出0～32Hz频段的信号，对于每一段低通滤波后的4秒长数据，采用小波变换进行五层分解，母小波选择与脑电信号时域形状较为相近的db5小波。小波分解的过程如图3所示，最后得到的脑电信号名称以及对应波段的名称如表2所示：

表2

以delta为例，患者某次发作前后的单通道时域图如图4a、4b、4c、5a、5b、5c、6a、6b、6c所示，为了重构0～32Hz频段的脑电信号，需要根据这四个频段的系数，分别重构出对应的信号，再将这四种信号重新叠加到一起，由此便得到0～32Hz的脑电信号数据。

步骤3：脑电信号的时频分析

对上述得到的0～32Hz频段的一维脑电信号，进行短时傅里叶变换(STFT)，得到对应脑电信号与时间、频率相关的信息。其中，STFT的变换参数如表3所示：

表3

由于STFT变换后，得到为时间、频率、以及与时间和频率相对应的幅值，为了采用深度卷积网络进行特征提取和分类，将这三维的数据映射为二维图像。其中，横坐标为时间，纵坐标为频率，并通过jet颜色映射(HSV的一种变形，以蓝色开始，红色结束)来与幅值进行对应，最终呈现出的为彩色图像。某位患者的某个通道的某次发作前、发作时、发作后的部分时频图分别如图7a、图7b、图7c、图8a、图8b、图8c、图9a、图9b、图9c所示。同时，为了方便后一步输入到深度卷积神经网络中，将时频图的图像大小调整为28x28。

步骤4：利用时频图训练深度卷积神经网络

本发明采用LeNet-5的深度卷积神经网络模型，具有两个卷积层和两个池化层，卷积层滤波器大小均为5x5，池化层滤波器大小也均为5x5。经典结构中的第一个卷积层具有6个滤波器，而第二个卷积层具有16个滤波器。本发明综合考虑算法的运行时间以及特征提取的复杂程度，在第一个卷积层设置了20个滤波器，第二个卷积层设置了50个滤波器，来对图像进行特征提取。模型结构如图10所示。该模型具体的结构参数如下：

1)输入层：

将输入的28x28的图像变为灰度图像，并将灰度值映射到0～255。

2)卷积层C1：

采用20种5x5的卷积核，对输入图像进行步长为1的卷积，得到20张24x24的特征图。

卷积层的总运算公式如式(1)所示：

其中，l为层数，M_j为第j个特征图,i为对应特征图内的窗口索引，为对应的卷积核，为对应的偏置，表示为卷积后的特征图。假设输入图片大小为rxc，滤波器大小为axb，步长为i，则输出矩阵大小dxe，输出矩阵大小计算公式如式(2)所示：

3)池化层S1：

对20张24x24的特征图进行下采样处理，卷积核为2x2，步长为2。下采样处理后，每张特征图变为12x12。

下采样处理公式如式(3)所示：

其中，down()表示下采样函数，称为下采样的乘子偏差，称为对应的附加偏差。即分别对每个特征图中nxn大小的邻域进行加权求和运算或取最大值等运算，最后乘以一个乘子偏差，加上一个附加偏差，最终得到的特征图大小为上一层特征图的1/n，根据本发明选取的网络模型的参数，下采样后特征图大小变为原来的1/2,也就是分辨率降低了1/2。

4)卷积层C3：

采用50种5x5的卷积核，对上述20张12x12的特征图进行步长为1的卷积，得到50张8x8的特征图。

5)池化层S4：

对50张8x8的特征图进行下采样处理，卷积核为2x2，步长为2，下采样处理后，使每张特征图变为4x4。

6)IP1层：

对50张4x4的特征图进行全连接，输出500维的数据。

7)ReLU层：

采用ReLU函数对IP1的输出进行激活操作，其中，ReLU的函数表达式如式(4)所示：

8)IP2层：

将上述得到的500维数据进行全连接方式的降维，得到2维的数据，其中输出结果中最大分量对应的位置就是输出的分类结果。

9)Softmax层：

采用Softmax+交叉熵函数进行分类，因为网络最终结果为二分类，因此Softmax的计算过程如式(5)所示：

Softmax的结果相当于输入图像被分到每个标签的概率分布，该函数为单调增函数，即输入值越大，输出也越大，输入图像属于该标签的概率也越大。

同时，对Softmax的结果计算交叉熵分类损失函数如式(6)所示：

其中，k为真实标签值，N为一个批量的大小。

误差传播算法算法选用误差反向传播算法(BP算法)，具体的网络结构如图11所示，增加了loss层和accuracy层。其中，网络的训练参数如表(4)所示：

表4

以患者01为例，根据其前两次发作的数据用于训练网络，并用其第三次发作的数据来对网络进行测试，对于患者的每个通道(共有23个通道，每个通道代表一个导联)均进行此操作，最终得到的准确率如表5所示：

表5

步骤5：选出患者特定的最优的五个通道，并计算权重

根据患者所有通道测试数据的识别结果，选出准确率排在前五的通道，并由式(7)计算这五个通道对应的权重：

其中，w_i为第i个通道的权重，A_i为第i个通道的准确率，此步骤计算出的权重，用于最终的癫痫发作识别，这五个通道，对患者具有特异性，计算的权重也各不相同。

步骤6：采用加权求和的方法结合最优的五个通道进行癫痫的识别

在通过上述步骤计算出最优的五个通道的权重后，采用式(8)进行癫痫发作状态标签的计算：

其中，Label为最终计算出的标签值，取值为0或者1，1代表癫痫发作，0代表非癫痫发作，label_i为第i个通道分类得到的标签值，多个通道的加权和与0之间的距离如果小于0.5，说明加权和与0更近一些，也就是更倾向于非癫痫发作状态，因此将最终结果判定为非癫痫发作状态，否则判定为癫痫发作状态。