CN114145749B - 基于优化模型的ecg信号有限新息率采样方法 - Google Patents
基于优化模型的ecg信号有限新息率采样方法 Download PDFInfo
- Publication number
- CN114145749B CN114145749B CN202111307289.8A CN202111307289A CN114145749B CN 114145749 B CN114145749 B CN 114145749B CN 202111307289 A CN202111307289 A CN 202111307289A CN 114145749 B CN114145749 B CN 114145749B
- Authority
- CN
- China
- Prior art keywords
- sampling
- signal
- ecg signal
- optimization
- particle
- 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.)
- Active
Links
- 238000005070 sampling Methods 0.000 title claims abstract description 164
- 238000005457 optimization Methods 0.000 title claims abstract description 69
- 238000000034 method Methods 0.000 title claims abstract description 57
- 239000002245 particle Substances 0.000 claims abstract description 91
- 238000001914 filtration Methods 0.000 claims abstract description 36
- 238000012952 Resampling Methods 0.000 claims description 18
- 230000008569 process Effects 0.000 claims description 12
- 230000004044 response Effects 0.000 claims description 6
- 238000001228 spectrum Methods 0.000 claims description 6
- 230000003595 spectral effect Effects 0.000 claims description 4
- 238000009826 distribution Methods 0.000 claims description 3
- 238000012804 iterative process Methods 0.000 claims description 3
- 238000005259 measurement Methods 0.000 claims description 3
- 238000009827 uniform distribution Methods 0.000 claims description 3
- 230000009467 reduction Effects 0.000 abstract description 5
- 230000008859 change Effects 0.000 description 6
- 238000012544 monitoring process Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 3
- 238000001208 nuclear magnetic resonance pulse sequence Methods 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 230000007774 longterm Effects 0.000 description 2
- 238000012806 monitoring device Methods 0.000 description 2
- 230000002028 premature Effects 0.000 description 2
- 230000002411 adverse Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005265 energy consumption Methods 0.000 description 1
- 230000036541 health Effects 0.000 description 1
- 230000005802 health problem Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
Classifications
-
- 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/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
-
- 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/7253—Details of waveform analysis characterised by using transforms
- A61B5/7257—Details of waveform analysis characterised by using transforms using Fourier transforms
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Molecular Biology (AREA)
- General Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Medical Informatics (AREA)
- Veterinary Medicine (AREA)
- Surgery (AREA)
- Animal Behavior & Ethology (AREA)
- Pathology (AREA)
- Public Health (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physiology (AREA)
- Psychiatry (AREA)
- Signal Processing (AREA)
- Cardiology (AREA)
- Mathematical Physics (AREA)
- Complex Calculations (AREA)
Abstract
一种基于优化模型的ECG信号有限新息率采样方法,适用于ECG信号FRI采样及信号重构,所提方法将ECG信号建模为一阶微分VPW模型,以利用FRI采样理论对其进行欠采样;然后使用双通道采样结构进行采样,采样通道一是一个FRI采样通道,用于获取信号的少量样本以进行参数估计;采样通道二为低速采样通道,用于采集少量的信号时域样本,以用于对ECG信号进行优化重构。为了提高ECG信号的重构精度,采用于粒子滤波的参数优化算法,以提高ECG信号参数估计的精度,从而提高信号重构的质量。为了提高参数估计的效率和优化的成功率,在基于粒子滤波的参数估计算法中引入了块坐标下降的思想,提升了优化的成功率和优化的效率。
Description
技术领域
本发明涉及心电图ECG信号处理技术领域,具体涉及一种ECG信号的欠采样方法。
背景技术
心电图(Electrocardiogram,ECG)信号是记录人体心跳活动的电信号,广泛应用于医学的人体监测领域。ECG信号是一种脉冲序列信号,如图1所示,可以视为由P,Q,R,S,T等脉冲组合而成。ECG信号是监测人体心脏健康的重要手段,对人体ECG信号进行长期监测,有利于及时有效地发现人体心脏的健康问题。近年来,ECG信号的监测设备出现了便捷化,小型化的发展趋势,因此对ECG信号进行欠采样以降低ECG信号监测设备长期监测产生的能量消耗具有非常重要的实际应用价值。
有限新息率(Finite Rate of Innovation,FRI)采样理论是一种针对参数稀疏信号的欠采样方法,可以对信号以大于等于信号的新息率的采样速率进行欠采样,并能够对信号进行精准重构。而信号的新息率通常远小于信号的奈奎斯特采样速率。ECG信号是一种脉冲序列信号,可以被建模为参数稀信号,因而可以根据FRI采样理论进行欠采样,以降低采样速率,从而减少采样的数据量。
ECG信号是一种脉冲序列信号,而脉冲序列信号是一种FRI信号。因此,可以根据FRI采样原理对ECG信号进行欠采样。首先,需要对ECG信号进行建模,将其建模为FRI信号模型。ECG信号的脉冲波形具有非对称的特点,因此更适合利用非对称脉冲序列模型对其进行FRI建模。目前常用的非对称脉冲序列模型主要有:Nagesh、Seelamantula等人以高斯脉冲及其导数组合为基函数构造的AP-FRI模型,Nagesh、Mulleti等人以高斯脉冲及其分数阶希尔伯特变换组合作为基函数的FrHT-FRI模型,Baechler等人以脉宽可变的VPW脉冲作为基函数的VPW-FRI模型。AP-FRI,FrHT-FRI信号模型的基函数具有脉冲固定的特点,与ECG信号的模型匹配误差较大,不适合对ECG信号进行建模。VPW-FRI模型的基函数具有脉宽可变的特点,能够适应ECG信号的脉宽变化的特点,但是由于VPW-FRI模型的基函数的变化速度不够迅速,很难适应ECG快速变化的波形,模型匹配误差仍然较大。一阶微分VPW模型是VPW-FRI模型的扩展,一阶微分VPW模型的基函数由VPW脉冲及其微分脉冲组合而成,能够适应ECG信号波形在波峰波谷快速变化的特点,从而减低模型匹配误差,提高信号重构的精度。
将ECG信号建模为FRI信号模型,并经过采样滤波之后,可以以大于等于ECG信号新息率的采样速率进行均匀采样,之后根据采样信号进行参数估计以重构ECG信号。但是,FRI信号模型与ECG信号并不是完美匹配的,一般都存在模型匹配误差。模型匹配误差会降低信号重构的精度。ECG信号的参数估计是一个频谱估计问题,参数估计的精度会影响ECG信号重构的质量。ECG信号FRI建模的模型不匹配以及采样过程中引入的噪声干扰都会对参数估计的精度以及信号重构的精度产生不利的影响,进而导致ECG信号重构质量的降低。因此,引入优化思想对ECG信号的参数估计过程进行优化,以提高信号重构的质量是一种非常有效的方法。综上所述,模型不匹配下的ECG信号FRI采样方法的研究是有意义的。
发明内容
为了克服已有技术的不足,针对ECG信号的FRI采样存在的模型不匹配和重构精度较低的问题,本发明提出了一种基于优化模型的ECG信号有限新息率采样方法,适用于ECG信号FRI采样及信号重构,所提方法将ECG信号建模为一阶微分VPW模型,以利用FRI采样理论对其进行欠采样。然后使用双通道采样结构进行采样,采样通道一是一个FRI采样通道,用于获取信号的少量样本以进行参数估计。采样通道二为低速采样通道,用于采集少量的信号时域样本,以用于对ECG信号进行优化重构。为了提高ECG信号的重构精度,本发明提出了一种基于粒子滤波的参数优化算法,以提高ECG信号参数估计的精度,从而提高信号重构的质量。为了提高参数估计的效率和优化的成功率,在基于粒子滤波的参数估计算法中引入了块坐标下降(block coordinate descent,BCD)的思想,提升了优化的成功率和优化的效率。
本发明解决其技术问题所采用的技术方案是:
一种基于优化模型的ECG信号有限新息率采样方法,包括以下步骤:
步骤1:根据FRI采样理论,将ECG信号建模为一阶微分VPW模型信号,并忽略模型匹配误差:
其中,s(t)是ECG信号,t表示时间,T是ECG信号s(t)的观测时长。K是基函数一阶微分VPW脉冲fk(t)的个数,k∈[1,K],K是一个正整数。是VPW脉冲的r阶微分脉冲,分别是第k个基函数的r阶微分脉冲/>的对称部分脉冲幅值,非对称部分脉冲幅值,脉宽和时延参数,一阶微分VPW模型的傅里叶系数如下所示:
其中S[m]是ECG信号的傅里叶级数系数;
步骤2:利用双通道采样结构对ECG信号进行低速均匀采样,过程如下:
采样通道一:对ECG信号s(t)进行采样核滤波,采样核为sinc采样核,得到滤波之后的信号g(t):
其中,Y(f)是滤波信号的频谱,S(f)是ECG信号的频谱,f表示频率,B是sinc采样核的带宽。h(t)是sinc采样核的冲激响应:
其中,H(f)是sinc采样核的频率响应,B是sinc采样核的带宽,通常大于等于信号的新息率ρ,滤波器的截止频率为fc=B/2,因此,通道一的采样速率fs1≥2fc。之后,以采样速率fs1对ECG信号进行低速均匀采样,得到FRI采样之后的信号y[n],n∈[0,N-1]是采样点,是采样通道一的采样点数;
采样通道二:对ECG信号进行低速均匀采样,采样速率为fs2=fs1,以获取ECG信号的少量时域样本s[n′],n′∈[0,N′-1]是采样通道二的采样点,是采样通道二的采样点数。采样通道二的样本用于辅助信号参数估计的优化;尽管当fs2小于ECG信号的奈奎斯特采样速率时会造成采样信号的频谱混叠,但是采样信号在时域上却不会受到影响;
步骤3:利用基于粒子滤波的参数优化算法对ECG信号的参数估计进行优化,过程如下:
步骤3.1:计算采样通道一获取的样本y[n]的傅里叶系数Y[m],m∈[0,M],
步骤3.2:开始迭代优化,由于ECG信号通常包含5个脉冲,即P,Q,R,S,T脉冲,因此,为了减少遍历的次数,令微分VPW模型的基函数个数取值范围为K=4:7;
步骤3.3:根据M+1个采样信号的傅里叶系数,利用零化滤波器算法估计出信号的所有参数
步骤3.4:利用估计的参数重构ECG信号,得到重构信号/>
步骤3.5:构造待优化的目标函数:令t=n′Ts2=n′/fs2,得到信号离散化的重构信号/>利用通道二获取的少量时域样本s[n]构造待优化的目标函数:
其中,dev是模型匹配误差,是待优化的变量集合;
步骤3.6:利用BCD方法的思想对待优化的变量集合X={xj},j∈[1,4]进行优化,其中 将估计出的ECG信号的参数作为待优化的目标函数自变量的初始值/>依次对xj进行优化,并利用优化得到的结果/>更新待优化的变量集合/>
步骤3.7:对变量集合X={xj},j∈[1,4]进行迭代优化:itra=1:10;
步骤3.8:依次优化xj:迭代次数j=1:4;
步骤3.9:开始对变量集合xj进行迭代优化:迭代次数n=1:5;
步骤3.10:初始化:将公式(5)所示的目标函数的向最优解的变化过程视为一个动态时变系统,系统的状态则是目标函数的解,利用粒子滤波的基本原理估计出系统状态的最小方差估计,即目标函数的最优解,将系统状态初始化为零化滤波器方法估计得到的参数
步骤3.11:粒子滤波迭代优化:i=1:I;
步骤3.12:采样:在系统状态的解空间随机采样Np个粒子,对于第m′∈[1,Np]个粒子,其采样满足:使用均匀分布/>代替/>ci=Λ/1.01i是一个随着迭代过程的进行而逐渐缩小的参数,Λ是待优化的变量X的取值范围,随着迭代的进行而不断地缩小采样区间,以使系统的状态粒子最终收敛于目标函数的全局最优解;
步骤3.13:更新全局最优解:将公式六所示的待优化的目标函数作为系统的测量函数,计算当前迭代次数下每个粒子的测量值之后,找到并保存测量值的最小值及其对应的粒子:
步骤3.14:更新粒子权重:将粒子的测量值/>与状态粒子xi-1的测量值yi-1进行比较:若是粒子的测量值/>大于状态粒子的测量值yi-1,则该粒子的权重设置为0;若是粒子的测量值/>小于等于状态粒子xi-1的测量值yi-1,则该粒子的权重由/>与yi-1的距离决定。具体的实现方法是将粒子的权重/>视为服从正态分布N(yi-1,σ2),其中σ2是样本的方差,则每个粒子的权重经过更新后由下列式子得到:
之后,对每个粒子的权值进行归一化处理:
步骤3.15:重采样:随着迭代的进行,有效粒子可能会出现大规模减少的现象,从而导致优化算法过早收敛于某一极值点,为了避免这种现象,当有效粒子的数目Neff小于一个阈值Nth=2N′/3时,必须启动重采样操作,有效粒子的数目通过以下式子计算:
当Neff<Nth时,启动重采样操作,独立重采样的方式进行重采样,经过重采样之后,权重较少的粒子则会减少,而权重较大的粒子则会保留;
步骤3.16:估计系统状态:经过对粒子的权重进行更新之后,可以估计出系统在第i个时刻的状态:
步骤3.17:判断:如果i≤I,则返回步骤3.11迭代优化;如果i>I,停止迭代,输出基于粒子滤波的参数估计算法对变量集合xj的优化结果及其对应的/>
步骤3.18:更新xj:当yj[n]趋于不变或者n>5时,停止迭代,否则返回步骤3.9继续迭代优化;
步骤3.19:如果j>4,则所有变量优化完成,变量X={xj},j∈[1,4]优化完成,对应的最优值y(itra)=y4(end),否则返回步骤3.8优化下一个变量集合;
步骤3.20:如果y(itra)的值趋于不变或者对变量X的迭代优化次数itra>10,则停止迭代,得到优化完成的变量集合X={xj},j∈[1,4],及其对应的最优值y(end),否则返回步骤3.7迭代优化;
步骤3.21:迭代结束之后,输出当前迭代次数K对应的最优值yopt(K-3)及其最优解Xopt(K-3,:):
步骤3.22:判断:如果K≤7,则返回步骤3.2迭代优化;如果K>7,停止迭代;
步骤4:从集合yopt找到模型匹配误差的最小值,以及对应的最优解:
步骤5:将基于粒子滤波的参数估计算法得到的最优基函数个数Kopt以及优化参数集合代入到一阶微分VPW模型当中,重构出ECG信号s′(t)。
进一步,所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.6将零化滤波器方法估计出的ECG信号的参数作为待优化的目标函数自变量的初始值。
所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.10将系统状态初始化为零化滤波器方法估计得到的参数。
本发明中,基于优化模型的ECG信号有限新息率采样方法,其系统结构如图2所示。首先将心电图ECG信号建模为一阶微分VPW模型,然后通过一个双通道采样框架对其进行采样。采样通道一是一个FRI采样通道,sinc采样核用于对ECG信号进行采样核滤波,然后由ADC(Analog to Digital Converter)进行低速采样,获取少量的信号样本。通道二是一个低速采样通道,由ADC低速采样少量信号的时域样本。根据两个通道的样本,本发明方法提出一种基于粒子滤波的参数优化算法来对信号参数进行优化估计,并最终重构ECG信号。
将ECG信号建模为一阶微分VPW模型,并针对ECG信号的欠采样与优化重构提出了一种双通道采样方法。采样通道一对ECG信号进行FRI采样,采样通道二对ECG信号进行低速均匀采样。为了提高信号参数估计和信号重构的精度,提出了一种基于粒子滤波的ECG信号参数优化算法来估计ECG信号的参数,从而精准重构出原信号。本发明利用零化滤波器算法估计信号的参数作为初始值,并利用采样通道二的时域样本辅助构建了待优化的目标函数。之后,将目标函数最优解的变化过程视为一个动态时变系统,将目标函数的最优解作为系统的状态,将优化问题转化为状态估计问题,从而可以利用粒子滤波的思想求解优化问题。最终,在优化完成之后,能够得到一个模型匹配误差最小化的重构信号。
本发明的有益效果主要表现在:提高了ECG信号重构的精度。
附图说明
图1是单个心跳周期的ECG信号波形图。
图2是一种基于优化模型的ECG信号有限新息率采样方法的系统机构图。
图3是VPW-FRI方案的重构结果图。
图4是本发明方法的重构结果图。
图5是不同ECG信号的重构结果图,其中,(a)是ECG信号1的重构结果,(b)是ECG信号2的重构结果,(c)是ECG信号3的重构结果,(d)是ECG信号6的重构结果。
具体实施方式
下面结合附图对本发明作进一步描述。
参照图1~图5,一种基于优化模型的ECG信号有限新息率采样方法,包括以下步骤:
步骤1:根据FRI采样理论,将ECG信号建模为一阶微分VPW模型信号,并忽略模型匹配误差:
其中,s(t)是ECG信号,t表示时间,T是ECG信号s(t)的观测时长。K是基函数一阶微分VPW脉冲fk(t)的个数,k∈[1,K],K是一个正整数,是VPW脉冲的r阶微分脉冲,分别是第k个基函数的r阶微分脉冲/>的对称部分脉冲幅值,非对称部分脉冲幅值,脉宽和时延参数,一阶微分VPW模型的傅里叶系数如下所示:
其中S[m]是ECG信号的傅里叶级数系数;
步骤2:利用双通道采样结构对ECG信号进行低速均匀采样,过程如下:
采样通道一:对ECG信号s(t)进行采样核滤波,采样核为sinc采样核,得到滤波之后的信号g(t):
其中,Y(f)是滤波信号的频谱,S(f)是ECG信号的频谱,f表示频率,B是sinc采样核的带宽。h(t)是sinc采样核的冲激响应:
其中,H(f)是sinc采样核的频率响应,B是sinc采样核的带宽,通常大于等于信号的新息率ρ,滤波器的截止频率为fc=B/2,因此,通道一的采样速率fs1≥2fc。之后,以采样速率fs1对ECG信号进行低速均匀采样,得到FRI采样之后的信号y[n],n∈[0,N-1]是采样点,是采样通道一的采样点数;
采样通道二:对ECG信号进行低速均匀采样,采样速率为fs2=fs1,以获取ECG信号的少量时域样本s[n′],n′∈[0,N′-1]是采样通道二的采样点,是采样通道二的采样点数。采样通道二的样本用于辅助信号参数估计的优化;尽管当fs2小于ECG信号的奈奎斯特采样速率时会造成采样信号的频谱混叠,但是采样信号在时域上却不会受到影响;
步骤3:利用基于粒子滤波的参数优化算法对ECG信号的参数估计进行优化,过程如下:
步骤3.1:计算采样通道一获取的样本y[n]的傅里叶系数Y[m],m∈[0,M],
步骤3.2:开始迭代优化,由于ECG信号通常包含5个脉冲,即P,Q,R,S,T脉冲,因此,为了减少遍历的次数,令微分VPW模型的基函数个数取值范围为K=4:7;
步骤3.3:根据M+1个采样信号的傅里叶系数,利用零化滤波器算法估计出信号的所有参数
步骤3.4:利用估计的参数重构ECG信号,得到重构信号/>
步骤3.5:构造待优化的目标函数:令t=n′Ts2=n′/fs2,得到信号离散化的重构信号/>利用通道二获取的少量时域样本s[n]构造待优化的目标函数:
其中,dev是模型匹配误差,是待优化的变量集合;
步骤3.6:利用BCD方法的思想对待优化的变量集合X={xj},j∈[1,4]进行优化,其中 将估计出的ECG信号的参数作为待优化的目标函数自变量的初始值/>依次对xj进行优化,并利用优化得到的结果/>更新待优化的变量集合/>
步骤3.7:对变量集合X={xj},j∈[1,4]进行迭代优化:itra=1:10;
步骤3.8:依次优化xj:迭代次数j=1:4;
步骤3.9:开始对变量集合xj进行迭代优化:迭代次数n=1:5;
步骤3.10:初始化:将公式(5)所示的目标函数的向最优解的变化过程视为一个动态时变系统,系统的状态则是目标函数的解,利用粒子滤波的基本原理估计出系统状态的最小方差估计,即目标函数的最优解,将系统状态初始化为零化滤波器方法估计得到的参数
步骤3.11:粒子滤波迭代优化:i=1:I;
步骤3.12:采样:在系统状态的解空间随机采样Np个粒子,对于第m′∈[1,Np]个粒子,其采样满足:使用均匀分布/>代替/>ci=Λ/1.01i是一个随着迭代过程的进行而逐渐缩小的参数,Λ是待优化的变量X的取值范围,随着迭代的进行而不断地缩小采样区间,以使系统的状态粒子最终收敛于目标函数的全局最优解;
步骤3.13:更新全局最优解:将公式六所示的待优化的目标函数作为系统的测量函数,计算当前迭代次数下每个粒子的测量值之后,找到并保存测量值的最小值及其对应的粒子:
步骤3.14:更新粒子权重:将粒子的测量值/>与状态粒子xi-1的测量值yi-1进行比较:若是粒子的测量值/>大于状态粒子的测量值yi-1,则该粒子的权重设置为0;若是粒子的测量值/>小于等于状态粒子xi-1的测量值yi-1,则该粒子的权重由/>与yi-1的距离决定。具体的实现方法是将粒子的权重/>视为服从正态分布N(yi-1,σ2),其中σ2是样本的方差,则每个粒子的权重经过更新后由下列式子得到:
之后,对每个粒子的权值进行归一化处理:
步骤3.15:重采样:随着迭代的进行,有效粒子可能会出现大规模减少的现象,从而导致优化算法过早收敛于某一极值点,为了避免这种现象,当有效粒子的数目Neff小于一个阈值Nth=2N′/3时,必须启动重采样操作,有效粒子的数目通过以下式子计算:
当Neff<Nth时,启动重采样操作,独立重采样的方式进行重采样,经过重采样之后,权重较少的粒子则会减少,而权重较大的粒子则会保留;
步骤3.16:估计系统状态:经过对粒子的权重进行更新之后,可以估计出系统在第i个时刻的状态:
步骤3.17:判断:如果i≤I,则返回步骤3.11迭代优化;如果i>I,停止迭代,输出基于粒子滤波的参数估计算法对变量集合xj的优化结果及其对应的/>
步骤3.18:更新xj:当yj[n]趋于不变或者n>5时,停止迭代,否则返回步骤3.9继续迭代优化;
步骤3.19:如果j>4,则所有变量优化完成,变量X={xj},j∈[1,4]优化完成,对应的最优值y(itra)=y4(end),否则返回步骤3.8优化下一个变量集合;
步骤3.20:如果y(itra)的值趋于不变或者对变量X的迭代优化次数itra>10,则停止迭代,得到优化完成的变量集合X={xj},j∈[1,4],及其对应的最优值y(end),否则返回步骤3.7迭代优化;
步骤3.21:迭代结束之后,输出当前迭代次数K对应的最优值yopt(K-3)及其最优解Xopt(K-3,:):
步骤3.22:判断:如果K≤7,则返回步骤3.2迭代优化;如果K>7,停止迭代;
步骤4:从集合yopt找到模型匹配误差的最小值,以及对应的最优解:
步骤5:将基于粒子滤波的参数估计算法得到的最优基函数个数Kopt以及优化参数集合代入到一阶微分VPW模型当中,重构出ECG信号s′(t)。
进一步,所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.6将零化滤波器方法估计出的ECG信号的参数作为待优化的目标函数自变量的初始值。
所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.10将系统状态初始化为零化滤波器方法估计得到的参数。
为了验证本发明的性能,将通过仿真实验来验证,并与VPW-FRI方案进行对比。VPW-FRI方案利用K个VPW脉冲对ECG信号进行建模,然后使用sinc采样核对ECG信号进行采样核滤波,并对滤波之后的ECG信号进行低速均匀采样。获得ECG信号的样本之后,计算ECG信号的傅里叶系数,并使用零化滤波器算法进行参数估计,并利用估计出的参数重构ECG信号。
本发明所提的方法将ECG信号建模为一阶微分VPW模型,以降低模型匹配误差。之后,使用双通道采样结构对ECG信号进行欠采样。其中,采样通道一与VPW-FRI方案的采样结构相同,先对ECG信号进行sinc采样核滤波,之后对滤波之后的ECG信号进行低速采样。采样通道二使用一个低速ADC对ECG信号进行均匀采样,获取ECG信号的少量时域样本。获取采样通道一和采样通道二获取的样本之后,利用本发明所提的基于粒子滤波的参数估计算法估计出ECG信号的参数,并完成信号重构。
为了衡量ECG信号的重构质量,使用重构信号的信噪比作为评价指标:
信噪比越大,意味着信号重构的效果越好,信号重构的质量越高。为了衡量ECG信号的欠采样程度,使用欠采样率作为评价指标:
系统的欠采样率越低,说明信号的欠采样程度更高。将sinc采样核的截止频率设为:fc=45Hz,VPW-FRI方案和采样通道一的采样速率设为:fs1=2fc=90Hz,采样通道二的采样速率设为:fs2=fs1=90Hz。
实验一:为了验证本发明所提方法的有效性,将实际的ECG信号分别建模为VPW模型和一阶微分VPW模型。之后,使用VPW-FRI方案和本发明所提方案对ECG信号进行欠采样并精准重构ECG信号。实验结果如图3,图4以及表1所示。图3展示了VPW-FRI方案的ECG信号重构的结果,VPW模型与ECG信号存在较大的模型匹配误差,导致ECG信号的重构效果不够理想,重构信号与原信号的模型匹配误差较大。图4展示了本发明所提方法的ECG信号重构的结果,一阶微分VPW模型的模型匹配误差明显更小,通过基于粒子滤波的参数估计算法对信号进行参数估计之后,重构信号与原始信号的模型匹配误差非常小,信号重构的质量非常高。从表1的数据可以看出,本发明所提方法的信号重构精度更高,但代价是ECG信号的欠采样程度降低。
表1。
为了进一步验证本发明所提方法的有效性,本次实验将利用不同的ECG信号来验证。实验结果如图5所示。本发明所提方法将ECG信号建模为一阶微分VPW模型,以降低模型匹配误差。双通道的采样结构虽然提高了采样速率,但是也提高了ECG信号重构的质量,使得模型匹配误差得到进一步的降低,大大提高了信号重构的精度。从图5的ECG信号重构结果中可以观察到,本发明所提方法的ECG信号重构效果非常好,重构信号与原始ECG信号具有非常高的匹配程度。
Claims (3)
1.一种基于优化模型的ECG信号有限新息率采样方法,其特征在于,所述方法包括以下步骤:
步骤1:根据FRI采样理论,将ECG信号建模为一阶微分VPW模型信号,并忽略模型匹配误差:
其中,s(t)是ECG信号,t表示时间,T是ECG信号s(t)的观测时长,K是基函数一阶微分VPW脉冲fk(t)的个数,k∈[1,K],K是一个正整数,是VPW脉冲的r阶微分脉冲,r=0,1,分别是第k个基函数的r阶微分脉冲/>的对称部分脉冲幅值,非对称部分脉冲幅值,脉宽和时延参数,一阶微分VPW模型的傅里叶系数如下所示:
其中S[m]是ECG信号的傅里叶级数系数;
步骤2:利用双通道采样结构对ECG信号进行低速均匀采样,过程如下:
采样通道一:对ECG信号s(t)进行采样核滤波,采样核为sinc采样核,得到滤波之后的信号g(t):
其中,Y(f)是滤波信号的频谱,S(f)是ECG信号的频谱,f表示频率,B是sinc采样核的带宽,h(t)是sinc采样核的冲激响应:
其中,H(f)是sinc采样核的频率响应,B是sinc采样核的带宽,通常大于等于信号的新息率ρ,滤波器的截止频率为fc=B/2,因此,通道一的采样速率fs1≥2fc,之后,以采样速率fs1对ECG信号进行低速均匀采样,得到FRI采样之后的信号y[n],n∈[0,N-1]是采样点,是采样通道一的采样点数;
采样通道二:对ECG信号进行低速均匀采样,采样速率为fs2=fs1,以获取ECG信号的少量时域样本s[n′],n′∈[0,N′-1]是采样通道二的采样点,是采样通道二的采样点数,采样通道二的样本用于辅助信号参数估计的优化;尽管当fs2小于ECG信号的奈奎斯特采样速率时会造成采样信号的频谱混叠,但是采样信号在时域上却不会受到影响;
步骤3:利用基于粒子滤波的参数优化算法对ECG信号的参数估计进行优化,过程如下:
步骤3.1:计算采样通道一获取的样本y[n]的傅里叶系数Y[m],m∈[0,M],
步骤3.2:开始迭代优化,由于ECG信号通常包含5个脉冲,即P,Q,R,S,T脉冲,因此,为了减少遍历的次数,令微分VPW模型的基函数个数取值范围为K=4:7;
步骤3.3:根据M+1个采样信号的傅里叶系数,利用零化滤波器算法估计出信号的所有参数
步骤3.4:利用估计的参数重构ECG信号,得到重构信号/>
步骤3.5:构造待优化的目标函数:令t=n′Ts2=n′/fs2,得到信号离散化的重构信号/>利用通道二获取的少量时域样本s[n]构造待优化的目标函数:
其中,dev是模型匹配误差,是待优化的变量集合;
步骤3.6:利用BCD方法的思想对待优化的变量集合X={xj},j∈[1,4]进行优化,其中 将估计出的ECG信号的参数作为待优化的目标函数自变量的初始值/>依次对xj进行优化,并利用优化得到的结果/>更新待优化的变量集合/>
步骤3.7:对变量集合X={xj},j∈[1,4]进行迭代优化:itra=1:10;
步骤3.8:依次优化xj:迭代次数j=1:4;
步骤3.9:开始对变量集合xj进行迭代优化:迭代次数n=1:5;
步骤3.10:初始化:将公式(5)所示的目标函数的向最优解的变化过程视为一个动态时变系统,系统的状态则是目标函数的解,利用粒子滤波的基本原理估计出系统状态的最小方差估计,即目标函数的最优解,将系统状态初始化为零化滤波器方法估计得到的参数
步骤3.11:粒子滤波迭代优化:i=1:I;
步骤3.12:采样:在系统状态的解空间随机采样Np个粒子,对于第m′∈[1,Np]个粒子,其采样满足:使用均匀分布/>代替/>ci=Λ/1.01i是一个随着迭代过程的进行而逐渐缩小的参数,Λ是待优化的变量X的取值范围,随着迭代的进行而不断地缩小采样区间,以使系统的状态粒子最终收敛于目标函数的全局最优解;
步骤3.13:更新全局最优解:将公式六所示的待优化的目标函数作为系统的测量函数,计算当前迭代次数下每个粒子的测量值之后,找到并保存测量值的最小值及其对应的粒子:
步骤3.14:更新粒子权重:将粒子的测量值/>与状态粒子xi-1的测量值yi-1进行比较:若是粒子的测量值/>大于状态粒子的测量值yi-1,则该粒子的权重设置为0;若是粒子的测量值/>小于等于状态粒子xi-1的测量值yi-1,则该粒子的权重由/>与yi-1的距离决定,实现方法是将粒子的权重/>视为服从正态分布N(yi-1,σ2),其中σ2是样本的方差,则每个粒子的权重经过更新后由下列式子得到:
之后,对每个粒子的权值进行归一化处理:
步骤3.15:重采样:随着迭代的进行,有效粒子可能会出现大规模减少的现象,从而导致优化算法过早收敛于某一极值点,为了避免这种现象,当有效粒子的数目Neff小于一个阈值Nth=2N′/3时,必须启动重采样操作,有效粒子的数目通过以下式子计算:
当Neff<Nth时,启动重采样操作,独立重采样的方式进行重采样,经过重采样之后,权重较少的粒子则会减少,而权重较大的粒子则会保留;
步骤3.16:估计系统状态:经过对粒子的权重进行更新之后,可以估计出系统在第i个时刻的状态:
步骤3.17:判断:如果i≤I,则返回步骤3.11迭代优化;如果i>I,停止迭代,输出基于粒子滤波的参数估计算法对变量集合xj的优化结果及其对应的/>
步骤3.18:更新当yj[n]趋于不变或者n>5时,停止迭代,否则返回步骤3.9继续迭代优化;
步骤3.19:如果j>4,则所有变量优化完成,变量X={xj},j∈[1,4]优化完成,对应的最优值y(itra)=y4(end),否则返回步骤3.8优化下一个变量集合;
步骤3.20:如果y(itra)的值趋于不变或者对变量X的迭代优化次数itra>10,则停止迭代,得到优化完成的变量集合X={xj},j∈[1,4],及其对应的最优值y(end),否则返回步骤3.7迭代优化;
步骤3.21:迭代结束之后,输出当前迭代次数K对应的最优值yopt(K-3)及其最优解Xopt(K-3,:):
步骤3.22:判断:如果K≤7,则返回步骤3.2迭代优化;如果K>7,停止迭代;
步骤4:从集合yopt找到模型匹配误差的最小值,以及对应的最优解:
步骤5:将基于粒子滤波的参数估计算法得到的最优基函数个数Kopt以及优化参数集合代入到一阶微分VPW模型当中,重构出ECG信号s′(t)。
2.如权利要求1所述的基于优化模型的ECG信号有限新息率采样方法,其特征在于,所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.6将零化滤波器方法估计出的ECG信号的参数作为待优化的目标函数自变量的初始值。
3.如权利要求1或2所述的基于优化模型的ECG信号有限新息率采样方法,其特征在于,所述步骤3中,所述的基于粒子滤波的参数估计算法的步骤3.10将系统状态初始化为零化滤波器方法估计得到的参数。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202111307289.8A CN114145749B (zh) | 2021-11-05 | 2021-11-05 | 基于优化模型的ecg信号有限新息率采样方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202111307289.8A CN114145749B (zh) | 2021-11-05 | 2021-11-05 | 基于优化模型的ecg信号有限新息率采样方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114145749A CN114145749A (zh) | 2022-03-08 |
CN114145749B true CN114145749B (zh) | 2024-04-05 |
Family
ID=80459009
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202111307289.8A Active CN114145749B (zh) | 2021-11-05 | 2021-11-05 | 基于优化模型的ecg信号有限新息率采样方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114145749B (zh) |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109889231A (zh) * | 2019-02-01 | 2019-06-14 | 哈尔滨工业大学 | 基于随机解调和有限新息率的脉冲串信号欠采样方法 |
CN110852888A (zh) * | 2019-10-18 | 2020-02-28 | 浙江工业大学 | 一种基于粒子滤波的证券投资组合优化方法 |
CN110944336A (zh) * | 2019-10-18 | 2020-03-31 | 浙江工业大学 | 一种基于有限新息率的时频谱感知方法 |
CN110992434A (zh) * | 2019-10-09 | 2020-04-10 | 浙江工业大学 | 一种基于有限新息率的emt图像重构方法 |
CN111820888A (zh) * | 2020-06-30 | 2020-10-27 | 浙江工业大学 | 基于一阶微分vpw模型的心电图ecg信号欠采样方法 |
CN112395546A (zh) * | 2020-11-27 | 2021-02-23 | 北京理工大学 | 一种基于线性正则域的有限新息率信号降采样与重构方法 |
-
2021
- 2021-11-05 CN CN202111307289.8A patent/CN114145749B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109889231A (zh) * | 2019-02-01 | 2019-06-14 | 哈尔滨工业大学 | 基于随机解调和有限新息率的脉冲串信号欠采样方法 |
CN110992434A (zh) * | 2019-10-09 | 2020-04-10 | 浙江工业大学 | 一种基于有限新息率的emt图像重构方法 |
CN110852888A (zh) * | 2019-10-18 | 2020-02-28 | 浙江工业大学 | 一种基于粒子滤波的证券投资组合优化方法 |
CN110944336A (zh) * | 2019-10-18 | 2020-03-31 | 浙江工业大学 | 一种基于有限新息率的时频谱感知方法 |
CN111820888A (zh) * | 2020-06-30 | 2020-10-27 | 浙江工业大学 | 基于一阶微分vpw模型的心电图ecg信号欠采样方法 |
CN112395546A (zh) * | 2020-11-27 | 2021-02-23 | 北京理工大学 | 一种基于线性正则域的有限新息率信号降采样与重构方法 |
Also Published As
Publication number | Publication date |
---|---|
CN114145749A (zh) | 2022-03-08 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN113820003B (zh) | 一种适用于桥梁振动监测的加速度实时重构动态位移方法 | |
CN108984474B (zh) | 一种非理想分段多项式信号的欠采样方法 | |
CN104101751B (zh) | 基于信息熵的数字存储示波器垂直分辨率提高方法 | |
CN107480391B (zh) | 基于数据驱动的近断层非平稳地震动模拟方法 | |
CN111814656A (zh) | 一种基于对抗生成网络的心电信号降噪方法 | |
CN113970420B (zh) | 一种基于深度学习的激波风洞测力信号频域分析方法 | |
CN111820888A (zh) | 基于一阶微分vpw模型的心电图ecg信号欠采样方法 | |
CN110490366A (zh) | 基于变分模态分解和迭代决策树的径流量预测方法 | |
CN104143031B (zh) | 一种基于小波多尺度分解的植被指数时序数据重构方法 | |
CN110313903B (zh) | 一种脉搏波频域特征参数提取方法及装置 | |
CN105137373A (zh) | 一种指数信号的去噪方法 | |
CN112362343A (zh) | 基于调频字典的齿轮箱变转速下分布型故障特征提取方法 | |
CN116049632B (zh) | 一种风电主轴轴承故障诊断方法、装置及应用 | |
CN112213560A (zh) | 一种基于z-adaline的高精度电网宽频信号测量方法 | |
CN114145749B (zh) | 基于优化模型的ecg信号有限新息率采样方法 | |
CN114813123A (zh) | 一种基于pso-vmd-mckd的滚动轴承微弱故障诊断方法 | |
CN101241145A (zh) | 一种针对示波器的波形采样和重建的方法 | |
CN108267311A (zh) | 一种基于张量分解的机械多维大数据处理方法 | |
CN112883318A (zh) | 相减策略的多频衰减信号参数估计算法 | |
CN108801296B (zh) | 基于误差模型迭代补偿的传感器频响函数计算方法 | |
CN106452693B (zh) | 一种基于双频点噪底能量分析的时钟相位抖动测量方法 | |
CN112129983B (zh) | 一种基于等时间间隔等效取样的波形恢复数据处理方法 | |
CN112468114B (zh) | 一种基于非理想sinc核的FRI采样系统及方法 | |
Cennamo et al. | Dynamic testing and diagnostics of digitizing signal analyzers | |
Jiang et al. | Time-frequency analysis—G (λ)-stationary processes |
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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |