一种自动识别和提取睡眠脑电中纺锤波的方法与流程

本发明涉及脑电信号
技术领域：
，特别涉及一种自动识别和提取睡眠脑电中纺锤波的方法。
背景技术：
：睡眠是一种由神经元的电活动的结构变换而引起的行为状态，它区别于清醒的最大特点在于对外围环境的可逆反应缺失。一个健康的成年人大约需要花费1/3的生命在睡眠上，但是，现如今随着科学技术的进步，生活节奏的加快和来自各方面压力的不断增大，失眠已经成为世界公认的问题。长期失眠容易导致幻想症、记忆障碍、语言障碍、注意力无法集中等问题。因此，近年来关于神经、精神系统疾病与睡眠的关系也得到了广泛的研究，通过睡眠脑电图的监测，可以为癫痫、脑血栓、脑缺氧等中枢神经系统疾病，以及抑郁、焦虑等精神系统疾病提供诊断和治疗等方面的信息。因此对睡眠的监测和研究具有重要的意义。睡眠效果的判定不仅要看睡眠时间的长短，更重要的是由睡眠深度来决定的。因此睡眠分期对睡眠研究有着重要意义。目前国际上将睡眠分为两种状态：非快速眼动期和快速眼动期，非快速眼动期又分为一、二、三期，其睡眠深度依次加深。睡眠进入非快速眼动二期的最重要的标志特征就是脑电图中出现纺锤波和K-复合波，而K-复合波的出现一般都伴随有纺锤波的发生。因此，纺锤波的识别和判定不仅为睡眠分期提供了重要的判别特征，对探索睡眠的生理机制、睡眠障碍、认知活动等均有着重要意义，而且在临床睡眠障碍疾病的检测和疗效评估中有很好的应用前景。纺锤波是一连串幅值大于10μv、频率介于12—16Hz的成串出现明显可辨的波形，持续时间一般大于0.5s，通常以在大脑中央区记录得到的信号幅值最大。传统纺锤波的检测方法主要依靠人工，不仅费时费力而且主观性强；随着对纺锤波的研究，逐渐出现了一些自动识别纺锤波的方法，主要包括：傅立叶变换、AR谱估计、神经网络技术等方法。傅立叶变换和AR谱估计虽可以突出了纺锤波的频域特性，但无法准确判定纺锤波所在的位置；同样神经网络也不能准确地判断纺锤波出现的位置。近两年，国外科研工作者提出了采用数字滤波结合短时傅立叶变换、Teager能量算子、MUSIC谱估计的方法检测纺锤波。但此类方法也都有其各自的缺点，短时傅立叶变换由于其具有时、频分辨率无法同时达到最优的缺点，容易造成波形起止位置判断不准确的问题；而MUSIC算法是一种空间角谱估计的方法，需要已知信号源子空间和噪声子空间的维度，但由于实际情况下信号子空间和噪声子空间相互渗透，故不能对信号源进行有效估计。以上因素均会导致纺锤波检测效果不佳、算法的通用性和鲁棒性较差。技术实现要素：为了克服上述方法的缺陷，本发明的目的在于提供一种自动识别和提取睡眠脑电中纺锤波的方法，本发明基于三重检测的思想，通过小波变换、Teager能量算子、平滑、0/1粗粒化等信号处理技术的有机结合，不仅很好地避免了对纺锤波的错判和漏判，有效地提高了算法的灵敏度和精度，而且去除了噪声的不良影响，提高了算法的抗噪能力和鲁棒性。为了达到上述目的，本发明通过以下具体步骤加以实现：一种自动识别和提取睡眠脑电中纺锤波的方法；包括以下步骤：步骤一、读取睡眠脑电信号x(t)，并对x(t)做0.5-30Hz的带通滤波，提取脑电的主要频段，避免高频噪声信号的干扰；随后对信号进行重新采样，由原来的200Hz降采样为128Hz；对降采样之后的信号进行分段，段长为30s，对每一段的数据利db4小波进行三层小波分解，依据纺锤波的主要频率分布范围，选择第3层的低频系数进行重构，得到预处理后脑电信号y(t)；步骤二、根据纺锤波频段主要分布在12—16Hz的特点，通过小波分解得到小波系数矩阵，根据小波系数矩阵判断小波能量在13—15Hz频段内是否达到最大，若该频段的小波能量达到最大，则初步判定为纺锤波，并根据小波系数矩阵的时间频率对应关系确定该波形出现的大体位置；步骤三、对预处理后脑电信号y(t)进行Teager能量算子TEO的处理，并提取Teager能量算子TEO的绝对值序列o(n)；具体为：Teager能量算子的定义为：z(n)＝y(n)2-y(n+1)y(n-1)o(n)＝|z(n)|其中，y(n)中第n个样本信号，z(n)代表Teager能量算子(TEO)的输出序列，o(n)为TEO的绝对值序列；步骤四、对Teager能量算子(TEO)的绝对值序列o(n)进行平滑处理，方法为：加矩形窗求平均值,窗长为64，步长为1，得到新的序列s(n)；步骤五、将平滑后的序列s(n)进行粗粒化转化为0/1序列p(n)；0/1粗粒化的方法为：p(n)=1,s(n)≥threshold0,s(n)<threshold]]>其中s(n)为平滑后TEO序列，threshold为设定的阈值，其定义为：threshold＝0.8×avg(s(n))，avg(s(n))为s(n)的平均值，p(n)为粗粒化转换后所得的新的0/1序列；步骤六、检测0/1序列p(n)中连续为1的持续时间，若持续时间大于0.5s且小于3s，则该位置的信号为纺锤波信号，并记录其起始和终止的位置；步骤七、对预处理后脑电信号y(t)再次分段，段长为3s，对每一段的数据计算中心频率SEF50，若中心频率大于10.7Hz，则认为该段数据中占据主要地位的波形为纺锤波；步骤八、对步骤二、步骤六和步骤七的三重检测结果进行对比，只有同时满足三重检测的要求才判定待检波形为纺锤波。综上所述，本发明对纺锤波的判定采用三重检测标准：首先，采用小波变换分解得到小波系数，判断小波能量是否在纺锤波频段能量达到最大；其次，采用小波分解和重构的方法提取纺锤波所在对应频段的重构信号，对重构信号采用Teager能量算子、平滑、0/1粗粒化等处理，并判断持续时间是否超过0.5秒；最后，计算重构信号的半功率值(SEF50)，判断SEF50是否大于阈值。只有在同时满足三重检验要求的情况下，待检波形才被判别为纺锤波。本发明能够准确地自动识别和判定睡眠脑电中的纺锤波，同时能够准确地确定纺锤波的起止位置，从而实现纺锤波的自动提取。本发明算法简便，抗干扰能力强。本发明为非快速眼动2期的识别和纺锤波的相关研究提供了新的方法和思路。附图说明图1是本发明算法总体流程图。图2是基于小波分析的能量分布图，图中显示纺锤波频段能量达到最大。图3是脑电原始信号和小波分解重构信号。图4是TEO绝对值序列。图5是平滑后的TEO序列。图6是0/1粗粒化后所得序列。图7是连续为1的持续时间大于0.5s的0/1序列。图8是小波重构信号的SEF50值。图9是检测出的纺锤波。图10是检测出的纺锤波的放大示意图。具体实施方式下面结合附图对本发明进一步进行说明。一种自动识别和提取睡眠脑电中纺锤波的方法，其实现的具体实施步骤参照图1：步骤一、读取睡眠脑电信号x(t)，并对x(t)做0.5-30Hz的带通滤波，提取脑电的主要频段，避免高频噪声信号的干扰；随后对信号进行重新采样，由原来的200Hz降采样为128Hz；对降采样之后的信号进行分段，由于睡眠分期一般以30s为一期，因此分段的段长为30s；根据纺锤波的频域特点：频率一般在12-16Hz,因此选用db4小波进行三层的小波分解，选择第3层的低频系数进行重构，得到8-16Hz频段的脑电信号，图3就是脑电原信号和小波重构后的信号。步骤二、根据纺锤波频段主要在12-16Hz的特点，同时为了与α波进行区分，需要判断小波能量是否在13-15Hz达到最大，若小波能量在13-15Hz频段范围内达到最大，则认为存在可能纺锤波，并可判定出波形出现的大体位置。图2为基于小波分析的能量分布图，图中显示纺锤波频段能量达到最大。步骤三、对小波重构信号进行Teager能量算子(TEO)的计算并取其绝对值，得到Teager能量算子(TEO)的绝对值序列；其中TEO的定义为：z(n)＝y(n)2-y(n+1)y(n-1)o(n)＝|z(n)|其中，y(n)中第n个样本值信号，z(n)代表Teager能量算子(TEO)的输出序列，o(n)为TEO的绝对值序列。图4为TEO绝对值序列。步骤四、为了能更好地检测出TEO的变化情况，需要对TEO序列进行平滑处理，即矩形窗为64，步长为1，为了得到和TEO序列同样长度的新序列，在加窗前需要在原始数据后补充63个值，补得值为最后一点的幅值。图5为平滑后的TEO序列。步骤五、将平滑后的TEO序列粗粒为0/1序列，粗粒化的标准为：其中s(n)为TEO平滑处理后序列，threshold为设定的阈值，其定义为：threshold＝0.8×avg(s(n))，avg(s(n))为s(n)的平均值，p(n)为粗粒化转换后的新的0/1序列，该阈值的设定是通过多次实验验证得到的。图6为0/1粗粒化转换后所得的结果。0/1粗粒化的方法为：p(n)=1,s(n)≥threshold0,s(n)<threshold]]>其中s(n)为平滑后TEO序列，threshold为设定的阈值，其定义为：threshold＝0.8×avg(s(n))，avg(s(n))为s(n)的平均值，p(n)为粗粒化转换后所得的新的0/1序列。步骤六、由于纺锤波的持续时间定义为0.5s-3s，因此检测0/1序列中连续为1的持续时间，若持续时间大于0.5s且小于3s，则信号为纺锤波，并记录其起始和终止的位置。图7为持续时间大于0.5s且小于3s的0/1序列结果。步骤七、根据纺锤波持续时间一般大于0.5s且小于3s的特性，对预处理之后的信号再次分段，段长为3s，对每一段的数据计算中心频率SEF50。由于纺锤波的频段范围和α波的频段范围有重叠，所以依据中心频率区分两种波：若中心频率大于10.7Hz，则认为该段数据中占据主要地位的为纺锤波而非α波；若中心频率小于10.7Hz，则认为该段数据中主要为α波而非纺锤波。图8为中心频率大于10.7Hz的数据段的中心频率值。步骤八、对步骤二、步骤六和步骤七的三重检测结果进行对比，只有同时满足三重检测的要求才判定待检波形为纺锤波；只要有任意一重检测不满足要求，则待测波形判定为非纺锤波。图9为算法最终检测出的纺锤波。图10为最终检测得到的纺锤波的放大示意图。当前第1页1 2 3