CN111643068B - 基于emd及其能量的心电信号去噪算法、设备和存储介质 - Google Patents

Abstract

本发明公开了基于EMD及其能量的心电信号去噪算法、设备和存储介质，该算法包括步骤：信号均值化处理及公差值的确定；对信号进行EMD分解，得到各阶IMF；根据EMD分解后各阶IMF的能量，计算能量曲线除边界点以外的第一个极大点和第一个极小点的所在阶数，结合所述公差值，判断需要去噪的IMF阶数变化点，若所述能量曲线的第一个极大点和第一个极小点均不存在，则预设所述IMF阶数变化点；对所述IMF阶数变化点之前的各阶IMF进行阈值去噪；将阈值去噪后的各阶IMF与未进行阈值去噪的IMF以及残差进行重构，生成去噪后的心电信号。本方案能够完整地保留大于阈值部分的信息，并且不会产生附加震荡和跳跃点，较好地保证了原始信号的平滑性，提高了去噪质量。

Description

基于EMD及其能量的心电信号去噪算法、设备和存储介质

技术领域

本发明涉及心电信号去噪技术领域，尤其是涉及一种基于EMD及其能量的心电信号去噪算法、设备和存储介质。

背景技术

心电信号(ECG)是人类较早研究并应用于医学临床的一种生理电信号。心电图信号在医学上被广泛运用于心脏病诊断。然而由于如人体的呼吸、运动、以及周边环境中的电磁干扰等等，心电图往往会包含大量的噪声。常见的噪声有工频干扰、肌电干扰和基线漂移。所有这些噪声都可能引起心电波形的畸变，掩盖对诊断有重要意义的微小特征。因此，消除心电信号中的噪声成为必要。

由于心电信号噪声源的种类复杂，对其进行去噪处理成为一个难题，在心电信号的自动处理中，这一困难更加突出。基于经验模式分解(EMD)的新信号分析方法，是将ECG信号分解为“n”个IMF分量和一个残差。低阶IMF分量往往对应于高频噪声。在传统的基于EMD的ECG去噪方法中，为了去除高频噪声分量，通常直接丢弃低阶IMF，然而低阶IMF中往往还包含了QRS波群等信息，因此导致了QRS波群等心电信息的丢失。

在去噪方法上，常用的有硬阈值法和软阈值法两种。传统硬阈值法在均方误差意义上优于软阈值法，但是硬阈值法会产生附加震荡和跳跃点，难以具备原始信号的平滑性。软阈值法得到的去噪后信号相对平滑得多，但是会压缩信号，产生一定的偏差，直接影响到重构的信号与真实信号的逼近程度。

发明内容

本发明旨在至少解决现有技术中存在的技术问题之一，第一方面，提出一种基于EMD及其能量的心电信号去噪算法，根据EMD分解后的各阶IMF能量，确定需要去噪的IMF阶数，与传统的直接丢弃低阶IMF的方法相比，有利于提高去噪质量。另外，对阈值去噪方法进行了改进，能有效去除噪声，改善信噪比，较为完整地保留原始信号特征。

根据本发明第一方面实施例的基于EMD及其能量的心电信号去噪算法，包括步骤：

S1、信号均值化处理及公差值的确定；

S2、对均值化处理后的信号进行EMD分解，得到各阶IMF；

S3、根据EMD分解后各阶IMF的能量，得到对应的能量曲线，计算该能量曲线除边界点以外的第一个极大点和第一个极小点的所在阶数，根据所述第一个极大点和第一个极小点的所在阶数及其与所述公差值的大小关系，判断需要去噪的IMF阶数变化点，若所述能量曲线的第一个极大点和第一个极小点均不存在，则预设所述IMF阶数变化点；

S4、对所述IMF阶数变化点之前的各阶IMF进行阈值去噪，所述IMF阶数变化点所在阶数的IMF及其之后的各阶IMF未进行阈值去噪；所述阈值去噪的方法如下：

其中d是需要进行阈值处理的系数，T是阈值，D是阈值处理后的系数；

S5、将阈值去噪后的各阶IMF与未进行阈值去噪的IMF以及残差进行重构，生成去噪后的心电信号。

根据本发明的一些实施例，设待求的所述IMF阶数变化点所在阶数为var_pots，所述能量曲线的第一个极大点所在阶数为max_pots(1)，第一个极小点所在阶数为min_pots(1)，步骤S3中，所述IMF阶数变化点采用以下方法确定：

S3a、当所述第一个极大点和第一个极小点同时存在时：

当max_pots(1)＜min_pots(1)：若max_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

当max_pots(1)＞min_pots(1)时：若min_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

S3b、当所述第一个极小点存在，但所述第一个极大点不存在时：若min_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

S3c、当所述第一个极大点存在，但所述第一个极小点不存在时：若max_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

S3d、当所述第一个极大点和第一个极小点均不存在时，定义为：var_pots＝3。

根据本发明的一些实施例，所述公差值不小于3。

根据本发明的一些实施例，所述公差值为3～6。

第二方面，提供一种基于EMD及其能量的心电信号去噪设备，包括至少一个处理器和与所述至少一个处理器通信连接的存储器，所述存储器上存储有可被所述至少一个处理器执行的指令，所述指令被所述至少一个处理器执行，以使所述至少一个处理器能够执行上述的基于EMD及其能量的心电信号去噪算法。

第三方面，提供一种计算机可读存储介质，计算机可读存储介质存储有计算机可执行指令，所述计算机可执行指令用于使计算机执行如上述的基于EMD及其能量的心电信号去噪算法。

第四方面，提供一种计算机程序产品，所述计算机程序产品包括存储在计算机可读存储介质上的计算机程序，所述计算机程序包括程序指令，当所述程序指令被计算机执行时，使计算机执行如上所述的基于EMD及其能量的心电信号去噪算法。

本发明实施例提供的一个或多个技术方案，至少具有如下有益效果：

根据EMD分解后的各阶IMF能量，确定需要去噪的IMF阶数，与传统的直接丢弃低阶IMF的方法相比，保留了QRS波群，并进行了有效去噪。阈值去噪过程中采用了改进阈值函数进行去噪，克服了采用传统软阈值函数及硬阈值函数的缺点，能够完整地保留大于阈值部分的信息，并且不会产生附加震荡和跳跃点，较好地保证了原始信号的平滑性，提高了去噪质量。

附图说明

图1是本发明实施例的基于EMD及其能量的心电信号去噪算法流程图。

图2是采用不同阈值处理方法对MIT-BIH数据库中114号记录进行去噪的效果对比图。

图3是图2中部分样本点的放大图。

图4是本发明实施例的基于EMD及其能量的心电信号去噪设备示意图。

附图标记：

101-处理器，102-存储器。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。

参照图1，根据本发明实施例的基于EMD及其能量的心电信号去噪算法，包括步骤：

S1、信号均值化处理及公差值的确定。

S2、对均值化处理后的信号进行EMD分解，得到各阶IMF。

S3、根据EMD分解后各阶IMF的能量，得到对应的能量曲线，计算该能量曲线除边界点以外的第一个极大点和第一个极小点的所在阶数，根据所述第一个极大点和第一个极小点的所在阶数及其与所述公差值的大小关系，判断需要去噪的IMF阶数变化点，若所述能量曲线的第一个极大点和第一个极小点均不存在，则预设该IMF阶数变化点。

S4、对IMF阶数变化点之前的各阶IMF进行阈值去噪，IMF阶数变化点所在阶数的IMF及其之后的各阶IMF未进行阈值去噪，阈值去噪的方法如下：

其中d是需要进行阈值处理的系数，T是阈值，D是阈值处理后的系数。

S5、将阈值去噪后的各阶IMF与未进行阈值去噪的IMF以及残差进行重构，生成去噪后的心电信号。

在对信号进行EMD分解之前，首先对信号进行零均值化处理，目的是去掉信号中的直流成分，使得零频信号幅值很小，在后续的功率谱分析中可以减少直流能量对低中频信号谱分析的影响。零均值化处理方法为：对输入信号X，进行X＝X-mean(X),mean(X)代表X的均值。

EMD分解后，根据求出的各阶IMF，计算各阶IMF分量包含的能量大小，计算公式为：E_j＝log_e((median(|IMF_j|)/0.6745)²)，其中E_j表示第j阶的IMF能量，|IMF_j|表示对第j阶IMF分量取绝对值，median()为对元素进行排序并取中值的函数。

根据计算出来的各阶IMF能量，可以得到对应的能量曲线，计算该能量曲线除边界点以外的第一个极大点和第一个极小点的所在阶数，并结合步骤S1确定的公差值，来判断需要去噪的IMF阶数变化点。

设待求的IMF阶数变化点所在阶数为var_pots，所述能量曲线的第一个极大点所在阶数为max_pots(1)，第一个极小点所在阶数为min_pots(1)，上述IMF阶数变化点可采用以下方法确定：

a、当第一个极大点和第一个极小点同时存在时：

当max_pots(1)＜min_pots(1)：若max_pots(1)大于公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

当max_pots(1)＞min_pots(1)时：若min_pots(1)大于公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

b、当第一个极小点存在，但所述第一个极大点不存在时：若min_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

c、当第一个极大点存在，但所述第一个极小点不存在时：若max_pots(1)大于公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

d、当第一个极大点和第一个极小点均不存在时，定义为：var_pots＝3。

公差值用于在后续的处理中判断需要去噪的IMF阶数。发明人发现，本实施例公差值不宜小于3，因为当公差值过小时，得到的去噪IMF过少，去噪效果不明显。对于一般的信号，公差值取3～6即可。当信号的噪声水平非常高时，可以适当增加公差值。

在阈值处理方法上，传统硬阈值法的原理是设定一个固定阈值，大于该阈值的系数保留不变，小于该阈值的系数直接置零，这种处理方法能够完整地保留大于阈值部分的信息，但信号会产生附加震荡和跳跃点，不具有原始信号的平滑性。软阈值法提出把大于阈值的系数进行相应的收缩，小于阈值的系数置零，处理后的效果是:在阈值附近的系数是连续的，得到的去噪后信号是平滑的，但是由于会压缩信号，会产生一定的偏差，直接影响到重构的信号与真实信号的逼近程度。

为了尽量完整保留原始信号特征，提出了一种新的阈值处理方法，由前述的阈值处理函数可以看出：

当d→±∞时，有

当d→±T时，有

即当d绝对值不断增大时，阈值处理后的系数将不断接近处理前的系数，当d绝对值趋近阈值时处理后的系数接近于0，说明该阈值处理函数在阈值点处是连续的，克服了采用传统软阈值函数处理后的系数与处理前系数存在较大差值的问题以及传统硬阈值函数不连续的缺点，能够完整地保留大于阈值部分的信息，并且不会产生附加震荡和跳跃点，得到的去噪后信号是平滑的，提高了去噪质量。

其中，阈值T的计算公式为本领域公知，计算公式为：

其中N为样本数，σ为噪声标准差，σ根据各阶IMF分别计算出来。

为证明本实施例的去噪效果，选取MIT-BIH数据库中114号记录进行去噪，去噪过程中，设定公差值为5，去噪效果如图2，图中“新阈值”代表本方案提出的阈值去噪的方法。将步骤S4的阈值去噪方法分别替换为传统软阈值、硬阈值方法，其余步骤不变，去噪效果比较见图2，图3是图2中信号样本点位于2650-2850的放大图。

可以看出，基于本实施例提出的去噪方法，采用改进阈值函数去噪后的信号，较为完整的保留了原始信号的波点、高峰等特征，去噪效果也较好；采用传统硬阈值去噪后的信号出现了更多的震荡，传统软阈值去噪后的信号则损失了部分高频信号，详见图3中椭圆形虚线框标记的部分。

另外对添加了20db白高斯噪声的MIT-BIH数据库的部分信号进行去噪，利用改善信噪比(SNR_imp)分析去噪信号的性能，研究重构信号的质量。改善信噪比的公式定义如下：

其中x(i)为原始信号，y(i)为添加噪声后的信号，s(i)为最终去噪并平滑后的信号。

将本实施例步骤S4的阈值去噪方法分别替换为传统软阈值、硬阈值方法，其余步骤不变，去噪效果如表1所示，可以看出，本方案改善信噪比较高。

表1

MIT-BIH	软阈值	硬阈值	本实施例
				106	1.1062	4.3969	4.7251
113	2.3292	5.3109	5.3258
				119	4.0147	5.2828	5.2915
平均值	2.4833	4.9968	5.1141

本发明实施例还提供基于EMD及其能量的心电信号去噪设备，包括至少一个处理器101，以及存储器102，图4以一个处理器101为例。处理器101和存储器102可以通过总线或者其他方式连接，图4以通过总线连接为例。

存储器102作为一种非暂态计算机可读存储介质，可用于存储计算机可执行程序、指令以及模块，如本发明实施例中的基于EMD及其能量的心电信号去噪算法对应的程序、指令以及模块。处理器101通过运行存储在存储器102中的程序、指令以及模块，实现上述方法实施例的基于EMD及其能量的心电信号去噪算法，例如，执行以上描述的图1中的方法步骤S1至S5。

存储器102可以包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能实现所需要的应用程序；存储数据区可存储所需要或所创建的数据等。存储器102可以包括高速随机存取存储器，还可以包括非暂态存储器，例如至少一个磁盘存储器件、闪存器件、或其他非暂态固态存储器件。

本发明实施例还提供了一种计算机可读存储介质，计算机可读存储介质存储有计算机可执行指令，该计算机可执行指令被一个或多个处理器执行，例如，被图4中的一个处理器101执行，从而执行上述方法实施例中的基于EMD及其能量的心电信号去噪算法，例如，执行以上描述的图1中的方法步骤S1至S5。

通过以上的实施方式的描述，本领域技术人员可以清楚地了解到各实施方式可借助软件加通用硬件平台的方式来实现。本领域技术人员可以理解实现上述实施例方法中的全部或部分流程是可以通过计算机程序来指令相关的硬件完成，程序可存储于一计算机可读取存储介质中，该程序在执行时，可包括如上述方法的实施例的流程。

尽管已经示出和描述了本发明的实施例，但是本发明不限于上述实施例，本领域的普通技术人员可以理解：在不脱离本发明的原理和宗旨的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由权利要求及其等同物限定。

Claims (5)

1.基于EMD及其能量的心电信号去噪算法，其特征在于，包括步骤：

S1、信号均值化处理及公差值的确定；

S2、对均值化处理后的信号进行EMD分解，得到各阶IMF；

S3、根据EMD分解后各阶IMF的能量，得到对应的能量曲线，计算该能量曲线除边界点以外的第一个极大点和第一个极小点的所在阶数，根据所述第一个极大点和第一个极小点的所在阶数及其与所述公差值的大小关系，判断需要去噪的IMF阶数变化点，若所述能量曲线的第一个极大点和第一个极小点均不存在，则预设所述IMF阶数变化点，其中，设待求的所述IMF阶数变化点所在阶数为var_pots，所述能量曲线的第一个极大点所在阶数为max_pots(1)，第一个极小点所在阶数为min_pots(1)，步骤S3中，所述IMF阶数变化点采用以下方法确定：

S3a、当所述第一个极大点和第一个极小点同时存在时：

当max_pots(1)＜min_pots(1)：若max_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

当max_pots(1)＞min_pots(1)时：若min_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

S3b、当所述第一个极小点存在，但所述第一个极大点不存在时：若min_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝min_pots(1)；或者，

S3c、当所述第一个极大点存在，但所述第一个极小点不存在时：若max_pots(1)大于所述公差值，计算公式为：

否则，计算公式为：var_pots＝max_pots(1)；或者，

S3d、当所述第一个极大点和第一个极小点均不存在时，定义为：var_pots＝3；

计算各阶IMF分量包含的能量大小为：

其中，E_j表示第j阶的IMF能量，|IMF_j|表示对第j阶IMF分量取绝对值，median()为对元素进行排序并取中值的函数；

S4、对所述IMF阶数变化点之前的各阶IMF进行阈值去噪，所述IMF阶数变化点所在阶数的IMF及其之后的各阶IMF未进行阈值去噪；所述阈值去噪的方法如下：

当d→±∞时，有

当d→±T时，有

其中d是需要进行阈值处理的系数，T是阈值，D是阈值处理后的系数；

S5、将阈值去噪后的各阶IMF与未进行阈值去噪的IMF以及残差进行重构，生成去噪后的心电信号。

2.根据权利要求1所述的基于EMD及其能量的心电信号去噪算法，其特征在于，所述公差值不小于3。

3.根据权利要求2所述的基于EMD及其能量的心电信号去噪算法，其特征在于，所述公差值为3～6。

4.一种基于EMD及其能量的心电信号去噪设备，其特征在于，包括至少一个处理器和与所述至少一个处理器通信连接的存储器，所述存储器上存储有可被所述至少一个处理器执行的指令，所述指令被所述至少一个处理器执行，以使所述至少一个处理器能够执行如权利要求1至3任一项所述的基于EMD及其能量的心电信号去噪算法。

5.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有计算机可执行指令，所述计算机可执行指令用于使计算机执行如权利要求1至3任一项所述的基于EMD及其能量的心电信号去噪算法。

Images