CN104983411A

CN104983411A - 动态脉率变异性信号复杂性测度的实时计算方法

Info

Publication number: CN104983411A
Application number: CN201510416062.5A
Authority: CN
Inventors: 张爱华; 丑永新; 漆宇晟
Original assignee: Lanzhou University of Technology
Current assignee: Lanzhou University of Technology
Priority date: 2015-07-16
Filing date: 2015-07-16
Publication date: 2015-10-21

Abstract

动态脉率变异性信号复杂性测度的实时计算方法，其目的为快速准确地实现脉率变异性信号复杂性测度计算，用于心血管疾病的在线监护与预警；其步骤为：(1)对动态脉搏信号检测与处理；(2)动态脉率变异性信号预处理与提取，记为PP(i)；(3)采用滑窗的方式更新动态脉率变异性信号采样点PP(N)；(4)重构PP(1)和PP(N)所在的时间序列；(5)符号化X(1)和X(N-m+1)序列；(6)产生{S₁(j)}和{S_N _-m+2(j)}编码的存储位置；(7)迭代更新熵值。

Description

动态脉率变异性信号复杂性测度的实时计算方法

技术领域

本发明涉及脉搏信号检测实时处理分析，具体涉及一种基于滑窗迭代基本尺度熵的动态脉率变异性信号复杂性测度实时计算方法。可用于动态脉率实时提取和实时复杂性计算。

背景技术

人体心血管***是一种复杂的动力学***，探测其动力学变化及区别不同时间段或不同生理状态下复杂性的变化，对心血管疾病的诊断和治疗十分重要。

心率变异性信号产生于心脏搏动周期的变化，蕴含着丰富的有关心血管***的生理和病理信息。但心率变异性信号从心电信号获取，繁杂的连线和电极使其在便携可穿戴式医疗仪器中应用受到限制。脉率变异性信号，也产生于心脏搏动周期的变化，与心率变异性信号相似，亦蕴含着大量有关心血管***的生理和病理信息。相比心电信号，脉率变异性信号从脉搏信号中获取，获取过程简单，可方便用于便携可穿戴式医疗仪器。动态脉率变异性信号从动态脉搏信号中提取，相比于离线信号分析，在线的实时分析对心血管疾病的监护和及时预警十分重要。

目前，对脉率变异性分析主要参考心率变异性的分析方法，常用的为时域、频域、时频域及非线性分析方法。脉率变异性信号为非平稳的时变信号，非线性方法可更加有效地提取信号的复杂性变化。一些非线性特征如样本熵、近似熵、符号序列熵、基本尺度熵等信息熵分析方法，可用来计算心率变异性信号复杂性。其中，基本尺度熵可有效地分析心率变异性信号，用于识别冠心病等疾病，取得了很好的效果。但基本尺度熵的运算量随着信号长度的增加而呈指数级增长，且运算中间变量占用***大量的内存，给动态脉率变异性信号复杂性测度的实时计算带来极大困难。

发明内容

本发明的目的是为快速准确地实现脉率变异性信号复杂性测度计算，用于心血管疾病的在线监护与预警。

本发明是动态脉率变异性信号复杂性测度实时计算方法，其步骤为：

（1）对动态脉搏信号检测与处理，动态脉搏信号检测与处理由光电脉搏传感器、脉搏信号检测模块、蓝牙模块和智能手机模块完成；

（2）动态脉搏信号预处理与动态脉率变异性信号提取，记为PP(i)，其中PP(i)为第i个采样值；

（3）采用滑窗的方式更新动态脉率变异性信号采样点PP(N)，其中PP(N)为第N个采样值；

（4）重构PP(1)和PP(N)所在的时间序列，分别记为X(1)和X(N-m+1)，其中m为时间序列长度；

（5）符号化X(1)和X(N-m+1) 序列，分别记为{S ₁(j)}和{S _N-m+2(j)}，其中j为相应序列符号化后第j点值；

（6）产生{S ₁(j)}和{S _N-m+2(j)}编码的存储位置；

（7）迭代更新熵值。

本发明的有益之处是：采用滑窗迭代基本尺度熵分析法是指采用滑窗迭代思想实现的基本尺度熵计算，相比于原基本尺度熵分析法，可以实现单采样点分析。在不影响计算精度的前提下，可提高算法运行速度，节约***内存。滑窗迭代对基本尺度熵的改进方法还可用于其他如符号序列熵等信息熵的改进。

心率变异性信号可以反映人体自主神经***活性,同时也可以用来评估交感神经与迷走神经的平衡性。通过滑窗迭代基本尺度熵分析法实时分析动态脉率变异性信号，得到可以反映心脏搏动变化的基本尺度熵。从而了解自主神经的***的状态，实现疾病的监测与预警。

附图说明

图1是本发明的脉搏信号检测与处理***框图，图2是本发明的滑窗迭代基本尺度熵分析的算法原理图。

具体实施方式

如图1、图2所示，本发明是动态脉率变异性信号复杂性测度实时计算方法，其步骤为：

（6）产生{S ₁(j)}和{S _N-m+2(j)}编码的存储位置；

（7）迭代更新熵值。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，述步骤（1）所述的动态脉搏信号的检测与处理由光电脉搏传感器、脉搏信号检测模块、蓝牙模块和智能手机模块完成。图1所示为动态脉搏信号的检测与处理***框图。实现动态脉搏信号的检测与处理，信号采样频率为250Hz。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（2）所述对采集动态脉搏信号进行预处理与提取，按如下步骤进行：

（1）通过截至频率为62.5Hz的整系数低通滤波器实时滤除肌电干扰和随机噪声；

（2）50Hz及其整数倍谐波陷波器去除基线漂移和工频干扰；

（3）对滤波后的信号采用动态差分阈值法提取动态脉率变异性信号。此处的脉率变异性信号为脉搏信号的主波间期,记为PP(i)，其中PP(i)为第i个采样值，单位为ms。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（3）所述采用滑窗的方式更新动态脉率变异性信号采样点PP(N)，其中PP(N)为第N个采样值，设定数据缓冲区,缓存提取的动态脉率变异性信号，记数据缓冲区长度为N个采样点。以数据缓冲区为窗口，通过滑窗的方式实现动态脉率变异性信号采样点的更新，在更新采样点的同时，通过迭代的方式计算基本尺度熵。按如下步骤进行：

（1）先在内存开辟空间，缓存新的采样点，将其记为PP(N+1)；

（2）剔除最早缓存的采样点PP(1)，高位的采样点向低位移动，PP(i)=PP(i+1)；

（3）再将预先缓存的PP(N+1)放在缓存区的最高位，即PP(N)=PP(N+1)。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（4）所述重构PP(1)和PP(N)所在的时间序列，其中m为时间序列长度，根据基本尺度熵的原理，需要将N个缓存的采样点重构产生N-m+1个长度为m的时间序列，每个时间序列代表一种心脏搏动模式。为了减少重构过程计算量和存储空间，本发明采用迭代方式实现熵值的计算。按如下步骤进行：

（1）将PP(1)所在的时间序列重构为：

(1)

（2）将PP(N)所在的时间序列重构为：

(2)。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（5）所述符号化X(1)和X(N-m+1) 序列，其中j为相应序列符号化后第j点值，按如下步骤进行：

（1）符号化X(1)序列。其中，采用下式对PP(1)符号化：

(3)

式中，μ ₁为时间序列X(1)的均值，α为常数，用来调节符号化边界。BS ₁为基本尺度，用来确定符号化边界。BS ₁为：

(4)

（2）符号化X(N-m+1)序列。过程同(1),记为{S_N-m+2(j)}，j=1, …, m。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（6）所述产生{S ₁(j)}和{S _N-m+2(j)}编码的存储位置。每一种不同的时间序列代表一种不同的心脏搏动模式，故包含4种符号且长度为m的时间序列可代表4^m种心脏搏动模式。统计每一种不同模式占整个N-m+1个模式的几率，用来计算整个时间序列的熵值，于是，需要更新{S ₁(j)}和{S _N-m+2(j)}的出现和个数。其中S ₁编码存储位置为：

(5)。

根据以上所述的动态脉率变异性信号复杂性测度实时计算方法，上述步骤（7）所述迭代更新熵值。

采用迭代的方式实现熵值的更新，故需从上一次计算结果中减去{S ₁(j)}的熵值，加上{S _N-m+2(j)}的熵值，即可实现所有采样点熵值的计算。每一次数据更新之后，{S ₁(j)}所代表的心脏搏动模式个数减1，{S _N-m+2(j)}所代表的模式个数加1。个数分别用n(h)和n(k)表示，如图2示。

更新完{S ₁(j)}和{S _N-m+2(j)}的个数之后，便可迭代更新熵值。记滑窗前熵值为BS′(m)，{S ₁(j)}的个数为n′(h)，出现几率为p′(h)，{S _N-m+2(j)}的个数为n′(k)，出现几率为p′(k)；滑窗后熵值为BS(m)，{S ₁(j)}的个数为n(h)，出现几率为p(h)，{S _N-m+2(j)}的个数为n(k)，出现几率为p(k)。初始BS′=0，迭代过程中n(h)=n′(h)-1，n(k)=n′(k)+1。根据两种模式个数的变化规律，以及熵值计算过程中对数的真数必须大于0；

按如下步骤进行：

（1）记滑窗前熵值为BS′(m)，{S ₁(j)}的个数为n′(h)，计算出现概率，记为p′(h)，{S _N-m+2(j)}的个数为n′(k)，计算出现概率，记为p′(k)；

（2）滑窗后熵值为BS(m)，{S ₁(j)}的个数为n(h)，计算出现概率，记为p(h)，{S _N-m+2(j)}的个数为n(k)，计算出现概率，记为p(k)；

（3）初始BS′=0，迭代过程中n(h)=n′(h)-1，n(k)=n′(k)+1；

（4）迭代更新。根据两种模式个数的变化规律，以及熵值计算过程中对数的真数必须大于0的限制，有以下四种迭代更新方式：

（a）n(h)>0，n(k)>1，表示{S ₁(j)}和{S _N-m+2(j)}所代表的模式在滑窗的前后都存在。于是，在计算熵值的过程中，不需要考虑对数真数为0的情况。由于：

(6)

(7)

式中，M=4^m。滑窗前后，只有{S ₁(j)}和{S _N-m+2(j)}所代表模式的个数发生了变化，其它模式不变，则-p(1)log₂ p(1)=-p′(1)log₂ p′(1)，…，-p(M)log₂ p(M)=-p′(M)log₂ p′(M)。于是，式（6）减式（7）可得：

(8)

进一步对式（8）化简：

(9)

特殊的，当滑窗前后{S ₁(j)}和{S _N-m+2(j)}代表同一种模式时。首先，{S ₁(j)}所代表模式的个数减1，n(h)=n′(h)-1；然后，{S _N-m+2(j)}所代表模式个数加1，n(k)=n′(k)+1=n(h)+1=n′(h)。则式（8）为：

(10)

即，滑窗前后熵值未改变。因为{S ₁(j)}和{S _N-m+2(j)}代表同一种模式时，滑窗前后各模式的个数没有改变。

或（b）n(h)=0，n(k)>1，表示{S ₁(j)}所代表的模式在滑窗之后消失。则n′(h)=1，p(h)=0，p′(h)=1/(N-m+1)。则由式（8）得：

(11)

或（c）n(h)>0，n(k)=1，表示{S _N-m+2(j)}所代表的模式在滑窗后消失。则n′(h)=n(h)+1，n(k)=1，n′(k)=0，p′(k)=0，p(k)=1/(N-m+1)。则由式（8）可得：

(12)

或（d） n(h)=0，n(k)=1，表示{S ₁(j)}所代表的模式消失，{S _N-m+2(j)}所代表的模式第一次出现，故总模式的个数不变。则n′(h)=1，n′(k)=0，p′(k)=p(h)=0，p′(h)=p(k)=1/(N-m+1)。由式（8）可得：

(13)

通过式（9）-式（13）可得到动态脉率变异性信号的基本尺度熵值，通过熵值的变化，实时地实现动态脉率变异性信号复杂性测度的计算。

Claims

1.动态脉率变异性信号复杂性测度实时计算方法，其特征在于，其步骤为：

（6）产生{S ₁(j)}和{S _N-m+2(j)}编码的存储位置；

（7）迭代更新熵值。

2.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（2）所述对采集动态脉搏信号进行预处理与提取，按如下步骤进行：

（2）50Hz及其整数倍谐波陷波器去除基线漂移和工频干扰；

（3）对滤波后的信号采用动态差分阈值法提取动态脉率变异性信号；

此处的脉率变异性信号为脉搏信号的主波间期，记为PP(i)，其中PP(i)为第i个采样值。

3.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（3）所述采用滑窗的方式更新动态脉率变异性信号采样点PP(N)，其中PP(N)为第N个采样值，按如下步骤进行：

4.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（4）所述重构PP(1)和PP(N)所在的时间序列，分别记为X(1)和X(N-m+1)，其中m为时间序列长度，按如下步骤进行：

（1）将PP(1)所在的时间序列重构为：

(1)

（2）将PP(N)所在的时间序列重构为：

(2)。

5.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（5）所述符号化X(1)和X(N-m+1) 序列，分别记为{S ₁(j)}和{S _N-m+2(j)}，其中j为相应序列符号化后第j点值，按如下步骤进行：

（1）符号化X(1)序列：

其中，采用下式对PP(1)符号化：

(3)

式中，μ ₁为时间序列X(1)的均值，α为常数，用来调节符号化边界；BS ₁为基本尺度，用来确定符号化边界；

BS ₁为：

(4)

（2）符号化X(N-m+1)序列：

过程同（1），记为{S_N-m+2(j)}，j=1, …, m。

6.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（6）所述产生{S ₁(j)}和{S _N-m+2(j)}编码的存储位置，其中S ₁编码存储位置为：

(5)。

7.根据权利要求1所述的动态脉率变异性信号复杂性测度实时计算方法，其特征在于步骤（7）所述迭代更新熵值，按如下步骤进行：

（3）初始BS′=0，迭代过程中n(h)=n′(h)-1，n(k)=n′(k)+1；

（4）迭代更新：根据两种模式个数的变化规律，以及熵值计算过程中对数的真数必须大于0的限制，有以下四种迭代更新方式：

（a） n(h)>0，n(k)>1，表示{S ₁(j)}和{S _N-m+2(j)}所表示的模式在滑窗的前后都存在；

熵值迭代公式为：

(6)

特殊的，当滑窗前后，{S ₁(j)}和{S _N-m+2(j)}代表同一种模式时；首先，{S ₁(j)}所代表模式的个数减1，n(h)=n′(h)-1；然后，{S _N-m+2(j)}所代表模式个数加1，n(k)=n′(k)+1=n(h)+1=n′(h)；熵值迭代公式为：

(7)

或（b）n(h)=0，n(k)>1，表示{S ₁(j)}所代表的模式在滑窗之后消失；则n′(h)=1，p(h)=0，p′(h)=1/(N-m+1)；熵值迭代公式为：

(8)

或（c） n(h)>0，n(k)=1，表示{S _N-m+2(j)}所代表的模式在滑窗后消失；则n′(h)=n(h)+1，n(k)=1，n′(k)=0，p′(k)=0，p(k)=1/(N-m+1)；熵值迭代公式为：

(9)

或（d）n(h)=0，n(k)=1，表示{S ₁(j)}所代表的模式消失，{S _N-m+2(j)}所代表的模式第一次出现，故总模式的个数不变；则n′(h)=1，n′(k)=0，p′(k)=p(h)=0，p′(h)=p(k)=1/(N-m+1)；熵值迭代公式为：

(10)

通过式（6）-式（10）可得到动态脉率变异性信号的基本尺度熵值，通过熵值的变化，实时地实现动态脉率变异性信号复杂性测度的计算。