CN110464337A - 一种基于极值能量分解法的心率变异性信号分析方法 - Google Patents

一种基于极值能量分解法的心率变异性信号分析方法 Download PDF

Info

Publication number
CN110464337A
CN110464337A CN201910841605.6A CN201910841605A CN110464337A CN 110464337 A CN110464337 A CN 110464337A CN 201910841605 A CN201910841605 A CN 201910841605A CN 110464337 A CN110464337 A CN 110464337A
Authority
CN
China
Prior art keywords
signal
extreme value
component
energy
mode function
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910841605.6A
Other languages
English (en)
Other versions
CN110464337B (zh
Inventor
周作建
宁新宝
王斌斌
姜晓东
王�华
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Jiangsu Huakang Information Technology Co Ltd
Original Assignee
Jiangsu Huakang Information Technology Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Jiangsu Huakang Information Technology Co Ltd filed Critical Jiangsu Huakang Information Technology Co Ltd
Priority to CN201910841605.6A priority Critical patent/CN110464337B/zh
Publication of CN110464337A publication Critical patent/CN110464337A/zh
Priority to PCT/CN2019/121416 priority patent/WO2021042591A1/zh
Priority to AU2019101756A priority patent/AU2019101756A4/en
Application granted granted Critical
Publication of CN110464337B publication Critical patent/CN110464337B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/0205—Simultaneously evaluating both cardiovascular conditions and different types of body conditions, e.g. heart and respiratory condition
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/024—Detecting, measuring or recording pulse rate or heart rate
    • A61B5/02405—Determining heart rate variability
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/024—Detecting, measuring or recording pulse rate or heart rate
    • A61B5/0245—Detecting, measuring or recording pulse rate or heart rate by using sensing means generating electric signals, i.e. ECG signals
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
    • A61B5/316—Modalities, i.e. specific diagnostic methods
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
    • A61B5/316—Modalities, i.e. specific diagnostic methods
    • A61B5/318—Heart-related electrical modalities, e.g. electrocardiography [ECG]
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7203—Signal processing specially adapted for physiological signals or for diagnostic purposes for noise prevention, reduction or removal
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235—Details of waveform analysis
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235—Details of waveform analysis
    • A61B5/725—Details of waveform analysis using specific filters therefor, e.g. Kalman or adaptive filters

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Surgery (AREA)
  • General Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Molecular Biology (AREA)
  • Physics & Mathematics (AREA)
  • Animal Behavior & Ethology (AREA)
  • Pathology (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Cardiology (AREA)
  • Physiology (AREA)
  • Signal Processing (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Psychiatry (AREA)
  • Pulmonology (AREA)
  • Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)

Abstract

本发明公开了一种基于极值能量分解法的心率变异性信号分析方法，包括获取给定时间和给定采样频率下的未知状态的ECG信号，去噪后得到RRI信号x(t)；将RRI信号x(t)作为原始信号，将原始信号x(t)分解为n个极值模态函数分量和一个余量，将原始信号x(t)分解得的n个极值模态函数分量，代表了原始信号不同频段的分量，根据n个极值模态函数分量判定该RRI信号是否为异常心率变异性信号。本发明采用极值能量分解方法分析RRI信号，将原始信号分解为多个分量，也就是极值分量函数，计算每一个分量的能量，得到其能量分布。

Description

一种基于极值能量分解法的心率变异性信号分析方法
技术领域
本发明涉及一种心电图信号分析，尤其涉及一种基于极值能量分解法的心率变异性信号分析方法。
背景技术
生理信号是由生命体多个***相互作用产生的，不同***作用的时间和强度不同，导致生理信号具有时间和空间上的复杂性。心率变异性(HRV)是指测量连续心动周期之间的时间变异数，准确地说，应该是测量连续出现的正常P-P间期之间的差异的变异数。但由于P波不如R波明显或P波顶端有时宽钝，所以对心率变异性信号的研究通常用与P-P间期相等的R-R间期信号(RRI)来代替。研究表明，HRV可做为植物神经***活动的无创性检测指标，尤其在判断某些心血管疾病的预后方面有着重要意义。
现有技术中研究心脏的长时调节节律(<1Hz))常用心律变异性信号(HRV信号)作为分析对象。大量的研究表明，人体HRV信号具有长时相关性和非线性动力学复杂性，并且年龄和疾病会导致动力学复杂性降低。对心率变异性(Heart Rate Varibility,HRV)信号的研究常用的是RR间期(Interbeat Intervals，RRI)信号，即连续RRI信号R波之间的时间间隔信号。
研究HRV信号的能量改变最常用的方法是功率谱分析(PSD)。PSD通过傅里叶变换将HRV信号的功率转化成频率的函数，研究不同频域范围的功率大小，通常HRV频谱分为高频(HF)、低频(LF)和极低频(VLF)等几个频段。LF/HF比值有重要的临床价值。心脏疾病会引起HRV功率谱的改变，比如心衰和心肌梗死引起归一化HF增高、LF和VLF降低。然而，PSD不是一种基于数据驱动的方法，并且对频域的分段比较粗糙，导致细节缺失，分割也不够灵活。
因此，亟待解决上述问题。
发明内容
发明目的：本发明的目的是提供一种可采用较少数据即可直观反映心电图能量分布和信号波动的真实规律的基于极值能量分解法的心率变异性信号分析方法。
技术方案：为实现以上目的，本发明公开了一种基于极值能量分解法的心率变异性信号分析方法，包括如下步骤：
(1)、获取给定时间和给定采样频率下的未知状态的ECG信号，然后对ECG信号进行去噪预处理，从中提取RRI信号，得到未知状态的RRI信号x(t)；
(2)、将RRI信号x(t)作为原始信号，求出原始信号的所有局部极值点，然后将原始信号的所有极大值点和所有极小值点采用样条曲线连起来分别形成上包络线emax和下包络线emin，得到上、下包络线的包络均值信号m(t)＝(emax+emin)/2；
(3)、将原始信号x(t)减去包络均值信号m(t)，得到h(t)＝x(t)-m(t)；然后判断h(t)是否满足极值模态函数的判定条件，如果不满足，将h(t)作为原始信号返回至步骤(3)，直到hk(t)满足极值模态函数的判定条件，则记c1(t)＝hk(t)，作为第一个极值模态函数分量；
(4)、将原始信号x(t)减去第一个极值模态函数分量c1(t)，得到余量r1(t)＝x(t)-c1(t)，然后判断hk(t)是否满足停止准则，如果不满足，将r1(t)作为新的原始序列x(t)，返回至步骤(2)和(3)；如果hk(t)满足停止准则但n＜8时，返回步骤(1)重新获取原始信号；如果hk(t)满足停止准则且n≥8时，得到第2、3、…、n个极值模态函数分量及余量rn(t)，于是将原始信号x(t)分解为n个极值模态函数分量和一个余量，即
(5)、对极值模态函数分量ci(t)，i＝1，2，…，n，进行频谱分析得到各极值模态函数分量的中心频率；
(6)、将原始信号x(t)分解得的n个极值模态函数分量，代表了原始信号不同频段的分量，然后计算其各个分量的能量
Ei＝∫|ci(t)|2dt，i＝1，2，…，n
将每一个能量值归一化，得到归一化的能量分布向量
pi＝Ei/E，i＝1，2，...，n
其中，第一个分量p1表示最高频段的能量，代表了信号在最高频段范围内能量分布的比例，最后一个分量pn表示信号在最低频段范围内能量分布的比例；根据归一化的能量分布向量绘制归一化能量分布图，其中横坐标表示分量层次，纵坐标表示归一化的能量分布向量值，曲线表示平均值，误差棒表示标准差；
(7)、选取第二个分量p2至第七个分量p7，计算能量差异值EDV，EDV＝(p2+p3+p4)-(p5+p6+p7)，当EDV≤第一标准值时，则判定该RRI信号为正常心率变异性信号，当第一标准值＜EDV＜第二标准值时，则判定该RRI信号为疑似异常心率变异性信号，当EDV≥第二标准值时，则判定该RRI信号为异常心率变异性信号。
其中，所述原始信号x(t)所需最少数据量N＝2n+1，其中n为分解出的极值模态函数分量的数量。
优选的，所述步骤(1)中去噪预处理的具体方法为：将ECG信号经过40Hz零相位FIR低通滤波器滤波消除高频噪声，然后经过中值滤波器去除基线漂移。
再者，所述步骤(3)中极值模态函数的判定条件为：(a)、在整个数据序列中，极值点的数量与过零点的数量相等或者相差一个；(b)、在任意时刻，上下包络线对于时间轴对称。
进一步，所述步骤(4)中hk(t)满足停止准则公式为：
ε表示筛选门限，取0.2～0.3之间。
优选的，所述步骤(7)中的第一标准值为-0.15，第二标准值为0.08。
有益效果：与现有技术相比，本发明具有以下显著优点：
本发明采用极值能量分解方法(Extremum Energy Decomposition,EED)分析RRI信号，将原始信号分解为多个分量，也就是极值分量函数，计算每一个分量的能量，得到其能量分布；本发明的可依据生RRI信号自身的波动规律将信号分解为从高频到低频的不同时间层次信号，对频段的分割较为细致；极值分解在所有层次上得到的数据长度相同，因而不会导致数据长度减小，从而使其可以用于短时间数据分析，即需要很少数据量即可分析得到准确结果；EED对于不同层次分量能量分析不容易受噪声的影响。
附图说明
图1为本发明中原始信号的示意图；
图2为本发明中原始信号求取包络线的示意图；
图3为本发明中原始信号的减去包络均值信号的示意图；
图4为本发明中得到第一个极值模态函数分量的示意图；
图5为本发明中极值能量分解法的流程示意图；
图6为本发明实施例1中RRI信号的EED分解示意图；
图7为本发明实施例1中RRI信号中归一化能量分布示意图。
具体实施方式
下面结合附图对本发明的技术方案作进一步说明。
如图1、图2、图3、图4和图5所示，本发明的一种基于极值能量分解法的心率变异性信号分析方法，包括如下步骤：
(1)、获取给定时间和给定采样频率下的未知状态的ECG信号，然后对ECG信号进行去噪预处理，从中提取RRI信号，得到未知状态的RRI信号x(t)；其中去噪预处理的具体方法为：由于ECG能量主要集中在0～40Hz，将ECG信号经过40Hz零相位FIR低通滤波器滤波消除高频噪声，然后经过中值滤波器去除基线漂移；
(2)、将RRI信号x(t)作为原始信号，原始信号x(t)所需最少数据量N＝2n+1，其中n为分解出的极值模态函数分量的数量；求出原始信号的所有局部极值点，然后将原始信号的所有极大值点和所有极小值点采用样条曲线连起来分别形成上包络线emax和下包络线emin，得到上、下包络线的包络均值信号m(t)＝(emax+emin)/2；
(3)、将原始信号x(t)减去包络均值信号m(t)，得到h(t)＝x(t)-m(t)；然后判断h(t)是否满足极值模态函数的判定条件，如果不满足，将h(t)作为原始信号返回至步骤(2)，直到hk(t)满足极值模态函数的判定条件，则记c1(t)＝hk(t)，作为第一个极值模态函数分量；其中极值模态函数的判定条件为：(a)、在整个数据序列中，极值点的数量与过零点的数量相等或者相差一个；(b)、在任意时刻，上下包络线对于时间轴对称；
(4)、将原始信号x(t)减去第一个极值模态函数分量c1(t)，得到余量r1(t)＝x(t)-c1(t)，然后判断hk(t)是否满足停止准则，如果不满足，将r1(t)作为新的原始序列x(t)，返回至步骤(2)和(3)；如果hk(t)满足停止准则但n＜8时，返回步骤(1)重新获取原始信号；如果hk(t)满足停止准则且n≥8时，得到第2、3、…、n个极值模态函数分量及余量rn(t)，于是将原始信号x(t)分解为n个极值模态函数分量和一个余量，即
其中hk(t)满足停止准则公式为：
ε表示筛选门限，取0.2～0.3之间；满足停止准则的极值模态分解则满足如下两个条件：(a)最后得到的极值模态函数分量cn(t)或者余量rn(t)小于预先设定的阈值；(b)残余信号rn(t)成为单调信号，不能从中再提取出极值模态函数信号；
(5)、对极值模态函数分量ci(t)，i＝1，2，…，n，进行频谱分析得到各极值模态函数分量的中心频率；
(6)、将原始信号x(t)分解得的n个极值模态函数分量，代表了原始信号不同频段的分量，然后计算其各个分量的能量
Ei＝∫|ci(t)|2dt，i＝1，2，…，n
将每一个能量值归一化，得到归一化的能量分布向量
pi＝Ei/E，i＝1，2，...，n
其中，第一个分量p1表示最高频段的能量，代表了信号在最高频段范围内能量分布的比例，最后一个分量pn表示信号在最低频段范围内能量分布的比例；根据归一化的能量分布向量绘制归一化能量分布图，其中横坐标表示分量层次，纵坐标表示归一化的能量分布向量值，曲线表示平均值，误差棒表示标准差；
(7)、选取第二个分量p2至第七个分量p7，计算能量差异值EDV，EDV＝(p2+p3+p4)-(p5+p6+p7)，当EDV≤第一标准值时，则判定该RRI信号为正常心率变异性信号，当第一标准值＜EDV＜第二标准值时，则判定该RRI信号为疑似异常心率变异性信号，当EDV≥第二标准值时，则判定该RRI信号为异常心率变异性信号；其中的第一标准值为-0.15，第二标准值为0.08。
本发明采用的极值能量分解方法(Extremum Energy Decomposition，EED)，其是基于极值模态函数的概念的方法，极值模态函数是同时满足下面两个条件的具有单一频率的一类信号，两个条件为：
(a)、在整个数据序列中，极值点(包括极大值和极小值)的数量与过零点的数量必须相等或者最多相差一个；
(b)在任意时刻，局部极大值点形成的上包络线与局部极小值点形成的下包络线的均值为零，也就是说局部上下包络线对于时间轴局部对称；
上面两个条件，条件(a)类似于高斯正态平稳过程对于传统窄带的要求，条件(b)保证了由极值模态函数计算得到的瞬时频率有物理意义。
本发明的极值模态函数分解终止的标准选择要适中，条件太严格，会导致最后几个极值模态函数分量失去意义；条件太宽松，会导致有用的分量丢失；在实际应用中，也可以根据需求设定需要分解的极值模态函数分量层数，当满足分解层数时即终止计算。
实施例1
将EED分析方法用于分析健康人和CHF患者ECG在不同层次下的能量分布。
健康人的一种基于极值能量分解法的心率变异性信号分析方法，包括如下步骤：
(1)、从physionet的RR间期数据库nsr2db获取健康人群的ECG信号；其中数据包含54个健康人(年龄28～76，平均61)，然后对ECG信号进行去噪预处理，从中提取RRI信号，得到RRI信号x(t)；其中去噪预处理的具体方法为：由于ECG能量主要集中在0～40Hz，将ECG信号经过40Hz零相位FIR低通滤波器滤波消除高频噪声，然后经过中值滤波器去除基线漂移；
(2)、将RRI信号x(t)作为原始信号，原始信号x(t)所需最少数据量N＝2n+1＝29，其中n为分解出的极值模态函数分量的数量，n＝8；求出原始信号的所有局部极值点，然后将原始信号的所有极大值点和所有极小值点采用样条曲线连起来分别形成上包络线emax和下包络线emin，得到上、下包络线的包络均值信号m(t)＝(emax+emin)/2；
(3)、将原始信号x(t)减去包络均值信号m(t)，得到h(t)＝x(t)-m(t)；然后判断h(t)是否满足极值模态函数的判定条件，如果不满足，将h(t)作为原始信号返回至步骤(2)，直到hk(t)满足极值模态函数的判定条件，则记c1(t)＝hk(t)，作为第一个极值模态函数分量；其中极值模态函数的判定条件为：(a)、在整个数据序列中，极值点的数量与过零点的数量相等或者相差一个；(b)、在任意时刻，上下包络线对于时间轴对称；
(5)、将原始信号x(t)减去第一个极值模态函数分量c1(t)，得到余量r1(t)＝x(t)-c1(t)，将r1(t)作为新的原始序列x(t)，返回至步骤(2)和(3)，得到第2、3、…、8个极值模态函数分量及余量r8(t)，于是将原始信号x(t)分解为8个极值模态函数分量和一个余量，即
得到如图6所示的健康人的极值模态分解示意图，可以看到分量1具有最高的频率，信号在最短的时间尺度上波动，随着分量序号增加，频率逐渐降低；
(6)、对极值模态函数分量ci(t)，i＝1，2，…，8，进行频谱分析得到各极值模态函数分量的中心频率，得到频域分析结果图；
(7)、将原始信号x(t)分解得的8个极值模态函数分量，代表了原始信号不同频段的分量，然后计算其各个分量的能量
Ei＝∫|ci(t)|2dt，i＝1，2，…，8
将每一个能量值归一化，得到归一化的能量分布向量
pi＝Ei/E，i＝1，2，...，8
其中，第一个分量p1表示最高频段的能量，代表了信号在最高频段范围内能量分布的比例，最后一个分量pn表示信号在最低频段范围内能量分布的比例。
(8)、选取第二个分量p2至第七个分量p7，计算能量差异值EDV，EDV＝(p2+p3+p4)-(p5+p6+p7)＝-0.2092±0.2940。
CHF患者的一种基于极值能量分解法的心率变异性信号分析方法，包括如下步骤：
(1)、从physionet的RR间期数据库chf2db数据库获取ECG信号，其中chfdb数据库包含29个充血性心衰(Congestive Heart Failure，CHF)患者(年龄34～79，平均55)，然后对ECG信号进行去噪预处理，从中提取RRI信号，得到RRI信号x(t)；其中去噪预处理的具体方法为：由于ECG能量主要集中在0～40Hz，将ECG信号经过40Hz零相位FIR低通滤波器滤波消除高频噪声，然后经过中值滤波器去除基线漂移；
(2)、将RRI信号x(t)作为原始信号，原始信号x(t)所需最少数据量N＝2n+1＝29，其中n为分解出的极值模态函数分量的数量，n＝8；求出原始信号的所有局部极值点，然后将原始信号的所有极大值点和所有极小值点采用样条曲线连起来分别形成上包络线emax和下包络线emin，得到上、下包络线的包络均值信号m(t)＝(emax+emin)/2；
(3)、将原始信号x(t)减去包络均值信号m(t)，得到h(t)＝x(t)-m(t)；然后判断h(t)是否满足极值模态函数的判定条件，如果不满足，将h(t)作为原始信号返回至步骤(2)，直到hk(t)满足极值模态函数的判定条件，则记c1(t)＝hk(t)，作为第一个极值模态函数分量；其中极值模态函数的判定条件为：(a)、在整个数据序列中，极值点的数量与过零点的数量相等或者相差一个；(b)、在任意时刻，上下包络线对于时间轴对称；
(4)、将原始信号x(t)减去第一个极值模态函数分量c1(t)，得到余量r1(t)＝x(t)-c1(t)，将r1(t)作为新的原始序列x(t)，返回至步骤(2)和(3)，得到第2、3、…、8个极值模态函数分量及余量r8(t)，于是将原始信号x(t)分解为8个极值模态函数分量和一个余量，即
(5)、对极值模态函数分量ci(t)，i＝1，2，…，8，进行频谱分析得到各极值模态函数分量的中心频率；
(6)、将原始信号x(t)分解得的8个极值模态函数分量，代表了原始信号不同频段的分量，然后计算其各个分量的能量
Ei＝∫|ci(t)|2dt，i＝1，2，…，8
将每一个能量值归一化，得到归一化的能量分布向量
pi＝Ei/E，i＝1，2，...，8
其中，第一个分量p1表示最高频段的能量，代表了信号在最高频段范围内能量分布的比例，最后一个分量pn表示信号在最低频段范围内能量分布的比例；根据健康人和CHF患者的归一化的能量分布向量绘制归一化能量分布图，其中横坐标表示分量层次，纵坐标表示归一化的能量分布向量值，曲线表示平均值，误差棒表示标准差；
(7)、选取第二个分量p2至第七个分量p7，计算能量差异值EDV，EDV＝(p2+p3+p4)-(p5+p6+p7)＝0.2642±0.4070。
本发明为了进一步分析了健康人和CHF患者HRV信号不同层次分量的平均中心频率，得到表1的结果。
表1健康人和CHF患者HRV信号不同分量层次的平均中心频率
可以得出随着分量层次增加HRV的中心频率逐渐降低。
传统的功率谱密度(Power Spectral Density，PSD)方法对频域分割的典型方式为：HF(0.15～0.4Hz)，LF(0.04～0.15Hz)，VLF(0.0033～0.04Hz)。本发明EED方法对分量层次1，频率高于HF；层次2在HF的频率范围；层次3、4在LF的范围；层次5～7在VLF的频率范围；层次8低于VLF。另外可以看到相同层次时CHF患者的频率略高于健康人，反应了心脏疾病对HRV波动节律的影响，同样层次下CHF患者HRV波动更快。
本发明将二个群体的HRV信号分解得的极值模态函数分量Ci(t)，计算获得归一化能量分布向量，并作EED曲线图，如图7所示。图7为健康人和CHF患者RR间期信号EED分析结果示意图，其中数据长度为10000点，曲线表示均值，误差棒表示标准差。曲线上方的符号*表示两组人群能量T检验p<0.01。在层次选择上，去除了层次1和高于7的层次，包括余量。层次1容易受到噪声的影响，导致能量较大的波动，引起结果标准差过大，并且其频率可高至几kHz以上，因而没有明确的生理意义；高于7的层次反应了信号的长时节律，非常容易受到外界环境的影响，并且其频率很低，生理意义不明。图7中在分量低层次(层次2、3)上，CHF患者的归一化能量值明显高于健康人，在分量高层次(层次>5)上，发生相反的变化，健康人能量高于CHF患者能量；健康人能量在层次2～5比较平稳，层次大于5时，能量缓慢上升，而CHF患者在层次2～4迅速下降，之后趋于稳定；在4个层次(2、3、6、7)上，二者能量有显著性区别(p<0.01)。为了对比，本发明增加了替代数据(Healthy Surrogate,CHF Surrogate)的EED分析结果，如图7所示，替代数据是通过将原始数据随机化的方法生成，替代数据的能量分布随着尺度增加单调下降，相比较CHF患者，在小时间尺度上能量更高，在长时间尺度上能量更低。
为了评估EED曲线在低层次和高层次分解上的能量分布差异，计算各人群的能量差异值EDV，高的EDV值表示RRI信号更高的分量低层次能量分布和更低的分量高层次能量分布。计算得到健康人、CHF患者以及他们替代数据的EDV值，如表2所示。
表2健康人和CHF患者的EDV值
*代表健康人和CHF患者T检验结果p<0.0001。
**代表替代数据与其原始数据T检验结果p<0.0001。
从表2可以看出，健康人和CHF患者EDV值具有显著的差别，CHF患者EDV值比健康人高很多。健康人的EDV值小于0，说明在分量高层次有更高的能量，预示着高层次分量有着更高的调节强度。两个人群的EDV值与他们的替代数据EDV值均有显著性差别，人体HRV的EDV明显小于随机序列。
心率变异性(HRV)是指测量连续心动周期之间的时间变异数，准确地说，应该是测量连续出现的正常P-P间期之间的差异的变异数。但由于P波不如R波明显或P波顶端有时宽钝，所以对心率变异性信号的研究通常用与P-P间期相等的R-R间期信号(RRI)来代替。研究表明，HRV可做为植物神经***活动的无创性检测指标，尤其在判断某些心血管疾病的预后方面有着重要意义。
在临床与实际应用中，本发明提出可采用512(n＝8)个连续RR间期的短时心搏信号用于上述的EED分析是有效的，并有一定的数据余量；由上述研究表明，EED分解达到7个分量层次(n＝7)已有很好的结果，则所需数据量最少可取N＝2n+1＝27+1＝256。

Claims (6)

1.一种基于极值能量分解法的心率变异性信号分析方法，其特征在于，包括如下步骤：
(1)、获取给定时间和给定采样频率下的未知状态的ECG信号，然后对ECG信号进行去噪预处理，从中提取RRI信号，得到未知状态的RRI信号x(t)；
(2)、将RRI信号x(t)作为原始信号，求出原始信号的所有局部极值点，然后将原始信号的所有极大值点和所有极小值点采用样条曲线连起来分别形成上包络线emax和下包络线emin，得到上、下包络线的包络均值信号m(t)＝(emax+emin)/2；
(3)、将原始信号x(t)减去包络均值信号m(t)，得到h(t)＝x(t)-m(t)；然后判断h(t)是否满足极值模态函数的判定条件，如果不满足，将h(t)作为原始信号返回至步骤(3)，直到hk(t)满足极值模态函数的判定条件，则记c1(t)＝hk(t)，作为第一个极值模态函数分量；
(4)、将原始信号x(t)减去第一个极值模态函数分量c1(t)，得到余量r1(t)＝x(t)-c1(t)，然后判断hk(t)是否满足停止准则，如果不满足，将r1(t)作为新的原始序列x(t)，返回至步骤(2)和(3)；如果hk(t)满足停止准则但n＜8时，返回步骤(1)重新获取原始信号；如果hk(t)满足停止准则且n≥8时，得到第2、3、…、n个极值模态函数分量及余量rn(t)，于是将原始信号x(t)分解为n个极值模态函数分量和一个余量，即
(5)、对极值模态函数分量ci(t)，i＝1，2，…，n，进行频谱分析得到各极值模态函数分量的中心频率；
(6)、将原始信号x(t)分解得的n个极值模态函数分量，代表了原始信号不同频段的分量，然后计算其各个分量的能量
Ei＝∫|ci(t)|2dt,i＝1，2，…，n
将每一个能量值归一化，得到归一化的能量分布向量
pi＝Ei/E,i＝1,2,…,n
其中，第一个分量p1表示最高频段的能量，代表了信号在最高频段范围内能量分布的比例，最后一个分量pn表示信号在最低频段范围内能量分布的比例；根据归一化的能量分布向量绘制归一化能量分布图，其中横坐标表示分量层次，纵坐标表示归一化的能量分布向量值，曲线表示平均值，误差棒表示标准差；
(7)、选取第二个分量p2至第七个分量p7，计算能量差异值EDV，EDV＝(p2+p3+p4)-(p5+p6+p7)，当EDV≤第一标准值时，则判定该RRI信号为正常心率变异性信号，当第一标准值＜EDV＜第二标准值时，则判定该RRI信号为疑似异常心率变异性信号，当EDV≥第二标准值时，则判定该RRI信号为异常心率变异性信号。
2.根据权利要求1的一种基于极值能量分解法的心率变异性信号分析方法，其特征在于：所述原始信号x(t)所需最少数据量N＝2n+1，其中n为分解出的极值模态函数分量的数量。
3.根据权利要求1的一种基于极值能量分解法的心率变异性信号分析方法，其特征在于：所述步骤(1)中去噪预处理的具体方法为：将ECG信号经过40Hz零相位FIR低通滤波器滤波消除高频噪声，然后经过中值滤波器去除基线漂移。
4.根据权利要求1的一种基于极值能量分解法的心率变异性信号分析方法，其特征在于：所述步骤(3)中极值模态函数的判定条件为：(a)、在整个数据序列中，极值点的数量与过零点的数量相等或者相差一个；(b)、在任意时刻，上下包络线对于时间轴对称。
5.根据权利要求1的一种基于极值能量分解法的心率变异性信号分析方法，其特征在于：所述步骤(4)中hk(t)满足停止准则公式为：
ε表示筛选门限，取0.2～0.3之间。
6.根据权利要求1的一种基于极值能量分解法的心率变异性信号分析方法，其特征在于：所述步骤(7)中的第一标准值为-0.15，第二标准值为0.08。
CN201910841605.6A 2019-09-06 2019-09-06 一种基于极值能量分解法的心率变异性信号分析方法 Active CN110464337B (zh)

Priority Applications (3)

Application Number Priority Date Filing Date Title
CN201910841605.6A CN110464337B (zh) 2019-09-06 2019-09-06 一种基于极值能量分解法的心率变异性信号分析方法
PCT/CN2019/121416 WO2021042591A1 (zh) 2019-09-06 2019-11-28 一种基于极值能量分解法的心率变异性信号分析方法
AU2019101756A AU2019101756A4 (en) 2019-09-06 2019-11-28 Method for analyzing heart rate variability signal based on extremum energy decomposition method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910841605.6A CN110464337B (zh) 2019-09-06 2019-09-06 一种基于极值能量分解法的心率变异性信号分析方法

Publications (2)

Publication Number Publication Date
CN110464337A true CN110464337A (zh) 2019-11-19
CN110464337B CN110464337B (zh) 2020-09-01

Family

ID=68515060

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910841605.6A Active CN110464337B (zh) 2019-09-06 2019-09-06 一种基于极值能量分解法的心率变异性信号分析方法

Country Status (3)

Country Link
CN (1) CN110464337B (zh)
AU (1) AU2019101756A4 (zh)
WO (1) WO2021042591A1 (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2021042591A1 (zh) * 2019-09-06 2021-03-11 江苏华康信息技术有限公司 一种基于极值能量分解法的心率变异性信号分析方法
CN116548928A (zh) * 2023-07-11 2023-08-08 西安浩阳志德医疗科技有限公司 一种基于互联网的护理服务***
CN117503153A (zh) * 2024-01-05 2024-02-06 北华大学 基于人工智能的患者术后康复评价方法

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116738153B (zh) * 2023-08-15 2023-10-17 中国科学院东北地理与农业生态研究所 基于光谱分析的有机肥利用效果评价方法
CN117390380B (zh) * 2023-12-12 2024-02-13 泰安金冠宏食品科技有限公司 一种油渣分离***中的数据分析方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080269628A1 (en) * 2007-04-25 2008-10-30 Siemens Medical Solutions Usa, Inc. Denoising and Artifact Rejection for Cardiac Signal in a Sensis System
CN101496716A (zh) * 2009-02-26 2009-08-05 周洪建 利用ecg信号检测睡眠呼吸暂停的测量方法
CN101642369A (zh) * 2008-08-04 2010-02-10 南京大学 自主神经功能生物反馈方法和***
CN108814579A (zh) * 2018-04-16 2018-11-16 西安交通大学 一种基于emd分解的心电、呼吸联合计算心率变异性的方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110464337B (zh) * 2019-09-06 2020-09-01 江苏华康信息技术有限公司 一种基于极值能量分解法的心率变异性信号分析方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080269628A1 (en) * 2007-04-25 2008-10-30 Siemens Medical Solutions Usa, Inc. Denoising and Artifact Rejection for Cardiac Signal in a Sensis System
CN101642369A (zh) * 2008-08-04 2010-02-10 南京大学 自主神经功能生物反馈方法和***
CN101496716A (zh) * 2009-02-26 2009-08-05 周洪建 利用ecg信号检测睡眠呼吸暂停的测量方法
CN108814579A (zh) * 2018-04-16 2018-11-16 西安交通大学 一种基于emd分解的心电、呼吸联合计算心率变异性的方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
庞春颖等: "一种基于改进经验模态分解的癫痫脑电识别新方法", 《中国生物医学工程学报》 *
董红生等: "基于Hilbert谱的心率变异信号时频分析方法", 《仪器仪表学报》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2021042591A1 (zh) * 2019-09-06 2021-03-11 江苏华康信息技术有限公司 一种基于极值能量分解法的心率变异性信号分析方法
CN116548928A (zh) * 2023-07-11 2023-08-08 西安浩阳志德医疗科技有限公司 一种基于互联网的护理服务***
CN116548928B (zh) * 2023-07-11 2023-09-08 西安浩阳志德医疗科技有限公司 一种基于互联网的护理服务***
CN117503153A (zh) * 2024-01-05 2024-02-06 北华大学 基于人工智能的患者术后康复评价方法
CN117503153B (zh) * 2024-01-05 2024-03-15 北华大学 基于人工智能的患者术后康复评价方法

Also Published As

Publication number Publication date
AU2019101756A4 (en) 2020-11-12
WO2021042591A1 (zh) 2021-03-11
CN110464337B (zh) 2020-09-01

Similar Documents

Publication Publication Date Title
CN110464337A (zh) 一种基于极值能量分解法的心率变异性信号分析方法
AU2019101755A4 (en) Method for quantitatively analyzing electrocardiogram signal based on extremum energy decomposition method
Acharya et al. Comprehensive analysis of cardiac health using heart rate signals
Li et al. On an automatic delineator for arterial blood pressure waveforms
CN105286815B (zh) 一种基于波形时域特征的脉搏波信号特征点检测方法
Podziemski et al. Fetal heart rate discovery: algorithm for detection of fetal heart rate from noisy, noninvasive fetal ECG recordings
CN108471961A (zh) 一种生理参数的计算方法及相应的医疗设备
Singh et al. Effects of RR segment duration on HRV spectrum estimation
CN109938719A (zh) 一种基于生理参数的驾驶员疲劳检测方法
CN110558959B (zh) 一种基于极值能量分解法的冥想训练的hrv信号分析方法
CN112656393A (zh) 一种心率变异性的检测方法及***
CA2888513A1 (en) Device and method for detecting and reporting of a stress condition of a person
Faust et al. Heart rate variability analysis for different age and gender
CN110558974B (zh) 一种基于极值能量分解法的心电图信号分析方法
Garcia-Gonzalez et al. Errors in the estimation of approximate entropy and other recurrence-plot-derived indices due to the finite resolution of RR time series
Nenova et al. An automated algorithm for fast pulse wave detection
Maciejewski et al. ECG parameter extraction and classification in noisy signals
Chen et al. A fast ECG diagnosis using frequency-based compressive neural network
Daqrouq et al. Arrhythmia detection using wavelet transform
Kovalchuk et al. Multifractal analysis of cardiac series and predictors of sudden cardiac death
Fahruzi et al. Screening of non-overlapping apnea and non-apnea from single lead ECG-apnea recordings using time-frequency approach
Daoud et al. HRV spectral estimation based on constrained Gaussian modeling in the nonstationary case
Schack et al. Adaptive methods of trend detection and their application in analysing biosignals
Patil et al. Discrimination between atrial fibrillation (AF) & normal sinus rhythm (NSR) using linear parameters
Castiglioni et al. How to check steady-state condition from cardiovascular time series

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB03 Change of inventor or designer information

Inventor after: Ning Xinbao

Inventor after: Zeng Peng

Inventor after: Zhou Zuojian

Inventor after: Jiang Xiaodong

Inventor after: Liu Guo

Inventor after: Wang Binbin

Inventor after: Wang Hua

Inventor before: Zhou Zuojian

Inventor before: Ning Xinbao

Inventor before: Wang Binbin

Inventor before: Jiang Xiaodong

Inventor before: Wang Hua

CB03 Change of inventor or designer information
GR01 Patent grant
GR01 Patent grant