CN107887031B - 基于同步筛选的脑区间耦合分析方法 - Google Patents

Abstract

本发明公开了一种基于同步筛选的脑区间耦合分析方法，本发明首先同步采集了不同握力输出时，32通道脑电信号EEG和12通道肌电信号EMG，为研究握力输出过程中，EEG运动区与感觉区的耦合关系，对大脑运动区及运动感觉区的导联C3，C4，CP5和CP6的脑电信号和肱桡肌信号进行分析；然后利用同步筛选算法将EEG与EMG的同步信息提取出来，得到EEG与EMG相关的数据；最后计算SSEM的符号传递熵，以确定其耦合关系。本发明可去除与运动的非相关数据，缩小数据规模。

Description

基于同步筛选的脑区间耦合分析方法

技术领域

本发明属于信号处理领域，涉及一种脑电、肌电信号的耦合性方法，特别涉及一种应用于脑区间状态分析方法。

背景技术

运动皮层的脑电信号(electroencephalogram，EEG)和肌电信号(electromyogra-phy,EMG)分别反映运动控制信息和肌肉对大脑控制意图的功能响应信息，脑肌电之间具有同步性。不同脑区拥有其特殊的功能，但是大脑作为一个有机的整体，各个脑区之间必然存在某种联系。大脑控制上肢运动时，运动神经系统通过神经的同步振荡来传递运动控制信息，相关的脑区和肌肉组织相关的周围神经之间也自动地发生了同步化，从而激发肢体的运动。大脑皮层与肌肉之间的同步振荡不仅从发出命令的大脑皮层传导至作为效应器的肌肉，同时也由肌肉反馈回大脑皮层，这些反馈信息必然会被大脑作为发出下一次控制信息的重要依据。近来，探究EEG与肢体运动之间的关系已成为人们研究的热点问题。Ponten等人利用同步似然性构建了癫痫患者脑电信号在癫痫发作不同阶段的脑功能网络，并通过对网络属性的分析发现癫痫在发病中和发病后，脑功能网络的特征路径长度会增大；Scott曾以奥托博克肌电假手为例，设计了电刺激感觉反馈系统，通过安装在假手食指上的应变式力矩传感器测量捏力大小，通过调制电刺激脉冲频率的方式将刺激频率与捏力的大小成比例对应，以此将捏力感觉反馈给使用者；美国国防先进研究项目局(DARPA)采用目标肌肉神经移植术重新连接被截肢位置和肌肉组织的神经信号，并作为一个“放大器”来控制假肢。同时，利用名为FINE(Flat Interface Nerve Electrode)技术直接向使用者的大脑反馈假手的感觉。FINE使神经元变得平缓，让部分神经元能够接触电流并提供反馈，恢复患者的触觉。然而在处理与运动相关的EEG时都面临的一个问题是，EEG中包含有大量与运动非相关的成分，这些非相关成分会使EEG数据规模变得庞大，湮没EEG中与运动相关的成分，这会使求解变得困难，同时增大计算复杂度。针对这一问题，本发明提出了一种基于同步筛选-符号传递熵的脑区间耦合分析方法，可以滤除EEG中与运动的非相关成分，去除非相关数据、降低运算复杂度的同时，保留了数据有效信息，可以快速准确的分析脑区间的耦合关系。

发明内容

为了能够快速准确的分析上肢握力输出时，大脑运动区与感觉区的耦合关系，本发明提出了一种基于同步筛选-符号传递熵的脑区间耦合分析方法。首先同步采集了不同握力输出时，32通道脑电信号(EEG)和12通道肌电信号(EMG)，为研究握力输出过程中，EEG运动区与感觉区的耦合关系，对大脑运动区及运动感觉区的导联C3，C4，CP5，CP6等多个通道的脑电信号和肱桡肌信号进行分析；然后利用同步筛选算法将EEG与EMG的同步信息提取出来，得到EEG与EMG相关的数据(Synchronous Screening of EEG Signals Based on EMGSignals,SSEM)；最后计算SSEM的符号传递熵，以确定其耦合关系。

本发明方法主要包括以下步骤：

(1).同步采集了在不同握力输出时，32通道脑电信号和12通道肌电信号。

(2).选取步骤(1)中获取的C3、C4、CP5、CP6通道脑电信号和肱桡肌肌电信号，对其进行小波阈值消噪处理，滤除心电、体位的影响并对EEG与EMG进行了时延分析。

(3).使用同步筛选算法，以EMG肱桡肌信号为基准数据，对EEG运动区C3、C4通道与感觉区CP5、CP6通道的信号进行同步筛选，得到EEG与EMG相关的数据SSEM。

同步筛选算法，具体算法如下：

对于任意两列时间序列：X₁＝{x₁,x₂,x₃,…}；Y₁＝{y1,y2,y3,…}；X₁、Y₁每一元素都有一个固定的标签，各元素在X₁或Y₁中的索引值为其对应的标签值。X₁,Y₁对应的标签序列为WX₁、WY₁,WX₁＝{1,2,3,…}；WY₁＝{1,2,3,…}。

X₂＝V(X₁) (1)

V(X)表示将X进行非递减的排序。在改变X₁中元素的位置得到非递减的X₂时，为了实现标签和元素的一一对应，标签序列中标签会发生与X₁中元素相同的位置改变，从而得到X₂对应的索引序列WX₂。

Y₂＝T(WX₂,Y₁) (2)

T(WX₂,Y₁)表示将WX₂中存放的标签序列当作Y₁新的索引序列，对Y₁进行重新排列。因为Y₁为非递减序列，若X₁、Y₁完全同步,Y₂必为非递减序列。当X₁、Y₁不完全同步时,Y₁中的非同步成分会导致无法将T(WX₂,Y₁)映射成非递减序列；根据这一特点，我们可通过Y₂去掉Y₁中的非同步成分；

若Y₂中某索引大的元素值等于索引小的元素值，则将该值赋值为0，从而得到新的序列Y₃。Y₃即：对于Y₂中的y_n、y_m,此处n、m为索引，若y_n＝y_m，n<m，则将y_m置0。然后将Y₃按照WY₁中的标签进行还原排序：

Y₄＝T(WY₁,Y₃) (4)

T(WY₁,Y₃)表示将WY₁中存放的标签序列当作Y₃新的索引，对Y₃进行重新排列。Y₄与原始时间序列Y₁具有相同的排序，但与X₁非同步的成分都已置0。将X₁中为0的元素全部去除，得到Y₅。Y₅即为Y₁基于X₁同步筛选后的结果。

(4).结合SSEM，计算在不同握力下大脑运动区C3、C4与感觉区CP5、CP6之间基于同步筛选的符号传递熵，分析脑区间的耦合性。

所述基于同步筛选的符号传递熵，具体算法如下：

同步筛选即为步骤(3)中提到的同步筛选算法。符号传递熵算法如下：

首先对时间序列进行符号化处理，根据时间序列的数值特点进行符号划分的静态法,公式如下：

式(5)中，S_i为划分好的符号集，x_i为时间序列X_t＝{x₁,x₂,x₃,…}中的数据点；B₁,B₂...B_N-1是一系列截断点。设符号划分数量为pieces，pieces＝N。pieces的大小可以任意的放大或缩小；将原始数据符号化后，计算的传递熵称为符号传递熵。

传递熵定义如式(6)：设X_t,Y_t为两个长度为n的时间序列，各自从t-1到t-p的历史记为

表示X与Y^-在X^-条件下的互信息。则Y_t对X_t的传递熵T_X→Y定义为X与Y^-在X^-下的互信息。

在相同的符号集下，计算T_X→Y越大，表明数据X传递到数据Y的信息量越多，X对Y的影响越大。

本发明与已有的诸多与运动相关的脑区间EEG耦合分析方法相比，具有如下特点：

探究动力系统之间的耦合关系时，传统的方法是先将时间序列符号化，再探求其耦合关系是一种重要的方法。对于两列含有大量非相关数据的时间序列中，非相关数据会使数据规模变大，会对耦合分析产生不利影响，使用传统的方法分析，计算复杂度往往很大，且难以寻找耦合关系。针对这一问题，本发明提出了基于同步筛选-符号传递熵的脑区间耦合分析方法，在探究耦合关系之前，先对数据进行同步筛选处理。同步筛选算法可去除与运动的非相关数据，缩小数据规模。这会降低运算复杂度，去除非相关数据对耦合分析的影响。

附图说明

图1为本发明的实施流程图。

图2为非同步成分示意图，图中标记的两个黑点为Y相对于X的非同步点。

图3-(a)、图3-(b)为左手发力，在不同握力输出下大脑运动区与感觉区间的符号传递熵的分析；

图3-(c)、图3-(d)为右手发力，在不同握力输出下大脑运动区与感觉区间的符号传递熵的分析。

具体实施方式

下面结合附图对本发明的实施例作详细说明：本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程。

其中图2为非同步成分示意图；

如图1所示，本实施例包括如下步骤：

步骤一，同步采集了不同握力输出时，32通道脑电信号(EEG)和12通道肌电信号(EMG)。具体的采集操作如下：

选取6名无任何病史的健康受试者，在安静的实验室内，受试者端坐于木制座椅上，将握力计放置于受试者手心。受试者在接收到信号之后，立即将握力计从零刻度握至目标刻度(握力计的目标刻度为5kg，10kg，20kg)，受试者持续输出握力，使握力计在目标刻度处保持5秒，受试者松开握力计，休息5秒，进行下一组实验。每名受试者的同一只手，在相同握力条件下进行10组实验。数据采集所用的是128导BrainAmp DC脑电采集系统(BrainProducts GmbH,Germany)。同步采集了上述实验范式下32通道脑电和12通道肌电信号。BrainAmp DC设备也用于采集左右手臂指浅屈肌(flexor digitorum superficialis,FDS)、肱桡肌(Brachioradialis muscle,BM)、桡侧腕屈肌(Radial wrist flexor,RWF)、尺侧腕屈肌(Ulnar wrist flexor,UWF)、肱二头肌(Musculus biceps brachii,MBB)和肱三头肌(Triceps)的肌电信号。在电极安放前，头皮已经清洗干净，皮肤表面也经过酒精清理，采样频率为1000Hz。

步骤二，选取步骤一中获取的C3,C4,CP5,CP6通道脑电信号和肱桡肌肌电信号，对其进行小波阈值消噪处理，滤除心电，体位，以及其它噪声的影响并针对EEG与EMG之间的shiyan2进行延时分析。

EEG->EMG,EMG->EEG存在一定的时延。研究发现EEG与EMG之间的信息流延迟时间为20-30ms。通过实验分析得出6名受试者的平均时延为：EEG->EMG的时延是21ms，EMG->EEG的时延是30ms。

步骤三，使用同步筛选算法，以EMG肱桡肌信号为基准数据，对EEG运动区(C3,C4通道)与感觉区(CP5,CP6通道)信号进行同步筛选，得到SSEM。

对于不同握力下大脑运动区(C3、C4通道)感觉区(CP5、CP6通道)的脑电信号，以左右手肱桡肌的数据为参照进行同步筛选，从而得到四组SSEM数据。EEG的数据点的个数均为4970个,由于SSEM较EEG去除了大量与运动无关的成分，SSEM的数据规模必定显著小于EEG。表1为左右手不同握力下,六组数据的SSEM数据点个数的均值与EEG数据点个数的比较。可以看出：同步筛选后得到的SSEM数据点大大减少，且随着握力的增加SSEM的数据点个数增多，这表明随着握力的增加，EEG与EMG的同步点增多，同步性增强，信息交互增强。

表1不同握力下原始EEG与SSEM数据点个数

记每个符号平均分得的数据点个数为pa，pa为数据点个数和pieces的比值,pa值取整数。pa值是对数据划分粗细的指标。显然pa值越小，平均每个符号代表的数据点个数越少，数据划分的越细。由表2可知，相同的pieces值下，SSEM的pa值远远小于原始EEG的。当计算与运动相关的脑区间符号化传递熵时，SSEM较EEG的数据点个数大大减小。符号化时，考虑到平均每个符号可代表的数据点个数。对于小规模的数据，用一个较小的pieces值就可实现大规模数据在较大pieces值下相同的划分粗细程度。数据的规模，符号集大小正是影响符号传递熵计算复杂度的关键所在。

表2输出握力10kg时原始EEG与SSEM在不同Pieces下pa值的比较

步骤四，结合SSEM计算在不同握力下运动区(C3,C4)与感觉区(CP5,CP6)之间基于同步筛选的符号传递熵。

如图3-(a)、图3-(b)、图3-(c)、图3-(d)所示，随着握力输出的增大，符号传递熵值也随之增大，即耦合关系增强；这是由于随着握力的增大，大脑运动区需要给肌肉神经传递更多的信息，同时肌肉也需要反馈更多的信息到大脑感觉区，这导致运动区与感觉区的信息流增大，耦合强度增强。感觉区传递到运动区的信息要多于运动区传递到感觉区的信息，并且这种差值随着握力的增大有增大的趋势；运动控制信号由大脑运动区发出，传递到肢体的肌肉神经，从而控制肢体运动，运动控制信号并不会由运动区传递到感觉区，因此我们测得的运动区到反馈区的传递熵并不是运动控制信息的信息量，而是一些其它信息，不是运动区激发的主要信息；感觉区的感觉信息作为整个运动系统的反馈信息传递到运动区，因此我们测得的感觉区到运动区的传递熵是运动系统的反馈信息量，正是感觉区激发出的主要信息。因此感觉区传递到运动区的信息要多于运动区传递到感觉区的信息。

Claims (1)

1.基于同步筛选的脑区间耦合分析方法，其特征在于，该方法包括如下步骤：

(1).同步采集了在不同握力输出时，32通道脑电信号EEG和12通道肌电信号EMG；

(2).选取步骤(1)中获取的C3、C4、CP5、CP6通道脑电信号和肱桡肌肌电信号，对其进行小波阈值消噪处理，滤除心电、体位的影响并对EEG与EMG进行了时延分析；

(3).使用同步筛选算法，以EMG肱桡肌信号为基准数据，对EEG运动区C3、C4通道与感觉区CP5、CP6通道的信号进行同步筛选，得到EEG与EMG相关的数据SSEM；

同步筛选算法，具体算法如下：

对于任意两列时间序列：X₁＝{x₁,x₂,x₃,…}；Y₁＝{y₁,y₂,y₃,…}；X₁、Y₁每一元素都有一个固定的标签，各元素在X₁或Y₁中的索引值为其对应的标签值；X₁,Y₁对应的标签序列为WX₁、WY₁,WX₁＝{1,2,3,…}；WY₁＝{1,2,3,…}；

X₂＝V(X₁) (1)

V(X)表示将X进行非递减的排序；在改变X₁中元素的位置得到非递减的X₂时，为了实现标签和元素的一一对应，标签序列中标签会发生与X₁中元素相同的位置改变，从而得到X₂对应的索引序列WX₂；

Y₂＝T(WX₂,Y₁) (2)

T(WX₂,Y₁)表示将WX₂中存放的标签序列当作Y₁新的索引序列，对Y₁进行重新排列；因为Y₁为非递减序列，若X₁、Y₁完全同步,Y₂必为非递减序列；当X₁、Y₁不完全同步时,Y₁中的非同步成分会导致无法将T(WX₂,Y₁)映射成非递减序列；根据这一特点，通过Y₂去掉Y₁中的非同步成分；

其中Y₃(i)表示Y₃中第i个索引元素值，Y₂(i)表示Y₂中第i个索引元素值；若Y₂中某索引大的元素值等于索引小的元素值，则将该值赋值为0，从而得到新的序列Y₃；Y₃即：对于Y₂中的y_n、y_m,此处n、m为索引，若y_n＝y_m，n<m，则将y_m置0；然后将Y₃按照WY₁中的标签进行还原排序：

Y₄＝T(WY₁,Y₃) (4)

T(WY₁,Y₃)表示将WY₁中存放的标签序列当作Y₃新的索引，对Y₃进行重新排列；Y₄与原始时间序列Y₁具有相同的排序，但与X₁非同步的成分都已置0；将X₁中为0的元素全部去除，得到Y₅；Y₅即为Y₁基于X₁同步筛选后的结果；

(4).结合SSEM，计算在不同握力下大脑运动区C3、C4与感觉区CP5、CP6之间基于同步筛选的符号传递熵，分析脑区间的耦合性；

所述基于同步筛选的符号传递熵，具体算法如下：

同步筛选即为步骤(3)中提到的同步筛选算法；符号传递熵算法如下：

首先对时间序列进行符号化处理，根据时间序列的数值特点进行符号划分的静态法,公式如下：

式(5)中，S_i为划分好的符号集，x_i为时间序列X_t＝{x₁,x₂,x₃,…}中的数据点；B₁,B₂...B_N-1是一系列截断点；设符号划分数量为pieces，pieces＝N；pieces的大小可以任意的放大或缩小；将原始数据符号化后，计算的传递熵称为符号传递熵；

传递熵定义如式(6)：设X_t,Y_t为两个长度为n的时间序列，各自从t-1到t-p的历史记为

I(X；Y^-|X^-)表示X与Y^-在X^-条件下的互信息；则Y_t对X_t的传递熵T_X→Y定义为X与Y^-在X^-下的互信息；

在相同的符号集下，计算T_X→Y越大，表明数据X传递到数据Y的信息量越多，X对Y的影响越大。

Images