CN114081507B - 结合谱聚类和事件相关性检测的脑功能区定位方法及系统 - Google Patents

Abstract

本发明公开了结合谱聚类和事件相关性检测的脑功能区定位方法及系统，包括：采集被检测者在静息状态下及任务状态下的第一及第二皮质脑电信号；对所有第一皮质脑电信号进行小波包分解，并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个聚类类簇；对所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，计算所有导联的事件相关性及k个聚类类簇的平均事件相关性，得到类簇导联所在区域即判定为该任务相关功能的区域；使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位。本发明的算法具有快速、准确、泛化性好等优点。

Description

结合谱聚类和事件相关性检测的脑功能区定位方法及系统

技术领域

本发明涉及脑电分析领域，具体涉及结合谱聚类和事件相关性检测的脑功能区定位方法及系统。

背景技术

大脑是神经系统的最高级部分，大脑的病变将会影响人类的记忆、语言、运动等多个功能，严重影响人的生存和生活质量，其中，癫痫、脑胶质瘤是最常见的脑部病变。目前，对于脑胶质瘤和难治性癫痫的首要治疗手段是手术切除。而在神经外科手术中，最为关键的步骤是脑功能区的定位。

皮质脑电(ECoG)是一种从大脑皮层表面测量电位的神经活动记录技术。相比于皮质电刺激，ECoG是被动记录，可以减轻患者痛苦，节约时间，且较为安全。相比于术前影像方法，ECoG定位准确性更高、速度更快。“个体化”的脑功能区定位算法使用每个患者不同导联的ECS结果作为标签，提取导联的不同频率段能量占比作为特征矢量，训练分类模型，然后对其他导联进行分类。此方法的准确率可以80％以上，但训练得到的模型无法泛化到其他人。所以，神经外科亟需一种准确、快速、泛化性好的ECoG术中脑功能区定位算法。

因此，针对临床需求和现有脑功能区定位方法的局限性，提出一种结合谱聚类和事件相关性检测的脑功能区定位方法及系统，该定位方法具有快速、准确、泛化性好等优点。

发明内容

为了克服现有技术存在的缺点与不足，本发明提供一种结合谱聚类和事件相关性检测的脑功能区定位方法及系统。

本发明采用如下技术方案：

一种结合谱聚类和事件相关性检测的脑功能区定位方法，包括：

采集数据步骤：采集被检测者在静息状态下的第一皮质脑电信号，及被测者在任务状态下的第二皮质脑电信号；

功能区导联聚类步骤：对所有第一皮质脑电信号进行小波包分解，并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个聚类类簇即将所有导联分为k组；

导联事件相关性检测步骤：将所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，计算所有导联的事件相关性及k个聚类类簇的平均事件相关性，并比较大小，平均事件相关性最大的类簇导联所在区域即判定为该任务相关功能的区域；

使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位。

进一步，所述静息状态是被检测者处于清醒、安静、放松的状态。

进一步，所述任务状态是被检测者处于手部自主运动、物品命名任务或者数数任务的状态。

进一步，采集被检测者在静息状态下的第一皮质脑电信号，及被测者在任务状态下的第二皮质脑电信号，其采集具体是将硬膜下矩阵电极贴合在被检测者的大脑皮质上采集皮质脑电，电极的另一端连接脑电图仪，采集完毕后，从脑电图仪中导出数据。

进一步，功能区导联聚类步骤还包括预处理步骤，所述预处理步骤具体为：首先进行共同平均参考，然后用自适应工频干扰抵消算法去除50Hz附近工频及其谐波，最后去除基线漂移。

进一步，功能区导联聚类步骤中的小波包分解，具体是采用coif3小波基，对皮质脑电信号进行8层小波包分解。

进一步，导联事件相关性检测步骤，具体为：

对所有导联皮质脑电信号进行预处理，所述预处理为使用自适应工频干扰抵消算法去除50Hz附近工频及其谐波，然后去除基线漂移。

对预处理后的所有导联信号使用小波分解，对于不同的任务提取特定单子频带的能量；所述小波分解具体是选定coif3小波基，对皮质脑电信号进行7层小波分解。对于手部自主运动任务，提取第5层细节系数能量；对于物品命名任务，提取第2层细节系数能量；对于数数任务，提取第3层细节系数能量。

计算所有导联的事件相关性以及其所在类簇的平均事件相关性；

平均事件相关性计算如下：

ERD＝(A-R)/R*100，

其中，ERD为事件相关性，A为任务态皮质脑电特定频带的能量，R为静息态皮质脑电特定频带的能量，i为聚类类簇号，N为聚类类簇号i中导联总数，n为聚类类簇号i中导联号，则ERD_i[n]为聚类类簇号i中第n号导联的事件相关性，M_i为第i号类簇的平均事件相关性。

比较各类簇平均事件相关性，平均事件相关性最大的类簇中导联所在区域即判定为该任务相关功能的区域；

进一步，使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位，具体为：

以不同颜色表示导联的功能区属性，并在图中标注每一种颜色代表的功能区属性；以各导联的任务相关性为三维高斯分布的高，电极直径为高斯分布的方差，电极圆点为高斯分布的均值，得到混合高斯分布，在脑电地形图中画出混合高斯分布的等高线，则可得到功能区在大脑皮质上的概率分布图。

一种实现所述脑功能区定位方法的系统，包括：

信号采集模块，用于采集被检测者在静息状态下的第一皮质脑电信号，及被测者在任务状态下的第二皮质脑电信号；

功能区导联聚类模块，用于对所有第一皮质脑电信号进行小波包分解，并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个聚类类簇；

导联事件相关性检测模块，用于将所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，计算所有导联的事件相关性及k个聚类类簇的平均事件相关性，并比较大小，平均事件相关性最大的类簇导联所在区域即判定为该任务相关功能的区域；

功能区地形区绘制模块，用于使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位。

本发明的有益效果：

本发明提出的一种结合谱聚类和事件相关性检测的ECoG术中脑功能区定位方法，不需要皮质电刺激作为标签训练出模型，融合了谱聚类在个体间泛化性能较好、定位边界较准确，事件相关检测可以大致定出功能区属性的优点。将本方法应用于真实的临床案例，运动区、语言区、命名区的分类准确率分别达到100％、80％、77.78％，证明本方法可以快速有效地实现大脑多功能区定位。

附图说明

图1是本发明的工作流程图；

图2是本发明的功能区导联聚类模块的工作流程图；

图3是本发明导联事件相关性检测模块的工作流程图；

图4是本发明实施例的功能区地形绘制图。

具体实施方式

下面结合实施例及附图，对本发明作进一步地详细说明，但本发明的实施方式不限于此。

实施例

如图1-图4所示，一种结合谱聚类和事件相关性检测的脑功能区定位方法,包括如下步骤：

S1信号采集：

使用硬膜下矩阵电极和脑电图仪采集3名脑肿瘤患者的静息态和任务态术中ECoG，出于对隐私的保护，本发明中使用A、B、C表示这三名患者。三名患者的数据详情如下表1所示，各导联所在功能区属性由皮质电刺激得到。

表1实验样本

S2采用功能区聚类算法：如图2所示，对所有导联静息状态下的皮质脑电信号进行小波包分解并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个类簇。

具体包括如下：

S2.1数据预处理

首先对患者的所有导联进行共同平均参考，然后用自适应工频干扰抵消算法去除50Hz附近工频及其谐波，最后去除基线漂移。

S2.2小波包分解

选定coif3小波基，对对预处理后各导联的数据进行8层小波包分解，然后重构β节律(14-32Hz)信号。

S2.3构建邻接矩阵

计算两两重构的β节律信号之间的似然同步性，构成邻接矩阵；

似然同步性的具体计算公式如下：

X_k,i＝(x_k,i,x_k,i+l,x_k,i+2l,...,x_k,i+(m-1)l) {1}

式{1}中X_k,i为内嵌向量，其中k表示信道数量(k＝1,2,…,M)，i表示信号记录的时间点(i＝1,2,…N)；式{2}中两个内嵌向量X_k,i和X_k,j之间的欧氏距离小于ε的概率为

通常研究中给定/>

概率为0.05；假定M信道信号中与内嵌向量X_k,i之间距离小于ε_k,i，并满足窗口条件(w₁＜|i-j|＜w₂)的内嵌向量X_k,j有个H_i,j，即有{3}；对于k道信号和另外M-1道信号i和j时刻的似然同步值S_k,i,j定义如{4}；最后对所有时间点j求平均，得到k道信号和其他M-1道信号在i时刻的似然同步值为{5}。

S2.4谱聚类算法

对邻接矩阵进行k类谱聚类，得到k个聚类类簇；k为功能区类别数。

S3采用事件相关检测算法：对于不同的任务，计算所有类簇的平均事件相关性，比较各类簇平均事件相关性的大小，平均事件相关性最大的类簇中导联所在区域即判定为该任务相关功能的区域；具体包括如下步骤：

S3.1数据预处理

对患者的所有导联用自适应工频干扰抵消算法去除50Hz附近工频及其谐波，然后去除基线漂移。

S3.2小波分解和提取特征

选定coif3小波基，对预处理后干净的皮质脑电信号进行8层小波分解。对于手部自主运动任务，分别提取任务态和静息态第5层细节系数能量；对于物品命名任务，分别提取任务态和静息态第2层细节系数能量；对于数数任务，分别提取任务态和静息态第3层细节系数能量。

S3.3计算各导联事件相关性

事件相关性计算如{6}：

ERD＝(A-R)/R*100 {6}

其中，ERD为事件相关性，A为任务态皮质脑电特定频带的能量，R为静息态皮质脑电特定频带的能量。

S3.4比较各类簇平均事件相关性

平均事件相关性计算如{7}：

其中，i为聚类类簇号，N为聚类类簇号i中导联总数，n为聚类类簇号i中导联号，则ERD_i[n]为聚类类簇号i中第n号导联的事件相关性，M_i为第i号类簇的平均事件相关性。

平均事件相关性最大的类簇中导联所在的大脑皮质区域即判定为该任务相关功能的区域。

S3.5重复S3.3～S3.4直至分析完所有任务。

S4算法相关性评价

由表2可以看出本方法分类得到的各导联的功能区属性与金标准皮质电刺激得出的基本相同，运动区分类准确率达到100％，语言区分类准确率达到80％，命名区分类准确率达到77.78％，证明本方法可以快速有效地实现大脑多功能区定位，有望在临床中得到广泛应用。

表2本方法结果与皮质电刺激对比

S5绘制功能地形图，效果图如图4所示。

以不同颜色表示导联的功能区属性，并在图中标注每一种颜色代表的功能区属性；以各导联的任务相关性为三维高斯分布的高，电极直径为高斯分布的方差，电极圆点为高斯分布的均值，得到混合高斯分布，在图中画出混合高斯分布的等高线，则可得到功能区在大脑皮质上的概率分布图。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受所述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (7)

1.一种结合谱聚类和事件相关性检测的脑功能区定位方法，其特征在于，包括：

采集被检测者在静息状态下的第一皮质脑电信号，及被测者在任务状态下的第二皮质脑电信号；

所述任务状态包括数数任务、手部自主运动任务及物品命名任务；

对所有第一皮质脑电信号进行小波包分解，并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个聚类类簇；

对所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，计算所有导联的事件相关性及k个聚类类簇的平均事件相关性，并比较k个聚类类簇的平均事件相关性的大小，平均事件相关性最大的类簇导联所在区域即判定为该任务相关功能的区域；

将所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，具体为：

选定coif3小波基，对预处理后的第一及第二皮质脑电信号进行八层小波分解，具体是：

对于手部自主运动任务，分别提取任务态和静息态第五层细节系数能量；

对于物品命名任务，分别提取任务态和静息态第二层细节系数能量；

对于数数任务，分别提取任务态和静息态第三层细节系数能量；

所述计算所有导联事件相关性，具体为：

ERD=(A-R)／R*100

其中，ERD为事件相关性，A为任务态皮质脑电特定频带的能量，R为静息态皮质脑电特定频带的能量；

使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位。

2.根据权利要求1所述的脑功能区定位方法，其特征在于，通过设置在硬膜下的矩阵电极和脑电图仪采集被测者的第一及第二皮质脑电信号。

3.根据权利要求1所述的脑功能区定位方法，其特征在于，所述静息状态包括被测者清醒、安静及放松的状态。

4.根据权利要求1所述的脑功能区定位方法，其特征在于，对所有第一皮质脑电信号进行小波包分解具体是对预处理后的各导联数据进行八层小波包分解，然后重构14 Hz-32Hz的β节律信号。

5.根据权利要求1所述的脑功能区定位方法，其特征在于，在所述对所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量之前还包括数据预处理步骤，所述数据预处理包括自适应工频干扰抵消算法去除50Hz附近工频及其谐波，然后去除基线漂移。

6.根据权利要求1-5任一项所述的脑功能区定位方法，其特征在于，所述k个聚类类簇的平均事件相关性，具体计算为：

其中，i为聚类类簇号，N为聚类类簇号i中导联总数，n为聚类类簇号i中导联号，则ERD_i［n］为聚类类簇号i中第n号导联的事件相关性，Mi为第i号类簇的平均事件相关性。

7.一种实现权利要求1-6任一项所述的脑功能区定位方法的系统，其特征在于，包括：

信号采集模块，用于采集被检测者在静息状态下的第一皮质脑电信号，及被测者在任务状态下的第二皮质脑电信号；

功能区导联聚类模块，用于对所有第一皮质脑电信号进行小波包分解，并重构β节律信号，计算导联间β节律信号的似然同步性构成邻接矩阵，采用谱聚类算法对此邻接矩阵进行聚类分析，得到k个聚类类簇；

导联事件相关性检测模块，用于将所有第一皮质脑电信号及第二皮质脑电信号进行小波分解并提取相应细节系数能量，计算所有导联的事件相关性及k个聚类类簇的平均事件相关性，并比较大小，平均事件相关性最大的类簇导联所在区域即判定为该任务相关功能的区域；

功能区地形区绘制模块，用于使用混合高斯函数画出脑电地形图，并标注所有导联的功能区属性，完成功能定位。

Images