CN112562701B - 心音信号双通道自适应降噪算法、装置、介质及设备 - Google Patents
心音信号双通道自适应降噪算法、装置、介质及设备 Download PDFInfo
- Publication number
- CN112562701B CN112562701B CN202011278315.4A CN202011278315A CN112562701B CN 112562701 B CN112562701 B CN 112562701B CN 202011278315 A CN202011278315 A CN 202011278315A CN 112562701 B CN112562701 B CN 112562701B
- Authority
- CN
- China
- Prior art keywords
- signal
- preprocessing
- filter
- noise reduction
- channel 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
- 230000009467 reduction Effects 0.000 title claims abstract description 66
- 230000005236 sound signal Effects 0.000 title claims abstract description 64
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 59
- 238000007781 pre-processing Methods 0.000 claims abstract description 112
- 230000003044 adaptive effect Effects 0.000 claims abstract description 104
- 238000012546 transfer Methods 0.000 claims abstract description 52
- 238000000034 method Methods 0.000 claims abstract description 33
- 238000005070 sampling Methods 0.000 claims abstract description 29
- 230000008569 process Effects 0.000 claims abstract description 23
- 238000001914 filtration Methods 0.000 claims abstract description 19
- 238000011946 reduction process Methods 0.000 claims abstract description 14
- 230000006870 function Effects 0.000 claims description 48
- 230000007613 environmental effect Effects 0.000 claims description 24
- 238000004364 calculation method Methods 0.000 claims description 10
- 238000005457 optimization Methods 0.000 claims description 8
- 230000000747 cardiac effect Effects 0.000 claims description 6
- 238000003860 storage Methods 0.000 claims description 5
- 238000004590 computer program Methods 0.000 claims description 4
- 230000009191 jumping Effects 0.000 claims description 4
- 230000000694 effects Effects 0.000 abstract description 9
- 239000013598 vector Substances 0.000 description 16
- 238000002555 auscultation Methods 0.000 description 9
- 238000010586 diagram Methods 0.000 description 7
- 230000002159 abnormal effect Effects 0.000 description 3
- 239000000654 additive Substances 0.000 description 3
- 230000000996 additive effect Effects 0.000 description 3
- 230000008859 change Effects 0.000 description 3
- 230000004083 survival effect Effects 0.000 description 3
- 238000005314 correlation function Methods 0.000 description 2
- 238000003745 diagnosis Methods 0.000 description 2
- 238000009826 distribution Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012805 post-processing Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 230000004044 response Effects 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 230000001133 acceleration Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000036772 blood pressure Effects 0.000 description 1
- 230000001684 chronic effect Effects 0.000 description 1
- 230000000052 comparative effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 230000002526 effect on cardiovascular system Effects 0.000 description 1
- 238000010884 ion-beam technique Methods 0.000 description 1
- 230000003137 locomotive effect Effects 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000007704 transition Effects 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
- G10L19/00—Speech or audio signals analysis-synthesis techniques for redundancy reduction, e.g. in vocoders; Coding or decoding of speech or audio signals, using source filter models or psychoacoustic analysis
- G10L19/04—Speech or audio signals analysis-synthesis techniques for redundancy reduction, e.g. in vocoders; Coding or decoding of speech or audio signals, using source filter models or psychoacoustic analysis using predictive techniques
- G10L19/26—Pre-filtering or post-filtering
- G10L19/265—Pre-filtering, e.g. high frequency emphasis prior to encoding
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B7/00—Instruments for auscultation
- A61B7/02—Stethoscopes
- A61B7/04—Electric stethoscopes
-
- 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
- G10L21/00—Speech or voice signal processing techniques to produce another audible or non-audible signal, e.g. visual or tactile, in order to modify its quality or its intelligibility
- G10L21/02—Speech enhancement, e.g. noise reduction or echo cancellation
- G10L21/0208—Noise filtering
- G10L21/0216—Noise filtering characterised by the method used for estimating noise
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
- Y02D30/00—Reducing energy consumption in communication networks
- Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks
Landscapes
- Engineering & Computer Science (AREA)
- Acoustics & Sound (AREA)
- Health & Medical Sciences (AREA)
- Physics & Mathematics (AREA)
- Computational Linguistics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Multimedia (AREA)
- Human Computer Interaction (AREA)
- Audiology, Speech & Language Pathology (AREA)
- Signal Processing (AREA)
- Biomedical Technology (AREA)
- General Health & Medical Sciences (AREA)
- Public Health (AREA)
- Veterinary Medicine (AREA)
- Animal Behavior & Ethology (AREA)
- Surgery (AREA)
- Molecular Biology (AREA)
- Medical Informatics (AREA)
- Heart & Thoracic Surgery (AREA)
- Quality & Reliability (AREA)
- Soundproofing, Sound Blocking, And Sound Damping (AREA)
Abstract
本发明提供了一种心音信号双通道自适应降噪算法、装置、介质及设备；其中算法包括预处理滤波器生成过程和主副通道信号预处理及降噪过程；预处理滤波器生成过程，是指在安静的环境下采样主通道信号，采用全极点模型对主通道信号作线性预测编码，求得全极点模型的传递函数，从而得到与全极点模型传递函数互为倒数的、预处理滤波器的传递函数；主副通道信号预处理及降噪过程是指：采样主通道信号和副通道信号；对主通道信号和副通道信号进行预处理，并对预处理后信号进行自适应滤波得到误差信号；之后采用传递函数与预处理滤波器传递函数互为倒数的滤波器对误差信号进行滤波处理，得到输出信号。本发明根据实测信号，以最小化心音信号方差为目标确定预处理滤波器的参数，所得预处理滤波器参数为方差最小化意义下的最优值，从而提升降噪效果。
Description
技术领域
本发明涉及医学测量和信号处理技术领域，尤其涉及一种心音信号双通道自适应降噪算法、装置、介质及设备。
背景技术
远程听诊可在心血管慢性病患者随诊上发挥重要作用，让患者在家享受医疗服务，而无需到大医院就诊，极大地减轻患者的负担，降低社会和经济成本。
心音听诊有效实施的前提在于获得高质量的心音信号。但远程听诊应用场合多样，不同强度水平和种类的环境噪声类型丰富；加之在远程听诊过程中，医生不了解患者所处环境的情况，难以判断所听到的异常声音为患者的额外心音或杂音还是环境噪声，容易出现误诊。且环境噪声极易掩盖微弱的异常心音信号，相当部分环境噪声的时频分布与异常心音的时频分布高度重叠，其降噪难度非常高。
双通道自适应滤波算法计算量小，易于硬件实现，且设计合理时可大幅提高信噪比，故常应用于语音信号降噪。其原理为：假设环境噪声为加性，主通道测量带噪心音，副通道采集环境噪音，副通道所测环境噪音经处理后用于部分抵消带噪心音中的噪音而实现降噪。
双通道自适应滤波算法通常采用最小均方(least mean square，LMS)、归一化最小均方(normalized least mean square，NLMS)算法或其改进算法优化滤波器的权值来实现。但是LMS或NLMS算法依赖于所测量信号以及噪声的统计特性，信号统计特性突变的情况下难以在快速收敛和低失调(misalignment)之间实现平衡。而信号统计特性突变在心音听诊过程中经常出现，特别是第一第二心音的幅值会周期性出现短时间内的大幅变化。变步长的LMS或NLMS算法可提高降噪性能，但对于第一第二心音这类幅值短时快速变化的信号，要确定合适的变步长策略并不容易。
与常用的变步长算法思路不同，在中国发明专利申请《用于采集体音信号的双麦克风自适应滤波算法及应用》(公开号：CN109545239A)中针对这种非平稳信号建立了一种双通道自适应滤波算法，以解决体音信号的降噪问题。双通道信号分别为主通道和副通道信号；主通道信号为带噪体音信号，副通道信号为环境噪音；对主、副通道信号作相同的预处理，对预处理后的主、副通道信号采用LMS或NLMS算法计算自适应滤波器权值并计算误差信号，以滤除主通道信号中的环境噪音；对误差信号作第一次低通滤波处理以复原体音信号，从而得到自适应滤波算法输出的体音信号。
该自适应降噪算法中所采用的主、副通道信号预处理方法可降低主通道信号的方差，提高降噪性能。但是在该专利中，其预处理方法所用的滤波器参数需要人为根据经验选定，且优选地由若干个一阶预加重环节串联而成；如此所得的滤波器效果会由于人为选定的滤波器参数不同而滤波效果有差异，不能确保滤波器取得最佳效果。更为理想的做法应为根据实际信号自动确定滤波器参数，而这仍是尚未解决的问题。
发明内容
为克服现有技术中的缺点与不足，本发明的目的在于提供一种心音信号双通道自适应降噪算法、装置、介质及设备；本发明根据实测信号，以最小化带噪心音信号方差为目标确定预处理滤波器的参数，所得预处理滤波器的参数为方差最小化意义下的最优值，从而提升降噪效果。
为了达到上述目的，本发明通过下述技术方案予以实现：一种心音信号双通道自适应降噪算法，其特征在于：经主通道采集带噪心音信号，以及副通道采集环境噪声；包括预处理滤波器生成过程和主副通道信号预处理及降噪过程；
所述预处理滤波器生成过程，是指：在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
所述主副通道信号预处理及降噪过程，是指：采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成过程得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)。
优选地，所述预处理滤波器生成过程，包括如下步骤：
SG1步，设定预处理滤波器的阶数N，其中N≥2，预处理滤波器的传递函数为H(z)＝1+a1z-1+...+aNz-N，其中ai，i∈{1，...，N}，为待定参数；
SG2步，在安静的环境下，以采样频率Fs测量主通道信号d′(k)，主通道信号d′(k)持续时间超过一次心动周期；
SG3步，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码(LinearPrediction Coding，LPC)，求得全极点模型的传递函数
SG4步，将bi，i∈{1，...，N}，代入到ai中，得到预处理滤波器的传递函数H(z)；
所述主副通道信号预处理及降噪过程，包括如下步骤：
SA1步，初始化当前时刻序号k＝0，自适应滤波器权值w(0，j)＝0，j∈{0，...，M-1}，其中M为自适应滤波器阶数，M＞N；
SA2步，以采样频率Fs分别测量当前时刻主通道信号d(k)和副通道信号n2(k)；
SA3步，计算当前时刻的预处理结果
和/>
若k＜N，则主通道预处理后信号
副通道预处理后信号
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若N≤k＜M，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；/>
若k≥M，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
SA4步，计算自适应滤波器输出y(k)：
SA5步，计算误差信号e(k)：
SA6步，更新自适应滤波器权值w(k+1，j)：
其中，ζ为一很小的正数以防止分母为0；η2为自适应滤波器权值优化的步长；
SA7步，采用传递函数为
的滤波器对误差信号e(k)进行滤波得到输出o(k)；
SA8步，判断自适应降噪终止指示变量：若自适应降噪终止指示变量为真，则结束自适应降噪算法，否则k＝k+1，跳至SA2步计算下一时刻自适应降噪算法的输出o(k)。
优选地，所述自适应滤波器阶数M的取值范围为：M∈[20，300]。
优选地，所述自适应滤波器权值优化的步长η2的取值范围为：η2∈[0.1，1]。
优选地，所述采样频率Fs≥4KHz。
优选地，所述SG2步中，主通道信号d′(k)持续时间的取值范围为[1s，3s]。
优选地，所述SG3步的线性预测编码中，采用Levinson-Durbin算法计算出全极点模型的传递函数G(z)。
本发明算法的技术原理是：
假设环境噪声为加性，主通道信号为带噪心音d(k)＝s(k)+n1(k)，其中s(k)为第k时刻纯净的心音信号，而n1(k)为第k时刻主通道所含的环境噪声信号；副通道信号为第k时刻环境噪音n2(k)，其经自适应滤波器加权后得到y(k)用于抵消带噪心音中的环境噪声n1(k)而实现降噪；其自适应滤波器权值w(k+1，j)的迭代算式分别如式(1)和式(2)所示：
w(k+1，j)＝w(k，j)+η1n2(k-j)[s(k)+n1(k)-y(k)] (1)
其中，e(k)＝d(k)-y(k)＝s(k)+n1(k)-y(k)，w(k+1，j)和w(k，j)分别为第k+1和第k时刻自适应滤波器的权值，n2(k-j)为第k-j时刻的副通道信号，η1和η2均为自适应滤波器权值优化的步长，ζ为一很小的正数以防止分母为0，M为自适应滤波器的阶数。
为方便起见，记自适应滤波器权值向量W(k)＝{w(k，0)，...，w(k，M-1)}。
假设心音信号s(k)和环境噪声n1(k)、n2(k)均为平稳信号，且假设步长η1和η2取值恰当使得自适应滤波器的权值迭代过程收敛；同时假设环境噪声n1(k)、n2(k)的互相关函数不随时间变化，则图3所示系统存在维纳解(最优解)，记为W*＝{w*(0)，...，w*(M-1)}，且权值向量W(k)可收敛到维纳解附近。但根据式(1)和(2)及e(k)的定义，图3所示系统进入稳态后，e(k)和权值向量W(k)依然会随s(k)和n1(k)、n2(k)的随机变化而波动。
进入稳态后e(k)的方差和权值向量偏差的范数(定义为||ΔW||2＝||W*-W(k)||2)均依赖于s(k)的方差。特别是在心音听诊中，第一心音和第二心音的幅值往往远高于环境噪声，因而导致偏差e(k)＝s(k)+n1(k)-y(k)在自适应滤波器权值收敛的过程中周期性跃变，并进而引起自适应滤波器权值周期性过调。减小s(k)的方差有助于降低稳态后e(k)的方差和权值向量偏差的范数||ΔW||2。
在本发明中，预处理滤波器参数通过对s(k)作线性预测编码求得，如此，可有效降低s(k)的方差：
对s(k)作线性预测编码可表示为
其中，参数bi，i∈{1，...，N}，通过最小化式(4)的性能指标J得到
其中，E{}表示求期望值。显然，使J取最小值的参数bi，i∈{1，...，N}，必能保证
此外，根据心音信号的特点，可假设E{s(k)}＝0，故
因此，对于心音降噪而言，E{s2(k)}和/>
即分别为s(k)和/>
的方差；/>
则表明线性预测编码有助于减小信号的方差。注意到
仅当s(k)与s(k-i)，i∈{1，...，N}，均线性无关时等号才成立，而当采样频率高于奈奎斯特频率时，采样时刻相近的心音信号之间存在显著的线性相关性，因此一般而言，/>
即线性预测编码有助于大幅减小信号的方差。
式(3)对应的传递函数即为
一种心音信号双通道自适应降噪装置，其特征在于，包括：
预处理滤波器生成模块，用于在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
以及主副通道信号预处理及降噪模块，用于采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成模块得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)。
一种存储介质，其特征在于，其中所述存储介质存储有计算机程序，所述计算机程序当被处理器执行时使所述处理器执行上述心音信号双通道自适应降噪算法。
一种计算设备，包括处理器以及用于存储处理器可执行程序的存储器，其特征在于，所述处理器执行存储器存储的程序时，实现上述心音信号双通道自适应降噪算法。
与现有技术相比，本发明具有如下优点与有益效果：
1、本发明根据实测信号，以最小化心音信号方差为目标确定预处理滤波器的参数，故所得预处理滤波器的参数为方差最小化意义下的最优值，从而提升降噪效果；
2、当预处理滤波器阶数较高时，根据经验确定较佳滤波器的参数难度很高；而本发明以最小化方差为目标确定最优滤波器参数，不受滤波器阶数影响；
3、本发明算法计算量小，响应速度快，对硬件设备的计算能力要求低，特别适用于小型可穿戴听诊设备和小型电子听诊器，同时本发明算法也适合在医院和家庭用电子听诊辅助诊疗系统中应用。
附图说明
图1为本发明算法的预处理滤波器生成过程流程图；
图2为本发明算法的主副通道信号预处理及降噪过程流程图；
图3为未采用主副通道信号预处理的自适应降噪算法的原理图；
图4为增加主副通道信号预处理的自适应降噪算法的原理图；
图5(a)和图5(b)分别为实施例一中在安静环境下所测得的主通道心音信号及其预处理后波形；
图6(a)和图6(b)分别为实施例二中在安静环境下所测得的主通道心音信号及其预处理后波形；
图7(a)至图7(c)分别为实施例一中带噪心音信号、自适应滤波后心音信号和相邻两时刻自适应滤波器权值向量差的范数||ΔW(k)||2＝||W(k+1)-W(k)||2随时间变化的曲线；
图8(a)至图8(c)分别为实施例三中带噪心音信号、自适应滤波后心音信号和相邻两时刻自适应滤波器权值向量差的范数||ΔW(k)||2＝||W(k+1)-W(k)||2随时间变化的曲线。
具体实施方式
下面结合附图与具体实施方式对本发明作进一步详细的描述。
实施例一
本实施例一种心音信号双通道自适应降噪算法，经主通道采集带噪心音信号，以及副通道采集环境噪声；包括预处理滤波器生成过程和主副通道信号预处理及降噪过程；
所述预处理滤波器生成过程，是指：在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
所述主副通道信号预处理及降噪过程，是指：采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成过程得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)。
具体地说，主通道信号和副通道信号的获取有多种方式，例如通过麦克风检测声音信号来实现、通过压电薄膜传感器或加速度传感器来检测振动或者位移信号来实现等。
预处理滤波器生成过程，如图1所示，包括如下步骤：
SG1步，设定预处理滤波器的阶数N，其中N≥2，预处理滤波器的传递函数为H(z)＝1+a1z-1+...+aNz-N，其中参数ai，i∈{1，...，N}，待定；预处理滤波器为阶数为N的有限冲击响应(finite impulse response，FIR)模型；
SG2步，在安静的环境下(环境噪声幅值较小)，以采样频率Fs测量主通道信号d′(k)，采样频率Fs优选为Fs≥4KHz，主通道信号d′(k)持续时间超过一次心动周期，优选为1s至3s；主通道信号d′(k)由心音信号s(k)和环境噪声n1(k)叠加而成；此时，环境噪声n1(k)幅值较小，且|s(k)|＞＞|n1(k)|，故d′(k)≈s(k)；
SG3步，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码(LinearPrediction Coding，LPC)，求得全极点模型的传递函数
优选采用Levinson-Durbin算法计算出全极点模型的传递函数G(z)；
SG4步，将bi，i∈{1，...，N}，代入到ai中，得到预处理滤波器的传递函数H(z)。
主副通道信号预处理及降噪过程，如图2所示，包括如下步骤：
SA1步，初始化当前时刻序号k＝0，自适应滤波器权值w(0，j)＝0，j∈{0，...，M-1}，其中M为自适应滤波器阶数，M＞N；M的取值范围优选为：M∈[20,300]；
SA2步，以采样频率Fs分别测量当前时刻主通道信号d(k)和副通道信号n2(k)；SA2步中采样频率与SG2步中采样频率相同，均为Fs；
SA3步，计算当前时刻的预处理结果
和/>
若k＜N，则主通道预处理后信号
副通道预处理后信号
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若N≤k＜M，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若k≥M，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
SA4步，计算自适应滤波器输出y(k)：
SA5步，计算误差信号e(k)：
SA6步，更新自适应滤波器权值w(k+1，j)：
其中，ζ为一很小的正数以防止分母为0，例如取ζ＝10-5；η2为自适应滤波器权值优化的步长；η2的取值范围优选为：η2∈[0.1，1]；
SA7步，采用传递函数为
的滤波器对误差信号e(k)进行滤波得到输出o(k)；
SA8步，判断自适应降噪终止指示变量：若自适应降噪终止指示变量为真，则结束自适应降噪算法，否则k＝k+1，跳至SA2步计算下一时刻自适应降噪算法的输出o(k)。
本发明算法的技术原理是：
图3所示为未采用主副通道信号预处理的自适应降噪算法的原理图，其原理为：假设环境噪声为加性，主通道信号为带噪心音d(k)＝s(k)+n1(k)，其中s(k)为第k时刻纯净的心音信号，而n1(k)为第k时刻主通道所含的环境噪声信号；副通道信号为第k时刻环境噪音n2(k)，其经自适应滤波器加权后得到y(k)用于抵消带噪心音中的环境噪声n1(k)而实现降噪；基于LMS或NLMS算法进行自适应降噪时，其自适应滤波器权值w(k+1，j)的迭代算式分别如式(1)和式(2)所示：
w(k+1，j)＝w(k，j)+η1n2(k-j)[s(k)+n1(k)-y(k)] (1)
其中，e(k)＝d(k)-y(k)＝s(k)+n1(k)-y(k)，w(k+1，j)和w(k，j)分别为第k+1和第k时刻自适应滤波器的权值，n2(k-j)为第k-j时刻的副通道信号，η1和η2均为自适应滤波器权值优化的步长，ζ为一很小的正数以防止分母为0，M为自适应滤波器的阶数。
为方便起见，记自适应滤波器权值向量W(k)＝{w(k，0)，...，w(k，M-1)}。
假设心音信号s(k)和环境噪声n1(k)、n2(k)均为平稳信号，且假设步长η1和η2取值恰当使得自适应滤波器的权值迭代过程收敛；同时假设环境噪声n1(k)、n2(k)的互相关函数不随时间变化，则图3所示系统存在维纳解(最优解)，记为W*＝{w*(0)，...，w*(M-1)}，且权值向量W(k)可收敛到维纳解附近。但根据式(1)和(2)及e(k)的定义，图3所示系统进入稳态后，e(k)和权值向量W(k)依然会随s(k)和n1(k)、n2(k)的随机变化而波动。
进入稳态后e(k)的方差和权值向量偏差的范数(定义为||ΔW||2＝||W*-W(k)||2)均依赖于s(k)的方差。特别是在心音听诊中，第一心音和第二心音的幅值往往远高于环境噪声，因而导致偏差e(k)＝s(k)+n1(k)-y(k)在自适应滤波器权值收敛的过程中周期性跃变，并进而引起自适应滤波器权值周期性过调。减小s(k)的方差有助于降低稳态后e(k)的方差和权值向量偏差的范数||ΔW||2。
中国发明专利申请《用于采集体音信号的双麦克风自适应滤波算法及应用》(公开号：CN109545239A)中所提出的双通道自适应降噪算法原理图如图4所示，其通过在主、副通道引入相同的预处理算法以减小s(k)的方差，从而达到降低稳态后e(k)的方差和权值向量偏差范数||ΔW||2的目的。但该方案中所用的预处理滤波器参数需要人为根据经验选定，且优选地由若干个一阶预加重环节串联而成；如此所得的滤波器效果会由于人为选定的滤波器参数不同而滤波效果有差异，不能确保滤波器取得最佳效果。
在本发明中，预处理滤波器参数通过对s(k)作线性预测编码求得，如此，可有效降低s(k)的方差：
对s(k)作线性预测编码可表示为
其中，参数bi，i∈{1，...，N}，通过最小化式(4)的性能指标J得到
其中，E{}表示求期望值。显然，使J取最小值的参数bi，i∈{1，...，N}，必能保证
此外，根据心音信号的特点，可假设E{s(k)}＝0，故
因此，对于心音降噪而言，E{s2(k)}和/>
即分别为s(k)和/>
的方差；/>
则表明线性预测编码有助于减小信号的方差。注意到
仅当s(k)与s(k-i)，i∈{1，...，N}，均线性无关时等号才成立，而当采样频率高于奈奎斯特频率时，采样时刻相近的心音信号之间存在显著的线性相关性，因此一般而言，/>
即线性预测编码有助于大幅减小信号的方差。
式(3)对应的传递函数即为
为实现上述心音信号双通道自适应降噪算法，本实施例提供一种心音信号双通道自适应降噪装置，包括：
预处理滤波器生成模块，用于在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
以及主副通道信号预处理及降噪模块，用于采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成模块得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)。
本实施例还提供一种存储介质，其中所述存储介质存储有计算机程序，所述计算机程序当被处理器执行时使所述处理器执行实施例一所述的心音信号双通道自适应降噪算法。
本实施例还提供一种计算设备，包括处理器以及用于存储处理器可执行程序的存储器，所述处理器执行存储器存储的程序时，实现实施例一所述的心音信号双通道自适应降噪算法。
下面以具体示例进行说明：
在预处理滤波器生成过程中，选定预处理滤波器的阶数N＝2；在安静的环境下，以Fs＝22050Hz为采样频率，通过主通道获取高信噪比的、持续时间大于1秒且小于3秒的四音律火车头奔马律心音信号d(k)；基于N＝2阶的全极点模型对d(k)作线性预测编码，采用现有的Levinson-Durbin算法求得该全极点模型的传递函数为
则得预处理滤波器传递函数为H(z)＝1-1.4463z-1+0.4517z-2；
主副通道信号预处理及降噪过程包括如下步骤：
SA1步，初始化当前时刻序号k＝0，自适应滤波器权值w(0，j)＝0，j∈{0，...，29}；
SA2步，以Fs＝22050Hz为采样频率，获取当前时刻的主通道带噪心音信号d(k)和副通道环境噪声信号n2(k)；
SA3步，
若k＜2，则取o(k)＝d(k)，副通道预处理后信号
w(k+1，j)＝w(k，j)，并转至SA8步；
若2≤k＜30，则取o(k)＝d(k)，副通道预处理后信号
同时令w(k+1，j)＝w(k，j)，并转至SA8步；
若k≥30，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
SA4步，计算自适应滤波器输出y(k)：
SA5步，计算误差信号e(k)：
SA6步，更新自适应滤波器权值w(k+1，j)：
其中，ζ＝10-5，η2＝0.25；
SA7步，采用传递函数为
的滤波器对e(k)进行滤波得到输出o(k)。
SA8步，判断自适应降噪终止指示变量：若自适应降噪终止指示变量为真，则结束自适应降噪算法；否则k＝k+1，跳至SA2步计算下一时刻自适应降噪算法的输出o(k)。
图5(a)和图5(b)分别为本实施例中在安静环境下所测得的主通道心音信号及其预处理后信号波形；其中，对主通道心音信号进行线性预测编码得
预处理滤波器传递函数为H(z)＝1-1.4463z-1+0.4517z-2；预处理前信号d(k)和预处理后信号/>
的方差分别为3.8724×10-4和2.2780×10-6；可见预处理大大减小了心音信号的方差。
实施例二
本实施例一种心音信号双通道自适应降噪算法，与实施例一的区别在于：本实施例中，在预处理滤波器生成过程中，选定预处理滤波器的阶数N＝5；所得全极点模型的传递函数为
所得预处理滤波器传递函数为H(z)＝1-1.3764z-1+0.4935z-2-0.5846z-3+0.4935z-4-0.0176z-5。
在SA3步中，
若k＜2，则取o(k)＝d(k)，副通道预处理后信号
w(k+1，j)＝w(k，j)，并转至SA8步；
若2≤k＜30，则取o(k)＝d(k)，副通道预处理后信号
同时令w(k+1，j)＝w(k，j)，并转至SA8步；
若k≥30，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
在SA7步中，采用传递函数为
的滤波器对e(k)进行滤波得到输出o(k)。
本实施例的其余步骤与实施例一相同。
图6(a)和图(b)分别为本实施例中在安静环境下所测得的主通道心音信号及其预处理后波形；其中，对主通道心音信号进行线性预测编码得
预处理滤波器传递函数为H(z)＝1-1.3764z-1+0.4935z-2-0.5846z-3+0.4935z-4-0.0176z-5；预处理前信号d(k)和预处理后信号/>
的方差分别为3.8724×10-4和1.7626×10-6；该实施例表明，本发明方法以最小化方差为目标确定最优预处理滤波器参数，不受预处理滤波器阶数影响。
实施例三
本实施例作为实施例一的对比例，目的是比较预处理的作用。本实施例仅实现心音信号双通道自适应降噪算法，与实施例一的区别在于：本实施例中，不包含预处理滤波器生成过程及在主副通道信号预处理及降噪过程中不对主副通道信号进行预处理；
在SA3步中，不管k取何值，均令
在SA7步中，采用传递函数为
的滤波器对e(k)进行滤波得到输出o(k)。
本实施例的其余步骤与实施例一相同。
图7(a)至图7(c)分别为实施例一中，带噪心音信号、自适应滤波后心音信号和相邻两时刻自适应滤波器权值向量差的范数||ΔW(k)||2＝||W(k+1)-W(k)||2随时间变化的曲线；由于维纳解W*通常未知，故以||ΔW(k)||2随时间的波动近似衡量||ΔW||2＝||W*-W(k)||2随时间的波动。
图8(a)至图8(c)分别为实施例三中，带噪心音信号、自适应滤波后心音信号和相邻两时刻自适应滤波器权值向量差的范数||ΔW(k)||2＝||W(k+1)-W(k)||2随时间变化的曲线。
图7中||ΔW(k)||2随时间变化相对平缓，其方差为8.8756×10-7，而图8中||ΔW(k)||2会随第一第二心音周期性过调，出现短时大幅跃变，其方差为3.2938×10-6。虽然图8中||ΔW(k)||2跃变后依然能随时间的推移而重新收敛至0附近，但在其收敛过程中，自适应降噪效果不佳。
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。
Claims (9)
1.一种心音信号双通道自适应降噪算法，其特征在于：经主通道采集带噪心音信号，以及副通道采集环境噪声；包括预处理滤波器生成过程和主副通道信号预处理及降噪过程；
所述预处理滤波器生成过程，是指：在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
所述主副通道信号预处理及降噪过程，是指：采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成过程得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)；
所述预处理滤波器生成过程，包括如下步骤：
SG1步，设定预处理滤波器的阶数N，其中N≥2，预处理滤波器的传递函数为H(z)＝1+a1z-1+…+aNz-N，其中ai，i∈{1，...，N}，为待定参数；
SG2步，在安静的环境下，以采样频率Fs测量主通道信号d′(k)，主通道信号d′(k)持续时间超过一次心动周期；
SG3步，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数
SG4步，将bi，i∈{1，...，N}，代入到ai中，得到预处理滤波器的传递函数H(z)；
所述主副通道信号预处理及降噪过程，包括如下步骤：
SA1步，初始化当前时刻序号k＝0，自适应滤波器权值w(0，j)＝0，j∈{0，...，M-1}，其中M为自适应滤波器阶数，M＞N；
SA2步，以采样频率Fs分别测量当前时刻主通道信号d(k)和副通道信号n2(k)；
SA3步，计算当前时刻的预处理结果
和/>
若k＜N，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若N≤k＜M，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若k≥M，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
SA4步，计算自适应滤波器输出y(k)：
SA5步，计算误差信号e(k)：
SA6步，更新自适应滤波器权值w(k+1，j)：
其中，ζ为一很小的正数以防止分母为0；η2为自适应滤波器权值优化的步长；
SA7步，采用传递函数为
的滤波器对误差信号e(k)进行滤波得到输出o(k)；
SA8步，判断自适应降噪终止指示变量：若自适应降噪终止指示变量为真，则结束自适应降噪算法，否则k＝k+1，跳至SA2步计算下一时刻自适应降噪算法的输出o(k)。
2.根据权利要求1所述的心音信号双通道自适应降噪算法，其特征在于：所述自适应滤波器阶数M的取值范围为：M∈[20，300]。
3.根据权利要求1所述的心音信号双通道自适应降噪算法，其特征在于：所述自适应滤波器权值优化的步长η2的取值范围为：η2∈[0.1，1]。
4.根据权利要求1所述的心音信号双通道自适应降噪算法，其特征在于：所述采样频率Fs≥4KHz。
5.根据权利要求1所述的心音信号双通道自适应降噪算法，其特征在于：所述SG2步中，主通道信号d′(k)持续时间的取值范围为[1s，3s]。
6.根据权利要求1所述的心音信号双通道自适应降噪算法，其特征在于：所述SG3步的线性预测编码中，采用Levinson-Durbin算法计算出全极点模型的传递函数G(z)。
7.一种心音信号双通道自适应降噪装置，其特征在于，包括：
预处理滤波器生成模块，用于在安静的环境下采样主通道信号d′(k)，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数G(z)，从而得到预处理滤波器的传递函数H(z)＝1/G(z)；
以及主副通道信号预处理及降噪模块，用于采样主通道信号d(k)和副通道信号n2(k)；对主通道信号d(k)和副通道信号n2(k)进行预处理；将主通道信号d(k)和副通道信号n2(k)预处理后的信号进行自适应滤波计算，得到误差信号e(k)；之后采用预处理滤波器生成模块得到的预处理滤波器传递函数的倒数1/H(z)对误差信号e(k)进行滤波处理，得到输出信号o(k)；
所述预处理滤波器生成模块，包括如下步骤：
SG1步，设定预处理滤波器的阶数N，其中N≥2，预处理滤波器的传递函数为H(z)＝1+a1z-1+…+aNz-N，其中ai，i∈{1，...，N}，为待定参数；
SG2步，在安静的环境下，以采样频率Fs测量主通道信号d′(k)，主通道信号d′(k)持续时间超过一次心动周期；
SG3步，采用N阶的全极点模型对主通道信号d′(k)作线性预测编码，求得全极点模型的传递函数
SG4步，将bi，i∈{1，...，N}，代入到ai中，得到预处理滤波器的传递函数H(z)；
所述主副通道信号预处理及降噪模块，包括如下步骤：
SA1步，初始化当前时刻序号k＝0，自适应滤波器权值w(0，j)＝0，j∈{0，...，M-1}，其中M为自适应滤波器阶数，M＞N；
SA2步，以采样频率Fs分别测量当前时刻主通道信号d(k)和副通道信号n2(k)；
SA3步，计算当前时刻的预处理结果
和/>
若k＜N，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若N≤k＜M，则主通道预处理后信号
副通道预处理后信号/>
自适应滤波器权值w(k+1，j)＝w(k，j)，输出信号o(k)＝d(k)，并转至SA8步；
若k≥M，则主通道预处理后信号
副通道预处理后信号/>
并转至SA4步；
SA4步，计算自适应滤波器输出y(k)：
SA5步，计算误差信号e(k)：
SA6步，更新自适应滤波器权值w(k+1，j)：
其中，ζ为一很小的正数以防止分母为0；η2为自适应滤波器权值优化的步长；
SA7步，采用传递函数为
的滤波器对误差信号e(k)进行滤波得到输出o(k)；
SA8步，判断自适应降噪终止指示变量：若自适应降噪终止指示变量为真，则结束自适应降噪算法，否则k＝k+1，跳至SA2步计算下一时刻自适应降噪算法的输出o(k)。
8.一种计算机可读存储介质，其特征在于，其中所述计算机可读存储介质存储有计算机程序，所述计算机程序当被处理器执行时使所述处理器执行权利要求1-6中任一项所述的心音信号双通道自适应降噪算法。
9.一种计算机设备，包括处理器以及用于存储处理器可执行程序的存储器，其特征在于，所述处理器执行存储器存储的程序时，实现权利要求1-6中任一项所述的心音信号双通道自适应降噪算法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011278315.4A CN112562701B (zh) | 2020-11-16 | 2020-11-16 | 心音信号双通道自适应降噪算法、装置、介质及设备 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011278315.4A CN112562701B (zh) | 2020-11-16 | 2020-11-16 | 心音信号双通道自适应降噪算法、装置、介质及设备 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112562701A CN112562701A (zh) | 2021-03-26 |
| CN112562701B true CN112562701B (zh) | 2023-03-28 |
Family
ID=75042480
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011278315.4A Active CN112562701B (zh) | 2020-11-16 | 2020-11-16 | 心音信号双通道自适应降噪算法、装置、介质及设备 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112562701B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113225046B (zh) * | 2021-04-14 | 2022-12-20 | 华南理工大学 | 一种用于数字听诊器的双通道自适应滤波器模型定阶方法 |
Family Cites Families (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5148488A (en) * | 1989-11-17 | 1992-09-15 | Nynex Corporation | Method and filter for enhancing a noisy speech signal |
| US6122610A (en) * | 1998-09-23 | 2000-09-19 | Verance Corporation | Noise suppression for low bitrate speech coder |
| US7139701B2 (en) * | 2004-06-30 | 2006-11-21 | Motorola, Inc. | Method for detecting and attenuating inhalation noise in a communication system |
| CN1815552B (zh) * | 2006-02-28 | 2010-05-12 | 安徽中科大讯飞信息科技有限公司 | 基于线谱频率及其阶间差分参数的频谱建模与语音增强方法 |
| EP2363853A1 (en) * | 2010-03-04 | 2011-09-07 | Österreichische Akademie der Wissenschaften | A method for estimating the clean spectrum of a signal |
| CN102429662B (zh) * | 2011-11-10 | 2014-04-09 | 大连理工大学 | 家庭环境中睡眠呼吸暂停综合征的筛查系统 |
| FR3023036A1 (fr) * | 2014-06-27 | 2016-01-01 | Orange | Re-echantillonnage par interpolation d'un signal audio pour un codage / decodage a bas retard |
| CN104157293B (zh) * | 2014-08-28 | 2017-04-05 | 福建师范大学福清分校 | 一种增强声环境中目标语音信号拾取的信号处理方法 |
| CN109545239B (zh) * | 2018-12-06 | 2021-11-05 | 华南理工大学 | 用于采集体音信号的双麦克风自适应滤波算法及应用 |
| CN111814515B (zh) * | 2019-04-11 | 2024-02-23 | 哈尔滨工业大学 | 基于改进变步长lms自适应的有源噪声对消方法 |
| CN110085247B (zh) * | 2019-05-06 | 2021-04-20 | 上海互问信息科技有限公司 | 一种针对复杂噪声环境的双麦克风降噪方法 |
-
2020
- 2020-11-16 CN CN202011278315.4A patent/CN112562701B/zh active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN112562701A (zh) | 2021-03-26 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109817209B (zh) | 一种基于双麦克风阵列的智能语音交互系统 | |
| KR102469516B1 (ko) | 마이크로폰 어레이 기반 타겟 음성 획득 방법 및 장치 | |
| CN109767783B (zh) | 语音增强方法、装置、设备及存储介质 | |
| CN109545239B (zh) | 用于采集体音信号的双麦克风自适应滤波算法及应用 | |
| EP2600344B1 (en) | Multi-input noise suppresion device, multi-input noise suppression method, program, and integrated circuit | |
| JP5153886B2 (ja) | 雑音抑圧装置および音声復号化装置 | |
| KR20220062598A (ko) | 오디오 신호 생성을 위한 시스템 및 방법 | |
| Wood | Motion artifact reduction for wearable photoplethysmogram sensors using micro accelerometers and Laguerre series adaptive filters | |
| JP6374120B2 (ja) | 発話の復元のためのシステムおよび方法 | |
| CN106999072B (zh) | 利用倒频谱平滑化和基于质量的动态信道选择的多信道心冲击描记器 | |
| CN112562701B (zh) | 心音信号双通道自适应降噪算法、装置、介质及设备 | |
| CN110349598A (zh) | 一种低信噪比环境下的端点检测方法 | |
| CN111599372B (zh) | 一种稳定的在线多通道语音去混响方法及系统 | |
| CN109905793A (zh) | 一种风噪声抑制方法及装置 | |
| CN111616695B (zh) | 一种心率获取方法、装置、系统和介质 | |
| CN108652609B (zh) | 一种心率获取方法、系统及可穿戴式设备 | |
| CN110111804B (zh) | 基于rls算法的自适应去混响方法 | |
| CN111883153B (zh) | 一种基于麦克风阵列的双端讲话状态检测方法及装置 | |
| Seltzer | Bridging the gap: Towards a unified framework for hands-free speech recognition using microphone arrays | |
| Paul et al. | Noise reduction for heart sounds using a modified minimum-mean squared error estimator with ECG gating | |
| Hong et al. | Dual-microphone noise reduction in car environments with determinant analysis of input correlation matrix | |
| CN114566179A (zh) | 一种时延可控的语音降噪方法 | |
| Fras et al. | Convolutive Weighted Multichannel Wiener Filter Front-end for Distant Automatic Speech Recognition in Reverberant Multispeaker Scenarios. | |
| Poleshenkov et al. | Research on dependences of speech pitch parameters on pulse and heartbeat signals | |
| CN114694675B (zh) | 一种基于麦克风阵列的广义旁瓣对消器及后置滤波算法 |
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 |