一种心冲击图信号的处理方法及系统

技术领域

本发明涉及生物医学信号处理技术领域，具体涉及一种心冲击图信号的处理方法及系统。

背景技术

大量临床研究证明，建立夜间连续基础生命体征，如静息心动周期、呼吸努力等数据，是实现心脑血管及呼吸类慢病患者居家管理的关键。然而，现有医疗器械或腕式穿戴设备在监测慢病患者夜间心动周期节律、呼吸努力变化过程，或因束缚性强、体验度差而无法长期居家使用，或因准确性不足无法开展医疗级分析与诊断。

心脏的周期性泵血过程中在人体重力方向将产生微弱且动态变化的力学信号，即心冲击图原始信号（Ballistocardiogram，以下简称BCG），BCG可反映心脏及心血管系统的工作状态。1939年Starr氏研究指出，正常的BCG应与心跳协调一致，具有重复性，且主要包含H、I、J、K、L和M、N波，其中H、I、J、K波的幅值最大，其组合类似于字母W的形状，常称为“W”组合波，如附图1所示。

BCG信号非常容易受到各种噪声和运动伪影的干扰。因此，对抗干扰的算法研究是提高BCG分析准确性的关键。目前常用的BCG心跳特征定位方法有：伪周期定位法、聚类自适应模板匹配法、信号熵值分析的方法来估计心脏跳动间隔。如现有专利文献CN201810104375.0公开了一种面向高血压病人监测的BCG信号分析方法及系统，是通过构建基于类关联规则的分类器，判断该BCG信号是否为高血压信号；CN201310393391.3公开了一种提取BCG信号特征的方法及系统，将累积残余熵算法引入到BCG信号的特征提取分类中，针对BCG信号的不同频段进行特征提取，可以分析得到哪一个频段对疾病的分类有突出贡献，并对分析结果进行分类，以达到准确提取BCG信号中有用的特征从而进行疾病分类诊断的目的。虽然这些方法可以定位出心跳间隔，无论是对于处理不同环境下得到的BCG信号的抗干扰能力，还是上述这些方法在嵌入式设备中进行运算所需要满足的算力，都不足以满足在实际应用中快速有效的定位BCG心跳间隔。

发明内容

针对BCG信号易受干扰，心跳周期定位不够精确等问题，本发明的目的旨在提供一种心冲击图信号的处理方法及系统。

为实现上述目的，本发明采用如下技术方案：

一种心冲击图信号的处理方法，包括以下步骤：

S1、进行BCG信号采集，并将采集到的所述BCG信号存入缓存内存中，得到BCG信号Ⅰ；

S2、对BCG信号Ⅰ进行预处理，通过低通滤波器消除BCG信号中的热噪声、工频干扰以及高频信号噪声，通过形态滤波器组去除BCG信号中低频信号的基线干扰，得到BCG信号Ⅱ；

S3、通过希尔伯特变换，从BCG信号Ⅱ中得到基于J峰特征的BCG信号包络；

S4、通过定位所述BCG信号包络的波峰得到的波峰序列，并计算出所述波峰序列中的最大值坐标，得到BCG信号Ⅱ的J峰坐标序列；

S5、通过所述J峰坐标序列计算出平均心搏周期，基于平均心搏周期构建多尺度BCG模板信号；

S6、在所述BCG信号Ⅰ获取截取出多个样本信号，将每个所述样本信号与所述多尺度BCG模板信号共同进行DTW算法运算，得到最终的BCG的J峰坐标序列，完成心冲击图信号的处理，实现心跳周期定位。

进一步地，所述步骤S1还包括：

S11、进行BCG信号采集，所述BCG信号的频率为1000Hz；

S12、建立容量为15S的缓存空间；

S13、将采集到的所述BCG信号存入所述缓存空间的缓存内存中，得到BCG信号Ⅰ；

S14、通过FIFO的方式实时更新BCG信号Ⅰ，每次更新1s的BCG信号Ⅰ。

进一步地，所述低通滤波器为零相位双线性Butterworth低通滤波器，所述低通滤波器的截止频率为20Hz。

进一步地，所述步骤S2还包括使用两组形态滤波器运算来提取低频基线干扰信号，具体步骤如下：

S21、第一组形态滤波器通过先开后闭运算得到信号MF

S22、把信号MF

S23、将低通滤波后的BCG信号Ⅰ减去基线信号MF(t)后得到的BCG信号为最终滤除低频基线干扰的BCG信号Ⅱ。

进一步地，所述步骤S3还包括：

S31、将所述BCG信号Ⅱ进行希尔伯特变换运算，得到相位位移为90度的BCG信号Ⅲ，具体计算公式如下：

其中

S32、计算所述BCG信号Ⅱ幅值得到所述BCG信号包络，具体公式如下：

其中，e(t)为BCG信号包络，

进一步地，所述步骤S4还包括：

S41、计算所述BCG信号包络出现的所有波峰坐标，组成BCG信号包络的波峰序列init_peak(n)，n为BCG信号包络的波峰序列中的坐标个数，n=1,2……i，

init_peak(n)=[init_peak

S42、计算所述波峰序列init_peak(n)中每个波峰到波谷的距离，得到BCG信号包络的垂直距离序列AP(n)；

S43、计算所述垂直距离序列AP(n)的平均幅值m_ap，将垂直距离小于平均幅值m_ap的波峰坐标从所述波峰序列init_peak(n)中删除，得到新的波峰序列init_peak2(n)；

init_peak2(n)=[init_peak2

S44、设计τ=50ms，以搜索区间为

[init_peak2

init_Jpeak(n)=[init_Jpeak

进一步地，所述步骤S5还包括：

S51、依次将所述J峰坐标序列中相邻坐标进行相减运算得到心跳间期序列init_interval(n);

S52、通过所述心跳间期序列init_interval(n)计算出，平均心跳间期mean_i；

S53、设定多尺度模板信号的长度，具体计算方法为：

M_Len=mean_i*d，M_Len为所述多尺度模板信号的长度，mean_i为平均心跳间期，d为平均心跳间期mean_i的倍数，d=1,2,...m，且满足最大长度mean_i*m小于BCG信号Ⅱ信号长度的一半。

进一步地，所述步骤S6还包括：

S61、设置DTW运算长度区间为l=[0 , 1*Δt , ……，(h-1)*Δt , ……1500]，Δt=10ms，h 为样本信号个数，h∈[0,150];

S62、在BCG信号s(t)上以init_Jpeak

S`

S63、对每个样本S`

S64、计算Min(D(h))得到的最小距离的坐标标记为BCG的J峰坐标；

S65、重复上述计算，遍历所有的init_JPeak(n)所得到的JPeak(n)为最终得到的BCG的J峰坐标序列。

一种心冲击图信号的处理系统，所述处理系统用于实现所述心冲击图信号的处理方法，包括BCG信号采集单元、预处理单元、波峰提取计算单元、多尺度BCG模板信号构建单元和DTW算法运算单元；

所述BCG信号采集单元用于采集和更新BCG信号；

所述预处理单元用于消除BCG信号中的热噪声、工频干扰以及高频信号噪声和低频信号的基线干扰；

所述波峰提取计算单元用于获得BCG信号的J峰坐标序列；

所述多尺度BCG模板信号构建单元用于通过所述J峰坐标序列计算出平均心搏周期和构建多尺度BCG模板信号；

所述DTW算法运算单元用于获取出多个样本信号，将每个所述样本信号与所述多尺度BCG模板信号共同进行DTW算法运算，得到最终的BCG的J峰坐标序列。

本发明的有益效果在于：

1、本发明采用数学形态滤波的方式有效的提取混叠呼吸基线的BCG信号，极大程度的保留了BCG信号在时域上的波形特征，提高BCG信号的信噪比；结合希尔伯特变换运算高效的提取BCG信号的信号包络，实现对BCG信号的心率特征范围的判断。

2、基于预估的心率特征范围构建多尺度的BCG信号模板，结合动态时间规整的匹配算法，提高对BCG信号J峰定位的准确度，以及有效扩大通过BCG信号计算有效心率的范围。

附图说明

图1为BCG信号特征波形图；

图2为本发明方法所提供的BCG心跳周期定位方法流程图；

图3为本发明方法中使用形态滤波器去除低频基线干扰信号的示意图；

图4为本发明方法形态学滤波器去除低频基线干扰信号设计流程图；

图5为本发明一种心冲击图信号的处理系统结构示意图。

具体实施方式

下面，结合附图以及具体实施方式，对本发明做进一步描述：

本发明提出一种心冲击图信号的处理方法，本方法包括以下步骤：

S1、进行BCG信号采集，并将采集到的所述BCG信号存入缓存内存中，得到BCG信号Ⅰ；

步骤S1还包括：

S11、进行BCG信号采集，所述BCG信号的频率为1000Hz；

S12、建立容量为15S的缓存空间；

S13、将采集到的所述BCG信号存入所述缓存空间的缓存内存中，得到BCG信号Ⅰ；

S14、通过FIFO的方式实时更新BCG信号Ⅰ，每次更新1s的BCG信号Ⅰ。

FIFO：简单理解为先进先出，在实际使用过程中描述为：在一个15s的信号缓存数组中，第一秒的数据存进数组位置0，存满15s信号后，第16s信号进来时，将第1s信号去除，第2s信号移到数组位置0，第16s信号位于原第15s信号的数组位置；通过这样的方式可以有效保持信号处理的实时性，节约嵌入式系统的内存使用。

S2：对BCG信号Ⅰ预处理滤波消除信号工频干扰；通过形态学滤波器去除0.5HZ以内的低频信号的基线干扰。

对缓存内存中的BCGⅠ信号进行预处理，包括使用截止频率为20Hz零相位双线性Butterworth低通滤波器消除BCG信号中的热噪声、工频干扰以及高频信号噪声，使用形态滤波器组合分离出BCG信号中0.5Hz以内低频基线干扰信号，最终得到去除高频热噪声和去除低频基线干扰的BCG信号，如图3所示。

预处理前的原始信号，是混叠呼吸运动和心跳BCG信号两种信号叠加起来的。预处理滤波：从硬件出来的信号会带有工频干扰，以及高频的热噪声信号，通过工频滤波去除50hz的工频信号，心冲击BCG信号的主要频率集中在0-20HZ，所以通过巴特沃兹低通滤波器来提取心冲击的主要信号频率，通过软件滤波的方式可以降低硬件滤波的成本，采用双线性零相位的巴特沃兹滤波可以消除IIR方式滤波的时延，同时滤波出来的信号结果更好。

本发明使用两组形态滤波器运算来提取低频基线干扰信号MF(t)，如图1所示，基于在1KHz采样下的有效的BCG的H、I、J、K特征宽度，设定形态学滤波器的结构函数g(m)的信号长度为200ms，数值为零的直线信号；第一组形态滤波器通过先开后闭运算得到信号MF

S3：通过希尔伯特变换（Hilbert transform，HT），从BCG信号Ⅱ中得到基于J峰特征的BCG信号包络，具体计算方式如下：

将15s的BCG信号s(t)进行HT运算，得到90度的相移的BCG信号

其中

根据如下公式计算BCG信号幅值得到BCG信号包络e(t):

其中，

即按照公式（1）可以得到90度的相移的BCG信号

通过计算提取BCG的信号包络，可以对15s缓存内的信号进行心跳间期的估计，得到15s信号的一个平均心跳间期，为后续的BCG信号J峰定位提供一个基准范围，为自适应识别心跳间期提供一个基础，包络信号的峰值，可以在后续进行BCG信号的J峰提取时，提供辅助定位的作用。

S4：通过S3得到的BCG信号包络e(t)，由于信号包络e(t)与BCG信号s(t)在时序上具有一致性，因此通过定位信号包络e(t)的波峰得到的波峰序列init_peak(n)，以50ms的搜索区间，在BCG信号s(t)上计算最大值的坐标，得到BCG的J峰坐标序列init_Jpeak(n)。后续可通过J峰坐标的间距反映心跳周期，

S4：BCG信号心跳特征J峰初始定位实现步骤如下：

S41、计算S3中信号e(t)出现的所有波峰坐标，组成BCG信号包络的波峰序列init_peak(n)，n为BCG信号包络的波峰序列中的坐标个数，n=1,2……i，

init_peak(n)=[init_peak

S42、计算init_peak(n)中每个波峰到波谷的距离得到BCG包络的垂直距离序列AP(n)；

S43、计算序列AP(n)的平均幅值m_ap，将垂直距离小于m_ap的波峰坐标从init_peak(n)中删除，得到新的波峰序列init_peak2(n)，init_peak2(n)=[init_peak2

S44、：设计τ=50ms，对BCG信号s(t)以搜索区间为

[init_peak2

S5：基于平均心搏周期构建多尺度BCG模板

构建BCG信号的多尺度模板信号

S51、将S4中得到的J峰坐标序列进行相邻坐标相减运算得到心跳间期序列init_interval(n)。

S52、计算得到平均心跳间期mean_i；平均心跳间期是指的将init_interval(n)序列中的间期计算平均值，得到平均心跳间期；

S53、以mean_i的倍数，设定多尺度模板信号的长度为M_Len=mean_i*d，M_Len为所述多尺度模板信号的长度，mean_i为平均心跳间期，d为平均心跳间期mean_i的倍数，d=1,2,...m，且满足最大长度mean_i*m小于BCG信号Ⅱ信号长度的一半；

M_Len的最大长度mean_i*d如果越长，最后能计算出来的间期数据就越少，意味着会损失一部分数据结果，同时M_Len的长度越长，需要的运算量也越大。故本发明中，将M_Len的最大长度定义在信号长度的一半。

S54、根据波峰序列init_Jpeak(n)在BCG信号s(t)中,根据J峰序列init_Jpeak(n)，n=0,1,2,...K-1，以每个J峰坐标init_Jpeak

S55、根据设定的多尺度模板长度M，可构建多尺度模板信号

S6：通过动态时间规整进行心跳周期定位。

本发明中使用DTW算法，是因为BCG信号本身是一个伪周期信号，我们可以通过DTW算法的特点来对伪周期信号进行峰值定位，当前周期峰值，与下一个周期的峰值间的距离是不同的，但是波形的特征是相似的，通过DTW算法，可以对这类信号进行更有效的自适应间期区间进行峰值定位。

根据动态时间规整（Dynamic Time Warping，DTW）的算法思路，设置DTW运算长度区间为l=[0 , 1*Δt , ……，(h-1)*Δt , ……1500]，Δt=10ms，h 为样本信号个数，h∈[0,150]，在BCG信号s(t)上以init_Jpeak

S`

在实际应用过程中，不断循环重复步骤S1-S6对BCG信号进行处理可快速有效的定位BCG心跳间隔。

如图5所示，一种心冲击图信号的处理系统，所述处理系统用于实现所述心冲击图信号的处理方法，包括BCG信号采集单元、预处理单元、波峰提取计算单元、多尺度BCG模板信号构建单元和DTW算法运算单元；

所述BCG信号采集单元用于采集和更新BCG信号；

所述预处理单元用于消除BCG信号中的热噪声、工频干扰以及高频信号噪声和低频信号的基线干扰；

所述波峰提取计算单元用于获得BCG信号的J峰坐标序列；

所述多尺度BCG模板信号构建单元用于通过所述J峰坐标序列计算出平均心搏周期和构建多尺度BCG模板信号；

所述DTW算法运算单元用于获取出多个样本信号，将每个所述样本信号与所述多尺度BCG模板信号共同进行DTW算法运算，得到最终的BCG的J峰坐标序列。

对本领域的技术人员来说，可根据以上描述的技术方案以及构思，做出其它各种相应的改变以及形变，而所有的这些改变以及形变都应该属于本发明权利要求的保护范围之内。