CN110263760B - 一种基于病理自适应的心音全变分滤波方法 - Google Patents
一种基于病理自适应的心音全变分滤波方法 Download PDFInfo
- Publication number
- CN110263760B CN110263760B CN201910587766.7A CN201910587766A CN110263760B CN 110263760 B CN110263760 B CN 110263760B CN 201910587766 A CN201910587766 A CN 201910587766A CN 110263760 B CN110263760 B CN 110263760B
- Authority
- CN
- China
- Prior art keywords
- heart sound
- equation
- sound signal
- noise
- signal
- 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
- 238000001914 filtration Methods 0.000 title claims abstract description 53
- 238000000034 method Methods 0.000 title claims abstract description 48
- 230000007170 pathology Effects 0.000 title claims abstract description 9
- 230000005236 sound signal Effects 0.000 claims abstract description 68
- 230000001575 pathological effect Effects 0.000 claims abstract description 9
- 238000012216 screening Methods 0.000 claims abstract description 9
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 8
- 239000000654 additive Substances 0.000 claims abstract description 4
- 230000000996 additive effect Effects 0.000 claims abstract description 4
- AYFVYJQAPQTCCC-GBXIJSLDSA-N L-threonine Chemical compound C[C@@H](O)[C@H](N)C(O)=O AYFVYJQAPQTCCC-GBXIJSLDSA-N 0.000 claims description 15
- 238000001228 spectrum Methods 0.000 claims description 7
- 230000008569 process Effects 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 claims description 5
- 238000010276 construction Methods 0.000 claims description 5
- 238000005070 sampling Methods 0.000 claims description 5
- 230000011218 segmentation Effects 0.000 claims description 5
- 239000000126 substance Substances 0.000 claims description 4
- 238000005457 optimization Methods 0.000 claims description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims description 2
- 238000010252 digital analysis Methods 0.000 abstract description 2
- 238000012545 processing Methods 0.000 abstract description 2
- 230000034994 death Effects 0.000 description 8
- 231100000517 death Toxicity 0.000 description 8
- 208000024172 Cardiovascular disease Diseases 0.000 description 7
- 208000019622 heart disease Diseases 0.000 description 6
- 238000004458 analytical method Methods 0.000 description 5
- 208000037656 Respiratory Sounds Diseases 0.000 description 3
- 238000003745 diagnosis Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 238000002555 auscultation Methods 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 201000010099 disease Diseases 0.000 description 2
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 2
- 230000002526 effect on cardiovascular system Effects 0.000 description 2
- 230000036541 health Effects 0.000 description 2
- 230000002265 prevention Effects 0.000 description 2
- 208000035211 Heart Murmurs Diseases 0.000 description 1
- 208000020128 Mitral stenosis Diseases 0.000 description 1
- 206010028980 Neoplasm Diseases 0.000 description 1
- 208000031662 Noncommunicable disease Diseases 0.000 description 1
- 230000001133 acceleration Effects 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 230000032683 aging Effects 0.000 description 1
- 238000010009 beating Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 201000011510 cancer Diseases 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 206010012601 diabetes mellitus Diseases 0.000 description 1
- 230000008034 disappearance Effects 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 230000004064 dysfunction Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000004927 fusion Effects 0.000 description 1
- 210000004072 lung Anatomy 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 210000004115 mitral valve Anatomy 0.000 description 1
- 208000006887 mitral valve stenosis Diseases 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 231100000915 pathological change Toxicity 0.000 description 1
- 230000036285 pathological change Effects 0.000 description 1
- 230000005180 public health Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 208000023504 respiratory system disease Diseases 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
- G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
- G10L25/48—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use
- G10L25/51—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use for comparison or discrimination
- G10L25/66—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use for comparison or discrimination for extracting parameters related to health condition
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
- G06F2218/04—Denoising
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
Landscapes
- Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Signal Processing (AREA)
- Public Health (AREA)
- General Health & Medical Sciences (AREA)
- Computational Linguistics (AREA)
- Epidemiology (AREA)
- Audiology, Speech & Language Pathology (AREA)
- Human Computer Interaction (AREA)
- Physics & Mathematics (AREA)
- Acoustics & Sound (AREA)
- Multimedia (AREA)
- Measuring Pulse, Heart Rate, Blood Pressure Or Blood Flow (AREA)
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
Abstract
本发明涉及一种基于病理自适应的心音全变分滤波方法，属于医学计算生物特征识别技术领域；该方法首先通过对输入的原始心音信号添加加性高斯白噪声，筛选得到满足要求的加噪信号，然后对上述加噪信号进行改进的经验模态分解，筛选得到最大相似度时的本征模函数，并对原始心音信号进行分段；使用基于线性系数的追踪演化算法计算原始心音信号不同段的信息熵，并将信息熵作为权重值构建心音信号的l1和限定l1正则化全变分方程，并使用改进的Split‑Bregman算法分别求解，最后使用滑动平均法，得到融合两次滤波结果的心音信号。对比现有技术，本发明能够对原始的心音信号进行更为精确的去噪处理，并且，能够实现保留病理信息，为心音图的数字化分析提供基础。
Description
技术领域
本发明涉及一种基于病理自适应的心音全变分滤波方法，尤其涉及一种融合心音不同的病理特征和全变分滤波的心音滤波方法，属于医学计算生物特征识别技术领域。
背景技术
长久以来，心血管疾病因其高发病率和高死亡率严重威胁着人类的健康和生命。据世界卫生组织(WHO)估计，目前世界上每年有3600万人死于心血管疾病、糖尿病、呼吸系统疾病和恶性肿瘤等非传染性疾病，占全球死亡总数的2/3，到2020年，数目要攀升到4400万。
《中国心血管病报告(2018)》显示，2016年中国心血管病死亡率仍居疾病死亡构成的首位，高于肿瘤及其他疾病，每5例死亡中就有2例死于心血管病，2016年农村心血管病死亡率为309.33/10万，其中心脏病死亡率为151.18/10万；城市心血管病死亡率为265.11/10万，其中心脏病死亡率为138.70/10万。2016年，农村心血管病死亡占全部死因的比率为45.50％，城市心血管死亡占全部死因的比率为43.16％。随着中国社会老龄化步伐的加速，这些数据均在快速增长，心血管疾病已成为我国的重大公共卫生问题之一。所以，心脏类疾病的预防、诊断和治疗尤为重要，是当今医学人员的巨大挑战和研究热点。
心音图(PCG)是将心脏活动过程中产生的心音(Heart Sound)及杂音，通过换能器、放大电路、检波等方法，以振动波的形式记录下来。心音图有助于确定心音或杂音发生的时间，区分易于混淆的声音；鉴别某些听诊杂音相似而病变迥异的心内疾病，将心音与杂音变为可长期保存、分析的图形，作为判断病情变化，手术或药物疗效的客观指标。
心音图仪将心音及心脏杂音描记下来供临床分析，有助于心脏病的病因诊断，了解心脏杂音的产生机制，并用于心脏病尤其是瓣膜功能障碍的早期预防与诊断。
但是，传统的频域降噪方法(巴特沃斯滤波、切比雪夫滤波等)存在如下问题：(1)去噪时，未考虑病理性杂音，如二尖瓣狭窄时存在的亢进、拍击状杂音，由于频率和幅值较高，会被误判为高频、脉冲噪声；当S3段幅值较高时，会被误诊为测量的噪声被滤除；(2)心音在移动采集时，由于心肺距离较近，二尖瓣区易受到左肺影响，心音与肺音发生混叠，由于二者频率相近，基于频域分析的去噪算法无法有效将肺音去除。
并且目前由于心音图的智能分析的准确性和可靠性并不能与心脏病专家的听诊分析相提并论，所以它在实际临床应用中的比例并不大。
发明内容
本发明的目的在于进一步提高现有心音信号的滤波的准确度以及不同病理信号在去噪时的有效保留，提出了一种基于病理自适应的心音全变分滤波方法。
本发明的核心思想为：对原始心音信号进行S1、S2、S3和S4段的划分；同时使用基于线性系数的心音信号信息熵追踪演化计算方式，计算不同心音段的权重值；在得到不同心音段的权重值后，使用全变分滤波的方式，对不同段的心音信号进行滤波，得到最佳的定位结果。
一种基于病理自适应的心音全变分滤波方法，包括如下步骤：
步骤A.对原始心音信号HSorigin(n)添加加性高斯白噪声wgn(n)，基于加噪信号HSnoise(n)的频率分布情况X(ω)和标准偏差S，调整高斯白噪声的标准偏差σ，根据设定的阈值thre_S，筛选得到满足要求的HSnoise(n)；
作为优选，所述调整高斯白噪声的标准偏差σ，根据设定的阈值thre_S，筛选得到满足要求的HSnoise(n)具体为通过公式(1)计算心音信号频谱的标准偏差S，若大于阈值thre_S，则选择此时的HSnoise(n)，进入下一步；若小于或等于阈值，则σ＝σ+k，调整HSnoise(n)；
其中，ω为心音信号的频率，ωmax为X(ω)的最大频率，
为心音信号频谱幅值的平均值，k为预设的能够使标准差变大的常数。
作为优选，所述thre_S为通过公式(2)确定：
其中，N为原始心音信号的采样点总数，||HSorigin(n)||1为HSorigin(n)的l1范数。
步骤B.使用改进的经验模态分解的方法对HSnoise(n)进行分解，得到不同层数下的本征模函数；然后计算不同本征模函数与HSorigin(n)的相似度，筛选得到最大相似度时的本征模函数IMFmax(n)，最后基于IMFmax(n)，根据心音信号划分规则，对原始心音信号HSorigin(n)进行划分，得到si(i＝1,2,3,4)；
作为优选，步骤B通过以下过程实现：
SB.0初始化r(n)＝HSnoise(n)，其中r(n)为心音信号的剩余分量；
SB.1根据公式(3)、(4)搜索HSnoise(n)中全部的局部极大值点high(d)和极小值点low(c)：
high(d)＝{r(j)|r(j)＞r(j-1),r(j)＞r(j+1)} (3)
low(c)＝{r(j)|r(j)＜r(j-1),r(j)＜r(j+1)} (4)
SB.2根据公式(5)中加权系数优化下的埃尔米特插值法，计算r(n)全部相邻极大值点high(d)和极小值点low(c)的上下包络线Hh(x),Hl(x)：
其中，ε是加权系数，取值范围为[0,1]，Hh(x)为区间[dp-1,dp)的上包络线，Hl(x)为区间[cq-1,cq)的下包络线，dp-1,dp分别为high(d)的第p-1,p个采样点索引，cq-1,cq分别为low(c)的第q-1,q个采样点索引，high(·)′和low(·)′分别表示r(n)在·处的梯度值；
SB.3根据公式(6)，计算本征模函数IMF(n)和剩余分量r(n)：
SB.4对剩余分量r(n)重复步骤SB.1-SB.3，直到剩余分量没有极大值点或极小值点，最后得到本征模函数IMFlev(n),(lev＝1,2,...,LEV)，其中，IMFlev(n)为第lev层的本征模函数；
SB.5根据公式(7)定义的相似度函数，计算不同层数本征模函数IMFlev(n)与原始心音信号HSorigin(n)的相似度alev：
其中，υ(0≤υ≤1)为权重系数，||IMFlev||和||HSorigin||分别表示IMFlev和HSorigin的l2范数，N表示原始心音信号的采样点总数；
SB.6找到alev最大时的IMFlev(n)，记为IMFmax(n)；
SB.7根据公式(8)定义的阈值thre，找到HSorigin(n)中全部大于阈值的点，对点根据心音S1、S3段的定义，将其划分为s1,s3；找到全部小于阈值的点，对点根据心音S2、S4段的定义，将其划分为s2,s4：
thre＝mul·||IMFmax(n)||p (8)
其中，mul为常数，||IMFmax(n)||p为IMFmax(n)的p范数。
步骤C.使用基于线性系数的信息熵追踪演化计算法，求心音si(i＝1,2,3,4)的信息熵，即权重值η1，η2，η3，η4；
作为优选，本步骤通过以下过程实现：
SC.1对心音信号si(i＝1,2,3,4)进行相空间重构，嵌入维数为m，时间延迟为τ，si的样本点总数为Ni，则重构的相空间为：
Y(ng)＝(si(ng),si(ng+τ),...,si(ng+(m-1)τ)),(ng＝1,2,...,Ni-(m-1)τ) (9)
其中，Y(ng)为第ng个相空间的序列；
SC.2根据式(10)，计算原始心音信号的信息熵ηi。
其中，Mi＝Ni-(m-1)τ，
sum_L(r)由公式(11)定义，表示第r个相点与其他相点的距离和：
ave_L表示平均相点距离，由公式(12)定义：
其中，Lgf表示Y(ng)和Y(nf)的欧式距离，ave_L表示平均相点距离；
L_DK表示演化距离，由公式(13)定义：
其中，Loc_DL表示在
范围内使L1f最短时的f值，为初始相点，P_DK为演化点。
步骤D.基于步骤C得到的各心音段权重ηi(i＝1,2,3,4)，构建基于心音不同段的l1正则化全变分方程，对方程求解，得到保留病理特征的四段心音的第1次滤波结果(b1 1,b2 1,b3 1,b4 1)；
作为优选，本步骤通过以下子步骤实现：
SD.1对于第i段心音信号si，建立(14)所示的原始心音信号si的l1正则化全变分方程：
其中，
为第1次滤波后的第i段信号，bi为第i段理想的干净心音信号，
二范数的平方，ηi为权重，Hl为原始心音信号的下包络线，Hh为原始心音信号的上包络线，Δbi为bi的梯度；
SD.2根据式(15)，将式(14)转换为无约束的求解问题：
其中，λ1,λ2为拉格朗日乘子；
SD.3初始化
使用式(16)建立的改进的Split Bregman迭代算法，求解方程(15)，当k＞K1或
时，迭代结束，
即为si滤波后的结果，其中K1,ζ1为常数；
SD.4重复SD.1-SD.3，直到所有的心音段均滤波完成，并将结果记为
表示第1次滤波后心音段i的滤波结果。
步骤E.基于步骤C得到的各心音段权重ηi(i＝1,2,3,4)，构建基于心音不同段的限定l1正则化全变分方程，对方程求解，得到第2次滤波后的心音段
最后将
与步骤D得到的
进行滑动平均，得到滤波后的心音
最后将心音段整合，得到原始心音滤波后的信号HSfilter(n)。
作为优选，本步骤通过以下子步骤实现：
SE.1对于第i段心音信号si，建立(17)所示的原始心音信号si的限定l1正则化全变分方程：
其中，r(Δbi)＝τmin(|Δbi|,θ)+(1-τ)log(1+|Δbi|/θ)，θ(θ＞0)为设定的阈值，τ(0≤τ≤1)为正则化方法的权重值，
为第2次滤波后第i段信号；
SE.2初始化
使用式(18)建立的改进的Split Bregman迭代算法，求解方程(17)，当k＞K2或
时，迭代结束，
即为si滤波后的结果，其中K2,ζ2为常数；
SE.3重复SE.1-SE.2，直到所有的心音段均滤波完成，并将结果记为
表示第2次滤波后心音段i的滤波结果；
SE.4使用(19)的滑动平均法对两次滤波后的结果进行融合：
其中，
表示滤波结果，ρ∈(0,1)，代表两种滤波方式的权重；
SE.5根据心音信号的分段规则，将
按顺序整合，得到原始心音滤波后的信号HSfilter(n)。
有益效果
本发明一种基于病理自适应的心音全变分滤波方法，对比已有技术，具有如下有益效果：能够对原始的心音信号进行更为精确的去噪处理，并且，能够实现保留病理信息，为心音图的数字化分析提供基础。
附图说明
图1是本发明一种基于病理自适应的心音全变分滤波方法流程图。
图2是本发明高斯白噪声幅值与标准差间的关系。
图3是添加高斯白噪声后的信号频谱图。
图4是心音信号的上下包络线。
图5为不同层数的本征模函数。
图6为心音分段的结果。
图7是第一次滤波后的结果。
图8是第二次滤波后的结果。
图9是两次滤波结果的融合。
图10是低通滤波的结果。
图11是巴特沃斯滤波的结果。
具体实施方式
下面，结合附图并举实施例，对本发明进行详细描述。
实施例1
本发明提出了一种基于病理自适应的心音全变分滤波方法，如图1所示。具体包括如下步骤：
步骤1.对原始心音信号添加加性高斯白噪声，再使用频谱分析的方法得到加噪信号的频率分布情况；
作为优选，步骤1又具体包括如下子步骤：
S1.A将高斯白噪声添加到心音信号中，得到加噪的信号HSnoise(n)，图2为添加的高斯白噪声与标准差的关系。图2中“δ”为高斯白噪声标准偏差值，“Sample points”为样本点个数，z轴为高斯白噪声的幅值。随着标准偏差值的增大，高斯白噪声的幅值的不确定区间也在不断增大，因此对于心音信号，如果标准偏差过小，则高斯白噪声的幅值过小，无法有效地降低混叠现象；如果标准偏差过大，则高斯白噪声的幅值不确定性区间增大，在经验模态分解时容易将原始心音信号的特征丢失，所以高斯白噪声应选择合适的标准偏差。
S1.B对HSnoise(n)进行离散傅里叶变换，得到频率分布X(ω)。图3为添加σ＝0.1的高斯白噪声后的频率分布图X(ω)，横坐标为ω，纵坐标为幅值，由频谱图可知，由于心音信号为低频信号，即使添加了噪声，能量仍主要集中在0～200Hz，高于200Hz的能量较少。
S1.C通过公式(1)计算心音信号频谱的标准偏差S，若大于阈值thre_S，则选择此时的HSnoise(n)，进入下一步；若小于或等于阈值，则σ＝σ+k，调整HSnoise(n)；为了使信号能量分布均匀，为在进行经验模态分解时可以有效地防止模态的混叠，本实例中计算得到thre_S＝0.3，k＝0.1。重复S1.A-S1.C，确定得到σ＝0.3。
步骤2.使用改进的经验模态分解的方法对HSnoise(n)进行分解，得到不同层数下的本征模函数；然后计算不同本征模函数与HSorigin(n)的相似度，筛选得到最大相似度时的本征模函数IMFmax(n)，最后基于IMFmax(n)，根据心音信号划分规则，对原始心音信号HSorigin(n)进行划分，得到si(i＝1,2,3,4)；
作为优选，步骤2又具体包括如下子步骤：
S2.O初始化r(n)＝HSnoise(n)，其中r(n)为心音信号的剩余分量；
S2.A根据公式(3)、(4)搜索HSnoise(n)中全部的极大值点high(k)和极小值点low(m)；
S2.B根据公式(5)中加权系数优化下的埃尔米特插值法，计算r(n)全部相邻极大值点和极小值点的上下包络线Hh(x),Hl(x)。图4显示了计算得到信号的上下包络线，其中，HSnoise(n)由实线表示，上包络线由虚线表示，下包络线由星号表示，由图4可知，上下包络线整体比较光滑，在极值点处不会产生梯度消失的现象，可以有效地反映信号幅度变化；
S2.C根据公式(6)，计算本征模函数IMF(n)和剩余分量r(n)；
S2.D对剩余分量r(n)重复步骤SB.1-SB.3，直到剩余分量没有极大值点或极小值点，最后得到本征模函数IMFlev(n),(lev＝1,2,...,LEV)，其中，IMFlev(n)为第lev层的本征模函数。图5为分解层数1～8时的本征模函数，横坐标为采样点索引，纵坐标为本征模函数幅值，在分解层数为8时，信号达到要求。
S2.E根据公式(7)定义的相似度函数，υ＝0.5，计算不同层数本征模函数IMFlev(n)与原始心音信号HSorigin(n)的相似度alev，计算结果如下：
a1＝2.46×10-4,a2＝3.91×10-3,a3＝5.19×10-4,a4＝1.11×10-3
a5＝2.31×10-3,a6＝5.83×10-4,a7＝1.43×10-3,a8＝2.74×10-4
S2.F由上述结果可知，第2层IMF与原始心音信号的相似度最大，此时，a2＝3.91×10-3，最大时的IMF2(n)，记为IMFmax(n)；
S2.G根据公式(8)定义的阈值thre，mul＝0.03找到全部大于阈值的点，对点根据心音S1、S3段的定义，将其划分为s1,s3；找到全部小于阈值的点，对点根据心音S2、S4段的定义，将其划分为s2,s4，图6为分段结果，虚线为分割线，实线为心音信号；
步骤3.使用基于线性系数的心音信号信息熵追踪演化计算法，求心音S1、S2、S3和S4段的权重值η1,η2,η3,η4；
作为优选，步骤3又包括如下步骤：
S3.A根据公式(9)，对心音信号si进行相空间重构，嵌入维数为m＝5，时间延迟为τ＝3，si的样本点总数为Ni＝200，则重构的相空间为：
Y(ng)＝(si(ng),si(ng+3),...,si(ng+(5-1)×3)),(ng＝1,2,...,200-(5-1)×3)
其中，Y(ng)为第ng个相空间的序列。
S3.B根据式(12)，计算心音信号信息熵η1＝0.062,η2＝0.057,η3＝0.036,η4＝0.053。
步骤4.基于步骤3得到的各心音段ηi(i＝1,2,3,4)，构建l1正则化全变分方程，然后使用改进的Slip-Bregman迭代方法，对上述方程求解，得到保留病理特征的四段心音的第1次滤波结果
作为优选，步骤4又具体包括如下子步骤：
S4.A根据公式(14)，建立原始心音信号si的全变分方程；
S4.B将式(14)转换为无约束的求解问题，如式(15)所示
S4.C初始化k＝0,b0＝0,m0＝0,d0＝0，使用式(16)建立的Split Bregman迭代算法，求解方程(15)；当k＞K1或
时，迭代结束，
即为滤波后的心音，其中K1＝100,ζ1＝0.001。
S4.D重复S4.A-S4.C，直到所有的心音段均滤波完成，并将结果记为
表示第1次滤波后心音段i的滤波结果。滤波结果见图7，使用此滤波方式，能够有效地保留细节纹理信息，但是仍有部分噪声无法有效地滤除。
步骤5.构建基于心音不同段的限定l1正则化全变分方程，求解后得到滤波后的心音段，然后使用滑动平均法得到滤波后的心音；
步骤5又具体包括如下子步骤：
S5.Aθ＝0.7，ω＝0.3，建立(19)所示的原始心音信号si的全变分方程；
S5.B初始化k＝0,b0＝0,m0＝0,d0＝0，使用式(18)建立的改进的Split Bregman迭代算法，求解方程(17)，当k＞K2或
时，迭代结束，
即为si滤波后的结果，其中K2＝20,ζ2＝0.0005；
S5.C重复SE.1-SE.2，直到所有的心音段均滤波完成，并将结果记为
滤波结果见图8，使用此滤波方式，能够有效地保留信号的概略信息，去除了全部噪声，但是会导致信号的部分失真和细节纹理特征的丢失；
S5.D令ρ＝0.4使用(19)的滑动平均法对两次滤波后的结果进行融合；
S5.E根据心音信号的分段规则，将
按顺序整合，得到原始心音滤波后的信号HSfilter(n)，滤波结果见图9，信噪比为55.5；图10与图11分别为低通滤波和巴特沃斯滤波结果，信噪比分别为0.9和26.4。由滤波结果可知，低通滤波和巴特沃斯滤波等常见的滤波器，会丢失细节纹理信息，且滤波后的信噪比低于本文提出的方法；使用本申请中的滤波方式，不仅能够有效地保留细节纹理信息，而且能够有效地保留信号的概略信息，去除噪声。
需要说明的是，本说明书所述的仅为本发明的较佳实施例而已，以上实施例仅用于说明本发明的技术方案而非对本发明的限制。凡本领域技术人员依本发明的构思通过逻辑分析、推理或者有限的实验可以得到的技术方案，皆应在本发明的范围之内。
Claims (6)
1.一种基于病理自适应的心音全变分滤波方法，其特征在于：包括以下步骤：
步骤A.对原始心音信号HSorigin(n)添加加性高斯白噪声，基于加噪信号HSnoise(n)的频率分布情况X(ω)和心音信号频谱的标准偏差，调整高斯白噪声的标准偏差σ，根据设定的阈值thre_S，筛选得到满足要求的HSnoise(n)；
步骤B.使用改进的经验模态分解的方法对HSnoise(n)进行分解，得到不同层数下的本征模函数；然后计算不同本征模函数与HSorigin(n)的相似度，筛选得到最大相似度时的本征模函数IMFmax(n)，最后基于IMFmax(n)，根据心音信号划分规则，对原始心音信号HSorigin(n)进行划分，得到si(i＝1,2,3,4)；
步骤C.使用基于线性系数的信息熵追踪演化计算法，求心音si(i＝1,2,3,4)的信息熵，即权重值η1，η2，η3，η4；过程如下：
SC.1对心音信号si(i＝1,2,3,4)进行相空间重构，嵌入维数为m，时间延迟为τ，si的样本点总数为Ni，则重构的相空间为：
Y(ng)＝(si(ng),si(ng+τ),...,si(ng+(m-1)τ)),(ng＝1,2,...,Ni-(m-1)τ) (9)
其中，Y(ng)为第ng个相空间的序列；
SC.2根据式(10)，计算原始心音信号的信息熵ηi；
其中，Mi＝Ni-(m-1)τ，
sum_L(r)由公式(11)定义，表示第r个相点与其他相点的距离和：
ave_L表示平均相点距离，由公式(12)定义：
其中，Lgf表示Y(ng)和Y(nf)的欧式距离，ave_L表示平均相点距离；
L_DK表示演化距离，由公式(13)定义：
其中，Loc_DL表示在
范围内使L1f最短时的f值，为初始相点，L1f表示g＝1时的Lgf，P_DK为演化点；
步骤D.基于步骤C得到的各心音段权重ηi(i＝1,2,3,4)，构建基于心音不同段的l1正则化全变分方程，对方程求解，得到保留病理特征的四段心音的第1次滤波结果
步骤E.基于步骤C得到的各心音段权重ηi(i＝1,2,3,4)，构建基于心音不同段的限定l1正则化全变分方程，对方程求解，得到第2次滤波后的心音段
最后将
与步骤D得到的
进行滑动平均，得到滤波后的心音
最后将心音段整合，得到原始心音滤波后的信号HSfilter(n)；所述限定l1正则化全变分方程为：
其中，r(Δbi)＝τmin(|Δbi|,θ)+(1-τ)log(1+|Δbi|/θ)，θ(θ＞0)为设定的阈值，Δbi为bi的梯度；τ(0≤τ≤1)为正则化方法的权重值，
为第2次滤波后第i段信号；‖算式‖2表示算式的2范数。
2.根据权利要求1所述的方法，其特征在于：所述调整高斯白噪声的标准偏差σ，根据设定的阈值thre_S，筛选得到满足要求的HSnoise(n)具体为通过公式(1)计算心音信号频谱的标准偏差S，若大于阈值thre_S，则选择此时的HSnoise(n)，进入下一步；若小于或等于阈值，则σ＝σ+k，调整HSnoise(n)；
其中，ω为心音信号的频率，ωmax为X(ω)的最大频率，
为心音信号频谱幅值的平均值，k为预设的能够使标准差变大的常数。
3.根据权利要求1所述的方法，其特征在于：所述thre_S为通过公式(2)确定：
其中，N为原始心音信号的采样点总数，||HSorigin(n)||1为HSorigin(n)的l1范数。
4.根据权利要求1所述的方法，其特征在于：步骤B通过以下过程实现：
SB.0初始化r(n)＝HSnoise(n)，其中r(n)为心音信号的剩余分量；
SB.1根据公式(3)、(4)搜索HSnoise(n)中全部的局部极大值点high(d)和极小值点low(c)：
high(d)＝{r(j)|r(j)＞r(j-1),r(j)＞r(j+1)} (3)
low(c)＝{r(j)|r(j)＜r(j-1),r(j)＜r(j+1)} (4)
SB.2根据公式(5)中加权系数优化下的埃尔米特插值法，计算r(n)全部相邻极大值点high(d)和极小值点low(c)的上下包络线Hh(x),Hl(x)：
其中，ε是加权系数，取值范围为[0,1]，Hh(x)为区间[dp-1,dp)的上包络线，Hl(x)为区间[cq-1,cq)的下包络线，dp-1,dp分别为high(d)的第p-1,p个采样点索引，cq-1,cq分别为low(c)的第q-1,q个采样点索引，high(·)′和low(·)′分别表示r(n)在·处的梯度值；
SB.3根据公式(6)，计算本征模函数IMF(n)和剩余分量r(n)：
SB.4对剩余分量r(n)重复步骤SB.1-SB.3，直到剩余分量没有极大值点或极小值点，最后得到本征模函数IMFlev(n),(lev＝1,2,...,LEV)，其中，IMFlev(n)为第lev层的本征模函数；
SB.5根据公式(7)定义的相似度函数，计算不同层数本征模函数IMFlev(n)与原始心音信号HSorigin(n)的相似度alev：
其中，υ(0≤υ≤1)为权重系数，||IMFlev||和||HSorigin||分别表示IMFlev和HSorigin的l2范数，N表示原始心音信号的采样点总数；
SB.6找到alev最大时的IMFlev(n)，记为IMFmax(n)；
SB.7根据公式(8)定义的阈值thre，找到HSorigin(n)中全部大于阈值的点，对点根据心音S1、S3段的定义，将其划分为s1,s3；找到全部小于阈值的点，对点根据心音S2、S4段的定义，将其划分为s2,s4：
thre＝mul·||IMFmax(n)||p (8)
其中，mul为常数，||IMFmax(n)||p为IMFmax(n)的p范数。
5.根据权利要求1所述的方法，其特征在于：所述步骤D通过以下过程实现：
SD.1对于第i段心音信号si，建立(14)所示的原始心音信号si的l1正则化全变分方程：
其中，
为第1次滤波后的第i段信号，bi为第i段理想的干净心音信号，
二范数的平方，ηi为权重，Hl为原始心音信号的下包络线，Hh为原始心音信号的上包络线，Δbi为bi的梯度；
SD.2根据式(15)，将式(14)转换为无约束的求解问题：
其中，λ1,λ2为拉格朗日乘子；
SD.3初始化k＝0,
使用式(16)建立的改进的Split Bregman迭代算法，求解方程(15)，当k＞K1或
时，迭代结束，
即为si滤波后的结果，其中K1,ζ1为常数；
其中，‖算式‖2表示算式的2范数；‖算式‖1表示算式的1范数；
SD.4重复SD.1-SD.3，直到所有的心音段均滤波完成，并将结果记为
表示第1次滤波后心音段i的滤波结果。
6.根据权利要求1-5任一所述的方法，其特征在于：所述步骤E通过以下过程实现：
SE.1对于第i段心音信号si，建立(17)所示的原始心音信号si的限定l1正则化全变分方程：
其中，r(Δbi)＝τmin(|Δbi|,θ)+(1-τ)log(1+|Δbi|/θ)，θ(θ＞0)为设定的阈值，τ(0≤τ≤1)为正则化方法的权重值，
为第2次滤波后第i段信号；
SE.2初始化k＝0,
使用式(18)建立的改进的Split Bregman迭代算法，求解方程(17)，当k＞K2或
时，迭代结束，
即为si滤波后的结果，其中K2,ζ2为常数；
其中，‖算式‖2表示算式的2范数；‖算式‖1表示算式的1范数；
SE.3重复SE.1-SE.2，直到所有的心音段均滤波完成，并将结果记为
表示第2次滤波后心音段i的滤波结果；
SE.4使用(19)的滑动平均法对两次滤波后的结果进行融合：
其中，
表示滤波结果，ρ∈(0,1),1-ρ∈(0,1)，代表两种滤波方式的权重；
SE.5根据心音信号的分段规则，将
按顺序整合，得到原始心音滤波后的信号HSfilter(n)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910587766.7A CN110263760B (zh) | 2019-07-02 | 2019-07-02 | 一种基于病理自适应的心音全变分滤波方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910587766.7A CN110263760B (zh) | 2019-07-02 | 2019-07-02 | 一种基于病理自适应的心音全变分滤波方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110263760A CN110263760A (zh) | 2019-09-20 |
| CN110263760B true CN110263760B (zh) | 2020-09-18 |
Family
ID=67923711
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910587766.7A Active CN110263760B (zh) | 2019-07-02 | 2019-07-02 | 一种基于病理自适应的心音全变分滤波方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110263760B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116702015B (zh) * | 2023-08-09 | 2023-11-28 | 中国科学技术大学 | 信号处理方法、装置、无人机及计算机可读存储介质 |
-
2019
- 2019-07-02 CN CN201910587766.7A patent/CN110263760B/zh active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN110263760A (zh) | 2019-09-20 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Fahad et al. | Microscopic abnormality classification of cardiac murmurs using ANFIS and HMM | |
| Dokur et al. | Heart sound classification using wavelet transform and incremental self-organizing map | |
| CN108470156B (zh) | 一种心音信号分类识别方法 | |
| CN111358455B (zh) | 一种多数据源的血压预测方法和装置 | |
| CN110970042B (zh) | 一种电子听诊器的肺部啰音人工智能实时分类方法、系统、装置及可读存储介质 | |
| CN111368627A (zh) | Cnn结合改进频率小波切片变换的心音分类方法及系统 | |
| CN114469124B (zh) | 一种运动过程中异常心电信号的识别方法 | |
| Avanzato et al. | Heart sound multiclass analysis based on raw data and convolutional neural network | |
| CN112307959B (zh) | 一种用于心电信号分析的小波去噪方法 | |
| CN111317499A (zh) | 一种基于小波技术的心音信号处理方法 | |
| Abdollahpur et al. | Cycle selection and neuro-voting system for classifying heart sound recordings | |
| CN112580486A (zh) | 一种基于雷达微多普勒信号分离的人体行为分类方法 | |
| Sujadevi et al. | A hybrid method for fundamental heart sound segmentation using group-sparsity denoising and variational mode decomposition | |
| CN110263760B (zh) | 一种基于病理自适应的心音全变分滤波方法 | |
| Arı et al. | A lung sound classification system based on data augmenting using ELM-wavelet-AE | |
| Celik | CovidCoughNet: A new method based on convolutional neural networks and deep feature extraction using pitch-shifting data augmentation for covid-19 detection from cough, breath, and voice signals | |
| Guermoui et al. | Heart sounds analysis using wavelets responses and support vector machines | |
| CN111938691B (zh) | 一种基础心音识别方法及设备 | |
| Morshed et al. | Automated heart valve disorder detection based on PDF modeling of formant variation pattern in PCG signal | |
| Patidar et al. | Classification of heart disorders based on tunable-Q wavelet transform of cardiac sound signals | |
| Kumar et al. | An adaptive approach to abnormal heart sound segmentation | |
| Shafi et al. | Validity-guided fuzzy clustering evaluation for neural network-based time-frequency reassignment | |
| Kaushal et al. | Deep Autoencoder Neural Networks for Heart Sound Classification | |
| Ge et al. | Detection of pulmonary arterial hypertension associated with congenital heart disease based on time–frequency domain and deep learning features | |
| Fahad et al. | Phonocardiogram heartbeat segmentation and autoregressive modeling for person identification |
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 |