基于蝶腭神经节刺激的脑电信号特征提取方法

本发明属于生物医学信息处理领域，尤其涉及到一种基于蝶腭神经节刺激的脑电信号特征提取方法。

背景技术：

蝶腭神经节是人体最大的颅外副交感神经节，其分支广泛分布在眶部、口腔及鼻腔。研究表明，对蝶腭神经节的神经刺激能够改善原发性头痛以及鼻炎病情。传统的蝶腭神经节刺激技术常用电刺激以及针刺刺激等有创刺激方法，存在肌肉痉挛，有创部位感染等副作用。蝶腭神经节刺激技术对非过敏性鼻炎、偏头痛等部分炎症类疾病有良好的疗效；采用静息态功能磁共振成像对该神经刺激技术进行研究，发现刺激后大脑边缘系统和初级感觉系统的内在功能活性显著提升。脑电活动是大脑皮层自发产生的节律性电位变化。用脑电图仪在头皮表面记录到的脑电活动，称为脑电图。研究证实，位于头皮的脑电图可记录来自深部脑区的电信号。

然而，上述现有技术还存在一些问题：1、现有技术的蝶腭神经节刺激技术存在肌肉痉挛，有创部位感染等副作用；2、缺乏有效的特征提取方式，难以评估刺激的可靠性。

技术实现要素：

有鉴于此，为了克服现有技术神经刺激带来的不良影响等问题，提供了一种采用超声波刺激方式的基于蝶腭神经节刺激的脑电信号特征提取方法。

本发明的技术解决方案是，提供了一种基于蝶腭神经节刺激的脑电信号特征提取方法，包括以下步骤：

步骤1：利用超声波通过两侧鼻翼实现对蝶腭神经的刺激；

步骤2：在头皮上放置脑电采集电极，利用脑电采集装置记录不同蝶腭神经节刺激阶段的脑电信号；

步骤3：进行脑电信号预处理操作，对经预处理后的脑电信号进行脑电相对功率特征提取和脑电有效连接特征提取；

步骤4：根据所述脑电相对功率特征和脑电有效连接特征评估超声波对蝶腭神经的刺激效果。

可选的，所述脑电信号预处理操作，包括：对脑电信号进行带通滤波，通带频率为0.5-80hz；对脑电信号进行带阻滤波，陷波频率为50hz；脑电信号参考电极重建为平均参考；采用伪影子空间重建技术处理脑电信号，去除不良通道，截除不良数据；采用独立成分分析技术处理脑电信号，去除眼动伪迹、肌电伪迹；将脑电信号按照频段分为5个信号节律，其中包括delta节律、theta节律、alpha节律、beta节律和gamma节律。

可选的，步骤4中所述脑电相对功率特征提取具体为：

采用welch功率谱估计算法，计算得出脑电信号的功率谱密度值；

采用以下公式计算脑电信号的相对功率特征：

其中，welch(fs，fe)表示频率从fs到fe之间功率密度值的总和。

可选的，所述welch功率谱估计算法采用4秒时长的汉宁窗函数，重叠率为50％。

可选的，步骤4中所述脑电有效连接特征提取具体为：

建立滑动窗口的多变量自回归模型，公式为：

其中，x(t)＝[x1(t)，x2(t)，…，xi(t)]^t为电极数量为i的脑电信号时间序列，a(n)为i×i系数矩阵，e(t)为脑电信号的白噪声余项，p为多变量自回归模型阶数；

采用频域转换算法将所述滑动窗口的多变量自回归模型转换至频域，公式为：

其中，h(f)为传递矩阵，f为频率；

采用直接定向传递函数计算脑电信号的有效连接特征，公式为：

其中，s(f)＝x(f)x(f)^*，k表示电极的总数，i、j表示电极序号。

可选的，所述多变量自回归模型阶数由汉南奎因信息准则算法确定。

可选的，步骤1的具体实现方法为：

在鼻翼两侧处涂抹超声耦合剂；将超声换能器分别贴敷于两侧鼻翼；启动超声换能器施加设定超声波频率、设定脉冲重复频率、设定脉冲长度、设定刺激时长、设定强度的超声刺激；

设定超声波频率为150khz至700khz；所述步骤3中所述设定脉冲重复频率为20hz至2000hz；所述设定脉冲长度为0.1ms至100ms；所述设定刺激时长为1分钟至15分钟；所述设定强度为0.1w/cm2至15w/cm2。

可选的，步骤2中，所述不同蝶腭神经节刺激阶段的脑电信号包括：蝶腭神经节刺激开始前30分钟记录到的脑电信号；蝶腭神经节刺激同时记录到的脑电信号；蝶腭神经节刺激结束后30分钟记录到的脑电信号。

本发明与现有技术相比，具有如下优点：本发明采用无创方法，避免了传统电刺激及植入式电刺激法所产生的损伤和副作用；本发明所提供的蝶腭神经刺激方法，操作较传统方法来说操作更为简单方便，成本更低，不需要很强的专业背景；本发明所提供的蝶腭神经刺激方法，可以与功能磁共振成像、脑电检测和肌电检测兼容；本发明所提供的蝶腭神经刺激方法是首次使用超声换能技术对蝶腭神经节的鼻部神经分布进行神经刺激，并通过脑电信号的特征提取方法是首次采用脑电信号记录和分析方法探究基于鼻部的超声刺激对于神经的调控效应。

附图说明

图1为本发明的原理示意图。

图2为10-10国际标准脑电导联系统图。

图3为蝶腭神经节刺激前后脑电相对功率在不同频段不同脑区的均值方差对照图。

图4为部分脑电导联(f3,f4,c3,c4,t7,t8,o1,o2)脑电有效连接在不同频段的色块图。

具体实施方式

以下结合附图对本发明的优选实施例进行详细描述，但本发明并不仅仅限于这些实施例。本发明涵盖任何在本发明的精神和范围上做的替代、修改、等效方法以及方案。

为了使公众对本发明有彻底的了解，在以下本发明优选实施例中详细说明了具体的细节，而对本领域技术人员来说没有这些细节的描述也可以完全理解本发明。

如图1所示，本发明的一种蝶腭神经节刺激方法，在实施蝶腭神经节刺激前，先将在被测试者鼻翼位置涂抹超声耦合剂，后使用超声换能器探头贴附鼻翼对被测试者进行超声刺激。本发明提供了一种基于蝶腭神经节刺激的脑电信号特征提取方法，包括以下步骤：

步骤1：利用超声波通过两侧鼻翼实现对蝶腭神经的刺激；

步骤2：在头皮上放置脑电采集电极，利用脑电采集装置记录不同蝶腭神经节刺激阶段的脑电信号；

步骤3：进行脑电信号预处理操作，对经预处理后的脑电信号进行脑电相对功率特征提取和脑电有效连接特征提取；

步骤4：根据所述脑电相对功率特征和脑电有效连接特征评估超声波对蝶腭神经的刺激效果。

所述脑电信号预处理操作，包括：对脑电信号进行带通滤波，通带频率为0.5-80hz；对脑电信号进行带阻滤波，陷波频率为50hz；脑电信号参考电极重建为平均参考；采用伪影子空间重建技术处理脑电信号，去除不良通道，截除不良数据；采用独立成分分析技术处理脑电信号，去除眼动伪迹、肌电伪迹；将脑电信号按照频段分为5个信号节律，其中包括delta节律、theta节律、alpha节律、beta节律和gamma节律。

步骤4中所述脑电相对功率特征提取具体为：

采用welch功率谱估计算法，计算得出脑电信号的功率谱密度值；

采用以下公式计算脑电信号的相对功率特征：

其中，welch(fs，fe)表示频率从fs到fe之间功率密度值的总和。

所述welch功率谱估计算法采用4秒时长的汉宁窗函数，重叠率为50％。

步骤4中所述脑电有效连接特征提取具体为：

建立滑动窗口的多变量自回归模型，公式为：

其中，x(t)＝[x1(t)，x2(t)，…，xi(t)]^t为电极数量为i的脑电信号时间序列，a(n)为i×i系数矩阵，e(t)为脑电信号的白噪声余项，p为多变量自回归模型阶数；

采用频域转换算法将所述滑动窗口的多变量自回归模型转换至频域，公式为：

其中，h(f)为传递矩阵，f为频率；

采用直接定向传递函数计算脑电信号的有效连接特征，公式为：

其中，s(f)＝x(f)x(f)^*，k表示电极的总数，i、j表示电极序号。

所述多变量自回归模型阶数由汉南奎因信息准则算法确定。

步骤1的具体实现方法为：

在鼻翼两侧处涂抹超声耦合剂；将超声换能器分别贴敷于两侧鼻翼；启动超声换能器施加设定超声波频率、设定脉冲重复频率、设定脉冲长度、设定刺激时长、设定强度的超声刺激；

设定超声波频率为150khz至700khz；所述步骤3中所述设定脉冲重复频率为20hz至2000hz；所述设定脉冲长度为0.1ms至100ms；所述设定刺激时长为1分钟至15分钟；所述设定强度为0.1w/cm2至15w/cm2。

步骤2中，所述不同蝶腭神经节刺激阶段的脑电信号包括：蝶腭神经节刺激开始前30分钟记录到的脑电信号；蝶腭神经节刺激同时记录到的脑电信号；蝶腭神经节刺激结束后30分钟记录到的脑电信号。

脑电信号采集电极分布遵循10-10国际标准导联系统，设置采样率为500赫兹，对大脑额叶、顶叶、颞叶、枕叶四个脑区进行脑电信号记录，具体电极位置如图2所示。其中设置cz为参考电极，fpz为地电极。

脑电信号记录在封闭的隔音室内进行，室内温度控制在24摄氏度。被测试者需清洁头发，在电极位置注入导电膏，确保记录过程中电极阻抗低于5k欧姆。

脑电信号记录被测试者需尽量避免面部及肢体活动，保持清醒放松。记录分为2个阶段：(1)刺激开始0.5小时前，记录5分钟；(2)刺激结束0.5小时后，记录5分钟。

对脑电信号样本进行预处理。步骤如下：(1)带通滤波0.5～45hz，去除基线漂移以及50hz工频干扰；(2)平均参考；(3)使用伪影子空间重建(asr)，去除不良通道和数据段；(4)使用独立成分分析(ica)，分解去除眼动以及肌电等噪声成分。(5)将脑电信号分为4个频段：delta(1～4hz)，theta(4～8hz)，alpha(8～13hz)，beta(13～30hz)。

采用welch法求脑电信号的功率谱密度，选用4秒时长的hanning窗，50％重叠。

各个频段相对功率计算公式为

式中welch(fs，fe)表示频段从fs到fe之间功率密度值的总和，计算结果如图3所示。

采用基于格兰杰因果关系(gc)的直接定向传递函数(ddtf)对脑电信号在不同频段的有效连接值进行计算，首先建立滑动窗口的多变量自回归模型(amvar)。

amvar定义为：

式中x(t)＝[x1(t)，x2(t)，…，xi(t)]^t为电极数量为i的脑电信号时间序列，

a(n)为i×i系数矩阵，e(t)为白噪声余项，p为由汉南奎因信息准则(hq)确定的多变量自回归模型阶数。

使用快速傅立叶变换将多变量自回归模型转化至频域，

x(f)＝h(f)e(f)

其中为系统的传递矩阵，f为频率。

ddtf定义为：

其中，

式中k表示通道的数量。

pij(f)定义为：

其中，

s(f)＝x(f)x(f)^*。

最后得到各频段ddtf矩阵如图4所示。

结果表明，蝶腭神经节刺激后alpha频段脑电的相对功率和有效连接以及beta频段的相对功率都显著增强。相关fmri研究证实，对鼻部施加振动刺激可增强边缘系统以及初级感觉系统的内在功能活性，并且边缘系统活性与eeg信号alpha频段的功率以及有效连接性增强呈正相关，beta频段的功率增强则源于初级体感皮层以及运动皮层的激活。本发明应用蝶腭神经节刺激以及eeg信号采集分析的方法，从脑电角度印证蝶腭神经节刺激对于边缘系统以及初级感觉网络的内在功能活性提升作用。本发明的应用可为蝶腭神经节刺激作用的评估，以及蝶腭神经节刺激治疗相关疾病的疗效评估提供更为简便、低成本且有效的方案。以上所述的实施方式，并不构成对该技术方案保护范围的限定。任何在上述实施方式的精神和原则之内所作的修改、等同替换和改进等，均应包含在该技术方案的保护范围之内。