CN113576410B - 一种入睡过程动态连续分析方法 - Google Patents

Abstract

本发明公开一种入睡过程动态连续分析方法，应用于电子信息、生物医学工程领域，针对现有技术存在的无法连续、动态地分析睡眠建立过程的问题；本发明基于马尔科夫模型以及粒子滤波，能够实现对入睡过程的连续、动态的分析；并且本发明将脑电信号与专家分期结果相结合，能够很好的体现在入睡过程中脑电各成分功率的变化与入睡程度间的联系，且能够输出稳定可靠的结果，能够有效地体现入睡过程的整体趋势与局部细节。

Description

一种入睡过程动态连续分析方法

技术领域

本发明属于电子信息、生物医学工程领域，特别涉及一种入睡过程动态连续分析技术。

背景技术

在人的睡眠过程中利用生理信息，根据国际睡眠医学标准，可将不同睡眠时间段划分为WAKE、REM、N1、N2和N3期，共五个正常类别。专家们将30秒长度的Epoch作为睡眠状态可辨别的最小时间单位，根据该时段内的EEG，执行人工检查，然后对睡眠状态进行评估。依据临床手册，在这种判断方式下如果要对人的清醒-睡眠的状态进行分类的话，人的睡眠状态只能分为清醒和入睡两种，从时间分辨率的角度来讲，30秒的时间分辨率明显太粗糙，考虑到人从清醒到入睡其实并不是一个特别长的过程，其入睡过程全长普遍不超过20分钟，如此粗糙的时间分辨率会导致大量关于入睡过程的信息被遗漏。这对我们动态地分析人的睡眠建立过程造成了极大的阻碍，从清醒程度的分辨率来讲，在传统的睡眠分期方法下，一个人在同一时间的入睡状态只有清醒和睡着两种，这种二值化的分析使得我们完全无法得知关于睡眠建立过程更多的细节，对睡眠建立过程的分析也就无从谈起。

人类睡眠状态转换过程的分析对于基础神经科学以及睡眠医学都具有极其重要的意义，但是目前对睡眠时间段进行划分的传统标准并不能满足我们对于动态、连续地分析睡眠建立过程所提出的时间分辨率以及清醒程度分辨率的要求，一种连续、动态地分析睡眠建立过程的方法显得极有必要。

发明内容

为解决上述技术问题，本发明提出一种入睡过程动态连续分析方法，结合连续的脑电信号对分析结果进行约束，实现了对入睡过程的连续、动态的分析。

本发明采用的技术方案为：一种入睡过程动态连续分析方法，包括：

S1、载入观测数据，包括专家的睡眠分期结果与脑电数据；

S2、构建多个状态变量；

S3、根据构建的多个状态变量，基于马尔科夫模型，构造当前时刻的初始粒子群；

S4、将当前时刻的所有粒子输入观测模型，估计得到观测数据；

S5、计算观测数据与真实数据之间的误差，并根据误差赋予各粒子权重；

S6、基于各粒子权重对所有粒子进行重采样，将采样出的新粒子群的中心值作为当前状态变量的真实值来计算当前时刻的清醒程度；

S7、重复步骤S3至S6，计算出每个时刻的清醒程度。

步骤S1所述专家的睡眠分期结果包括：清醒期与睡眠期，将清醒期的值设置为0，睡眠期的值设置为1

步骤S2具体构建两个与脑电数据相关的状态变量，具体的：根据脑电数据的三个频段功率α、δ、θ，建立的两个状态变量记为

与

表示从脑电α频段的频谱功率估计觉醒程度的状态变量，

是从脑电δ与θ频段的频谱功率估计觉醒程度的状态变量。

当观测数据包括肌电数据，步骤S2所述状态变量还包括与肌电有关的状态变量。

步骤S2还包括计算当前时刻初始的清醒程度，具体的将构建的多个状态变量通过线性组合后输入一个S型函数，该S型函数的结果即为当前时刻初始的清醒程度。

步骤S5所述根据似然函数分别计算当前时刻初始清醒程度与当前时刻的专家睡眠分期结果之间的误差、各粒子的观测数据与真实数据之间的误差；最后得到的总似然函数的绝对值作为粒子权重。

总似然函数的绝对值越小的粒子权重越大。

步骤S6所述计算当前时刻的清醒程度，具体为：将步骤S6得到的多个状态变量的真实值通过线性组合后输入一个S型函数，该S型函数的结果即为当前时刻的清醒程度。

本发明的有益效果：基于马尔科夫模型以及粒子滤波，能够实现对入睡过程的连续、动态的分析；并且本发明将脑电信号与专家分期结果相结合，能够很好的体现在入睡过程中脑电各成分功率的变化与入睡程度间的联系，且能够输出稳定可靠的结果，能够有效地体现入睡过程的整体趋势与局部细节；采用本发明的方法可以实现对入睡过程的连续、动态分析，对于基础神经科学以及睡眠医学都具有极其重要的意义；本发明的包括以下优点：

1、本发明的最大创新便是实现了对入睡过程的连续、动态的分析；传统的睡眠分期方法由专家进行人工分期，每次对最短15s的脑电信号进行判断，也就是最高的时间分辨率只有15s；而且，这种睡眠分期的结果虽然将睡眠分成了若干个阶段(比如快速眼动期、N1期、N2期、N3期，还有清醒状态)，但从入睡分析的角度来说，在现有分期结果下人的入睡状态只有清醒和已入睡2种状态，也就是说入睡程度的分辨率太低，这样的结果对于分析入睡过程显然不能提供足够的信息；综上所述，传统睡眠分期方法的时间分辨率与入睡程度的分辨率均太低，无法满足对入睡过程进行分析的要求；而本发明结合了连续的脑电信号对分析结果进行约束，其时间分辨率可以高于1s，同时能够连续地显示入睡程度，实现了对入睡过程的连续、动态的分析，非常适合用于分析入睡过程；

2、本发明具有非常好的适用性，本发明的算法将脑电信号与专家分期结果相结合，能够很好的体现在入睡过程中脑电各成分功率的变化与入睡程度间的联系，且能够输出稳定可靠的结果，能够有效地体现入睡过程的整体趋势与局部细节；

3、本发明具有极好的拓展性，可根据实际情况非常方便地添加新的变量或者新的真实数据作为约束条件；如要添加新的变量(如肌电反映的估计觉醒程度的状态变量)只需添加新变量的马尔可夫模型并对公式(1)稍作修改即可；而如果要增加新的真实数据作为约束条件，则只需在计算观测模型以及计算似然函数时分别添加上对应新真实数据的观测函数以及似然函数即可。本发明极好的扩展性大大加宽了本发明的适用范围，是本发明的重要创新。

附图说明

图1为本发明的方法流程图。

具体实施方式

为便于本领域技术人员理解本发明的技术内容，下面结合附图对本发明内容进一步阐释。

本发明利用粒子滤波器的基本原理，以蒙特卡洛方法以及贝叶斯定理为基础，建立若干随机服从马尔可夫模型的状态散点粒子，通过一个观测函数来与实际观测值进行联系，并依据观测函数的结果与实际观测值误差对散点粒子进行重采样，以重采样得到的状态散点粒子来表示当前真实的入睡状态，即清醒程度。

如图1所示，本发明的方法包括以下步骤：

S1、对需要进行分析的观测数据进行处理，标准的观测数据应该包含：专家的睡眠分期结果(将清醒期的值设置为0，睡眠期的值设置为1)以及脑电数据(δ、θ和α频段的功率分别随时间变化的曲线，其中δ、θ和α为脑电的三个波段，即δ(1－3Hz)、θ(4－7Hz)、α(8－13Hz))。例如，入睡起始过程中，专家对睡眠判读后的分期结果将从0变为1，同时α频段功率减小，θ频段功率增加，δ频段功率增加，所以这些数据包含有能够帮助分析睡眠建立过程的特征。如果采用如上的数据构成，则观测数据会以4×N矩阵的形式输入到分析算法当中去(N为数据的长度，4行分别为专家的睡眠分期结果、脑电α频段功率、脑电θ频段功率、脑电δ频段功率)。

需要补充的是，目前的睡眠分期结果大多依靠专家对采集到的脑电数据分析得来，这是一个经验性的工作，根据脑电信号在不同睡眠阶段的不同波形以及幅度等特征，专家会依据R&K标准(将睡眠分为4个非快速眼动期与1个快速眼动期，当然也包括清醒状态)或者AASM标准(将睡眠分为3个非快速眼动期与1个快速眼动期，以及清醒状态)等标准将整个睡眠过程分到符合标准的睡眠期，从而实现对睡眠过程的分析。

由于本发明着重解决对入睡过程的分析问题，所以在本发明中只关注从清醒阶段到睡眠阶段的过程，所以才有了将清醒期的值设置为0，睡眠期的值设置为1的处理方式。值得提及的是，在实际应用中可以使用多种方式来改进这种0-1形式的处理方式，比如使用不同的标准来对专家的分期结果进行这种0-1映射，或者将专家的分期结果映射到0至1区间内更多的离散值上。对于这种处理产生的变量增加，我们只需在下面S4步骤时添加对应变量的似然函数即可。

S2、建立两个状态变量

与

(根据实际可建立更多的状态变量，比如加入肌电数据，即可再建立一个与肌电有关的状态变量

)，

是从脑电α频段的频谱功率估计觉醒程度的状态变量，

是从脑电δ与θ频段的频谱功率估计觉醒程度的状态变量，上述状态变量经过线性组合(即公式(1))后输入一个S形函数(sigmoid)(即公式(2))。该函数的结果即作为当前的清醒程度

此处应该注意的是以上两个状态变量并不是对应信号的功率幅度，根据粒子滤波的特点，可以将其看作t时刻由实际观测结果所隐含的在对应方面关于觉醒程度的真实状态。对于每个时刻t的状态变量，t+1时刻的状态变量为时刻t的状态变量加上一个随机误差后的结果，这是一个典型的马尔可夫模型，在算法的实际运行中，根据粒子滤波器的原理本发明会构造一群状态变量的散点粒子，以上述马尔可夫模型对下一刻粒子进行适度初始化，马尔可夫模型形式如式(3)与式(4)，其中

与

均为构造粒子模型所需引入的随机误差。

S3、对于每个时刻t的粒子，算法中将会把这些粒子输入观测模型当中去，观测模型亦是对状态变量的线性组合使用一个S形函数进行变换，其形式同式(5)至(6)，其目的是计算对于每个状态变量粒子所估计出的观测数据。建立观测模型的意义在于通过状态变量对该时刻实际观测到的数据进行估计，在下面的粒子滤波过程中，通过对估计结果与真实观测结果的差异比较，可以确定各个粒子的权重，根据权重对粒子进行随机采样，对于权重越大的粒子，本发明会赋予其更大的被采样概率。

其中

与

与

与

均分别为对应脑电频段α、δ、θ真实数据功率谱的最小值与最大值。

是一串散点粒子，其作用是加强粒子滤波器的效果，其不同时刻间的变化亦是一个马尔可夫模型，其值一般在0-1之间，亦可根据情况改变其取值范围。在计算出

(其分别代表对脑电α、δ、θ频段功率的估计值)后，下一步便是计算权重后按权重采样。

S4、在这一步骤中，计算各粒子所估计出的观测数据与真实数据之间的误差，根据差异，本发明会对各粒子赋予不同的权重，粒子权重计算公式(即似然函数)的形式同式(8)至(12)：

总似然函数为：

所得的loglikehood的绝对值大小代表各粒子的权重。b_t代表的是专家的睡眠分期结果在t时刻的值，

分别代表t时刻脑电α、δ、θ频段功率的真实值，σ_α、σ_Δ、σ_θ为归一化误差。

对于误差越小，即loglikehood的绝对值越大的粒子本发明赋予越大的权重，根据这个权重归一化后的结果对粒子重采样，则权重越大的粒子被重复采到的概率越大，在对粒子按照原来的数量重新采样后，得到了采样出的新粒子群，新粒子群的中心值被作为当前状态变量的真实值来计算当前的P_wake，并取P_wake的97.5％至2.5％的范围作为置信区间。

S5、根据马尔可夫模型计算出下一刻的初始粒子群，重复S3与S4的步骤计算出每个时刻t的清醒程度

在遍历所有的时刻t以后，输出总结果，总结果就包含了各个时刻t的清醒程度

由于通过S形函数计算出的

是个连续变化的值，于是就实现了对于睡眠建立过程的动态且连续的分析，通过分析

值的变化情况以及各个时刻散点的分布状况，即可科学地对睡眠建立过程进行分析。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。

Claims (4)

1.一种入睡过程动态连续分析方法，其特征在于，包括：

S1、载入观测数据，包括专家的睡眠分期结果与脑电数据；

S2、构建多个状态变量；步骤S2具体构建两个与脑电数据相关的状态变量，具体的：根据脑电数据的三个频段功率α、δ、θ，建立的两个状态变量记为

与

表示从脑电α频段的频谱功率估计觉醒程度的状态变量，

是从脑电δ与θ频段的频谱功率估计觉醒程度的状态变量；状态变量

与

经过公式(1)的线性组合后输入公式(2)的S形函数，S形函数的结果即作为当前的清醒程度

S3、根据构建的多个状态变量，基于马尔科夫模型，构造当前时刻的初始粒子群；对于每个时刻t的粒子，把这些粒子输入观测模型当中去，观测模型亦是对状态变量的线性组合使用一个S形函数进行变换，其形式同式(5)至(6)，其目的是计算对于每个状态变量粒子所估计出的观测数据；建立观测模型的意义在于通过状态变量对该时刻实际观测到的数据进行估计，在下面的粒子滤波过程中，通过对估计结果与真实观测结果的差异比较，确定各个粒子的权重，根据权重对粒子进行随机采样，对于权重越大的粒子，会赋予其更大的被采样概率；

其中，

分别代表对脑电α、δ、θ频段功率的估计值，

与

与

与

分别为对应脑电频段α、δ、θ真实数据功率谱的最小值与最大值，

是一串散点粒子；

S4、将当前时刻的所有粒子输入步骤S3中的观测模型，估计得到观测数据；

S5、计算观测数据与真实数据之间的误差，并根据误差赋予各粒子权重；粒子权重计算公式的形式同式(8)至(12)：

总似然函数为：

所得的loglikehood的绝对值大小代表各粒子的权重，b_t代表的是专家的睡眠分期结果在t时刻的值，

分别代表t时刻脑电α、δ、θ频段功率的真实值，σ_α、σ_Δ、σ_θ为归一化误差；

S6、基于各粒子权重对所有粒子进行重采样，将采样出的新粒子群的中心值作为当前状态变量的真实值来计算当前时刻的清醒程度；具体为：将步骤S6得到的多个状态变量的真实值通过线性组合后输入一个S型函数

该S型函数的结果即为当前时刻的清醒程度；

S7、重复步骤S3至S6，计算出每个时刻的清醒程度。

2.根据权利要求1所述的一种入睡过程动态连续分析方法，其特征在于，步骤S1所述专家的睡眠分期结果包括：清醒期与睡眠期，将清醒期的值设置为0，睡眠期的值设置为1。

3.根据权利要求2所述的一种入睡过程动态连续分析方法，其特征在于，当观测数据包括肌电数据，步骤S2所述状态变量还包括与肌电有关的状态变量。

4.根据权利要求3所述的一种入睡过程动态连续分析方法，其特征在于，总似然函数的绝对值越小的粒子权重越大。

Images