CN106798554A - 一种含噪imf分量及心电信号的去噪方法 - Google Patents

Abstract

本发明涉及一种含噪IMF分量及心电信号的去噪方法，将EEMD算法和新小波阈值去噪有机地结合在一起，相对于传统的小波去噪算法，克服了软阈值函数和硬阈值函数各自的缺点，并对传统阈值进行了改进，增加了阈值使用的灵活性和实用性。另外，采用了双阈值和单阈值的组合形式，实现了计算量的减少，效率的提高。

Description

一种含噪IMF分量及心电信号的去噪方法

技术领域

本发明属于心电信号处理技术领域，更具体的说，涉及一种基于EEMD与小波新阈值的含噪IMF分量心电信号去噪方法。

背景技术

在心电去噪方面，有着很多的方法，诸如：传统数字滤波器FIR(有限长单位冲激响应)滤波器和IIR(无限脉冲响应)滤波器、自适应滤波、中值滤波、小波变换、数学形态法、EEMD等等。在这些方法中，EEMD和小波变换的去噪效果较好，它们的原理都是通过分解-重构的方式来进行去噪处理的。小波变换的原理是：通过选定的小波基，利用Mallet算法对原始信号进行多尺度分解，得到各尺度的近似系数和小波系数，然后对含噪尺度中的小波系数进行阈值处理，之后再重构，即可达到去噪的目的。EEMD的原理是：根据原始信号本身进行分解，得到有限个IMF分量和一个趋势项，根据噪声的统计特性确定含噪的IMF分量并将其去除后重构，即可达到去噪的目的。但是小波变换存在着小波基的选取问题，所以为了克服小波基的选取问题。有人设计了基于EEMD与小波软阈值的心电信号去噪算法，将EEMD分解出来的含噪IMF分量进行小波软阈值处理后再重构，在心电去噪上取得了较好的效果。

现有的阈值函数主要由两种：硬阈值函数和软阈值函数，但它们各自都存在着相应的缺点。硬阈值的缺点是：在阈值点处不连续，这有可能会使重构的心电信号产生震荡；软阈值函数的缺点是：当小波系数大于阈值时，处理后的小波系数与原小波系数之间存在恒定的偏差，这在一定程度上会影响心电信号的精确度。

目前使用较多的有VisuShrink阈值、Rigrsure阈值、Sqtwolog阈值、Heursure阈值和Minimax阈值。在这些阈值中，VisuShrink阈值采用的是一种全局统一阈值，即该阈值应用于小波分解后各尺度的小波系数，因而误差会比较大，去噪效果会比较差；Rigrsure阈值是根据Stein无偏似然估计(SURE)的阈值，该阈值是专门针对软阈值函数得出的结论，故存在软阈值函数的缺陷，而且对于高信噪比的信号，该阈值抑制噪声的效果不明显；Sqtwolog阈值是一种固定的阈值，它是根据分解信号的长度来确定阈值，灵活性小且不适合信噪比高的信号；Heursure阈值又称启发式SURE阈值，它是结合Rigrsure阈值和Sqtwolog阈值，需要进行参数的对比，计算量较大、运算时间较长；Minimax阈值为最大最小阈值，它的准则是使所选用的阈值能够产生最小的极大方差，是一种固定阈值，该阈值的灵活性较小。

综上所述，在传统阈值的选取方面，都只选取了一个阈值作为参考，这样使得某些含有信息量的小波系数有可能因为低于给定阈值而被误认为是噪声而被去除，从而加大了噪声误检的可能性。同时，传统的阈值选取方式，包括上述的全局阈值和局部适应阈值都没有随着分解尺度的变化而有所改变，缺乏灵活性。

发明内容

本发明要解决的技术问题是针对ECG(心电图)信号的噪声，提出了一种基于EEMD(集合经验模分解)与小波新阈值的心电信号去噪方法。首先根据EEMD分解理论对心电信号进行分解，得到有限个IMF(本征模函数)分量和一个趋势项；其次根据噪声统计特性，确定出含噪的IMF分量；然后利用小波新阈值对含噪IMF分量进行降噪处理，最后将处理后的IMF分量与其余信号IMF分量进行重构，获得去噪的心电信号。

本发明的目的是为了克服硬阈值函数、软阈值函数的缺点，以及阈值选取问题，采用了两种不同的阈值组合方式对小波分解后的不同尺度做相应的去噪处理。(1)用双阈值的方式处理低尺度(如第一尺度)的噪声，能够减小噪声误检的可能性，而且方法简单，计算量较小；(2)用本发明采用新的单阈值函数和改进的阈值处理较低尺度(如第三尺度)的噪声，该阈值函数不仅能够克服硬阈值函数和软阈值函数的缺点，而且改进后的阈值能够随着分解尺度的变化而改变，增加了阈值的实用性和灵活性。

本发明解决上述问题的技术方案为：提供一种含噪IMF分量的去噪方法，包括如下步骤：

S10、对第一个含噪的IMF分量进行小波分解，分解N层，N为自然数；

S20、前K个层次的小波系数利用双阈值方法处理，K为自然数且K＜N；

S30、第K+1到N层中的小波系数利用单阈值方法处理；

S40、阈值处理后进行小波重构，得到第一个去噪的IMF分量；

S50、重复上述S1～S4四个步骤，将所有含噪的IMF分量都进行阈值去噪处理。

在本发明提供的含噪IMF分量的去噪方法中，所述步骤S20包括如下步骤：

S21、小波分解后得到第i层的小波系数di，i≤K且i为自然数；

S22、求取第i层小波系数di的平均值g，计算公式如下：

其中，n为信号长度，di为小波系数；

S23、根据该平均值计算出第i层小波系数中的上阈值Ht、下阈值Lt，计算公式如下：

Ht＝g+(Max-g)·θ

Lt＝g-(g-Min)·θ，其中，Max和Min分别为第i层小波系数中的最大值和最小值，θ为上下阈值选取参数，其取值范围为(0,1)；

S24、通过上阈值Ht、下阈值Lt对di的小波系数进行处理，保留上下阈值之间的小波系数，达到去噪的目的；

S25、重复上述S21～S24四个步骤，直到前K个尺度都进行双阈值去噪处理。

在本发明提供的含噪IMF分量的去噪方法中，所述步骤S30包括如下步骤：

S31、小波分解后得到第i层的小波系数di，K+1≤i≤N且i为自然数；

S32、使用如下阈值函数处理第i层小波系数di：

其中，sgn()为符号函数，d_i为小波系数，α和β为调整阈值函数结构的参数；λ为改进阈值，计算公式为：

其中，N为信号长度，j为分解尺度，σ为高斯白噪声的标准差，

S33、重复上述S21和S21两个步骤，直到第K+1～N个尺度都进行单阈值去噪处理。

本发明还提供一种心电信号去噪方法，包括如下步骤：

S100、将含噪心电信号进行EEMD分解；

S200、根据噪声统计特性判定含噪IMF和信号IMF的分界点序号；

S300、使用如前述任意一种方法处理所述含噪IMF分量，得到小波阈值处理后的IMF分量；

S400、将所述小波阈值处理后的IMF分量与所述信号IMF分量重构，得到去噪的心电信号。

在本发明提供的心电信号去噪方法中，所述步骤S100包括如下步骤：

S101、将信号f加入均值为零、方差为常数的高斯白噪声s(t)，得到一个新的信号X(t)，即X(t)＝f(t)+s(t)；

S102、将信号X(t)进行EMD分解为有限个IMF分量和趋势项，即其中j表示第j个IMF分量；

S103、每次加入不同的白噪声，重复上述S101和S102两个步骤，得到其中i表示第i次加入白噪声；

S104、将上述结果进行总体平均运算，用来消除多次加入白噪声对真实IMF的影响，得到IMF分量，即其中s表示IMF分量的个数。

在本发明提供的心电信号去噪方法中，所述步骤S200包括如下步骤：

S201、计算各阶IMF与原始心电信号x(t_i)的相关系数R(x(t_i)，IMF_j(t_i))；它由下式确定：

其中，cov()代表协方差；

S202、含噪IMF和信号IMF的分界点序号k由下式确定：

其中，firstlocal表示第一个局部极小值点；

S203、根据分界点k，将其前面视为含噪IMF分量，其后面视为信号IMF分量。

实施本发明，具有如下有益效果：方法简单，易于实现，将EEMD算法和新小波阈值去噪有机的结合在一起，相对于传统的小波去噪算法，克服了软阈值函数和硬阈值函数各自的缺点，并对传统阈值进行了改进，增加了阈值使用的灵活性和实用性。另外还采用了双阈值和单阈值组合的形式，从而减小了计算量，提高了效率。

附图说明

下面将结合附图及实施例对本发明作进一步说明，附图中：

图1为本发明心电信号去噪方法较佳实施例的流程图；

图2为本发明含噪IMF分量的去噪方法较佳实施例的流程图；

图3为双阈值方法去噪的流程图；

图4为单阈值去噪方法的流程图；

图5为本发明心电信号去噪方法实验仿真结果图。

具体实施方式

下面将结合附图对本发明的实施例进行具体描述。

图1为本发明心电信号去噪方法较佳实施例的流程图，如图1所示，该心电信号去噪方法包括如下步骤：

S100、将含噪心电信号进行EEMD分解；

S200、根据噪声统计特性判定含噪IMF和信号IMF的分界点序号；

S300、双阈值和单阈值两种不同的阈值处理方式的组合对小波阈值去噪对所述含噪IMF分量进行处理，得到小波阈值处理后的IMF分量；

S400、将所述小波阈值处理后的IMF分量与所述信号IMF分量重构，得到去噪的心电信号。

先分步具体说明如下：

步骤1：将含噪心电信号f进行EEMD分解，具体步骤如下：

步骤1-1：将信号f加入均值为零、方差为常数的高斯白噪声s(t)，得到一个新的信号X(t)，即X(t)＝f(t)+s(t)；

步骤1-2：将信号X(t)进行EMD(经验模态分解)分解为有限个IMF分量和趋势项，即其中j表示第j个IMF分量；

步骤1-3：每次加入不同的白噪声，重复上述两个步骤，得到其中i表示第i次加入白噪声；

步骤1-4：将上述结果进行总体平均运算，用来消除多次加入白噪声对真实IMF的影响，得到IMF分量，即其中s表示IMF分量的个数。

步骤2：根据噪声统计特性判定含噪IMF和信号IMF的分界点序号k，即，计算各阶IMF与原始心电信号x(t_i)的相关系数R(x(t_i)，IMF_j(t_i))；它由下式确定：

其中，cov()代表协方差。含噪IMF和信号IMF的分界点序号k由下式确定：

其中，firstlocal表示第一个局部极小值点。

根据分界点k，将其前面视为含噪IMF分量，其后视为信号IMF分量。

步骤3：对含噪IMF分量进行小波分解并分别采用本发明的两种阈值方法对不同尺度的小波系数进行处理，图2为本发明含噪IMF分量的去噪方法较佳实施例的流程图，如图2所示，具体步骤包括：①对含噪的IMF₁-IMF_i分量进行小波分解，分解至N层；②利用本发明采用的双阈值方法对前K个层次的小波系数进行阈值处理；③利用本发明采用的单阈值方法对第K+1到N层的小波系数进行阈值处理；④阈值处理后进行小波重构，即可得到一个去噪的IMF_i分量；⑤重复上述四个步骤，将i个含噪的IMF分量都进行去噪处理；⑥将去噪后的i个IMF分量与其余信号的IMF分量进行重构，得到最终干净的心电信号。

步骤3-1：采用双阈值的方法进行去噪处理，图3为双阈值方法去噪的流程图，如图3所示，具体步骤包括：小波分解N个尺度后，对前K个尺度采用双阈值的方法进行去噪处理。具体步骤是：①小波分解后得到第一层的小波系数d1；②求取第一层小波系数d1的平均值；③根据该平均值计算出所一层小波系数中的上下阈值；④通过上下阈值对d1的小波系数进行处理，保留上下阈值之间的小波系数，从而达到去噪的目的。⑤重复上述四个步骤，直到前K个尺度结束。

小波系数d_i平均值g的计算公式：

其中，n为信号长度，d_i为小波系数；

上下阈值的计算公式：

Ht＝g+(Max-g)·θ (4)

Lt＝g-(g-Min)·θ (5)

其中，Ht为上阈值，Lt为下阈值，Max和Min分别为小波系数中的最大值和最小值，θ为上下阈值选取参数，其取值范围为(0,1)；

本文采用了双阈值的方法，同过上下阈值的限定，可以把阈值缩小在一个比较窄的范围之内，从而降低了噪声误检的可能性。另外，这种方法相比较与传统阈值，其计算量小，计算复杂度较低。

步骤3-2：小波分解N尺度后，对第K+1到N个尺度采用改进的单阈值进行处理，图4为单阈值去噪方法的流程图，如图4所示，采用的新的阈值函数为：

其中，sgn()为符号函数，d_i为小波系数，α和β为调整阈值函数结构的参数；λ为改进阈值，它的计算公式为：

其中，N为信号长度，j为分解尺度，σ为高斯白噪声的标准差，其中，d_i是第i尺度的小波系数；

本发明采用了一种新的阈值函数见(6)式，它主要是能够继承软阈值函数的平滑性，因此该函数具有高阶可导性；同时也要能够克服硬阈值的缺点，因此该函数具备在阈值点处连续的特性。

式(6)中，λ为阈值，参数α和β共同决定了阈值函数的处理与结构，其中，α∈[0,1]，β≥0。随着α和β在各自固定的取值区间内取值，该阈值函数也会产生不同的效果。

①当α＝0时，无论β取何值，该阈值函数就变成了硬阈值函数；

②当α∈(0,1]且β->+∞时，该阈值函数也变成了硬阈值函数；

③当α∈(0,1]且β＝0时，该阈值函数就变成了一种类似于软阈值的阈值函数。

另外在该函数中，参数α可以调节阈值函数对小波系数的压缩程度，弥补了传统软阈值函数在这方面的不足。

由此可见，通过改变参数α和β的取值，可以决定该阈值函数对小波系数的控制。

①β值的改变主要确定了阈值函数的趋向，是变成软阈值函数还是硬阈值函数，或是一种位于两个函数之间的类似软阈值函数；

②α值的改变主要决定了阈值函数对小波系数压缩程度。

综上所述，式(6)的阈值函数不仅在小波域内具有连续性，而且在|x|≥0时具有高阶可导的性质。通过α和β这两个参数的共同作用，该阈值函数就变成了硬阈值函数和软阈值函数的一种推广函数，这样该阈值函数不仅继承了软阈值函数的连续性、平滑性的优点，而且还克服了软阈值函数在处理过程中小波系数与原系数之间存在固有偏差的缺点。另外该阈值函数也具备了硬阈值函数的优点，同时也克服了它的缺点。这样使用起来更加方便、灵活，去噪效果更好。

在阈值选取上，本发明采用了两种阈值组合的方式对小波分解后的小波系数进行了处理。

1)由于小波分解后，低层次的噪声含量会比较大，若是选用单阈值的方式处理，可能由于某些含有信息量的小波系数低于给定阈值而被误认为是噪声去除，从而加大了噪声误检的可能性。故本发明采用了双阈值的方法：

①通过对分解层数的信号长度N，求取分解层小波系数的平均值g，见式(3)；

②由平均值以及参数θ的选取，分别求取该层次内的上下阈值Ht和Lt，见式(4)和式(5)，θ的取值范围为(0,1)，当噪声比较密集时，θ的取值应较小，一般在0.1-0.2；当噪声比较稀疏时，θ的取值应较大；从而使得阈值选取更为精确,去噪效果会更好。

由于双阈值的方法可以把阈值缩小在一个比较窄的范围之内，从而降低了噪声误检的可能性。在这里，为了达到最佳的去噪效果，参数θ取值为0.1。

2)对于较高的层数，则可以采用单阈值去噪的方法：

①基于本发明采用了一种新的阈值函数,见式(6)，它是通过α和β这两个参数的共同作用，使得该阈值函数变成了硬阈值函数和软阈值函数的一种推广函数。

②阈值选取方式是在统一阈值的基础上加以改进的。

②-1保留了传统统一阈值中的噪声标准差σ以及信号长度N；

②-2在分母上添加了分解层数j，使得阈值可以随着分解尺度而改变，分解尺度越大阈值就会相应的减小，这符合噪声的统计特性，随着分解层数的增加，含噪的小波系数分布也有所不同。

综上所述，新增加的分解层数j可以针对小波分解中不同层次的小波系数做相应的不同处理，从而增加了阈值的自适应性，减少了传统使用传统统一阈值所引起的偏差。

步骤4：将小波新阈值处理后的IMF分量与其余信号IMF分量重构，得到最终去噪的心电信号f’。将本发明心电信号去噪方法进行实验仿真，结果如图5所示。

本发明采用了双阈值和单阈值两种不同的阈值处理方式的组合对小波阈值去噪进行处理。采用双阈值的方法，降低了噪声误检的可能性。采用了一种新的阈值函数，它是硬阈值函数和软阈值函数的一种推广，在一定程度上既继承了它们的优点也克服了它们各自存在的缺点；同时还对传统的统一阈值加以改进，通过添加分解层数j，使得每层的阈值会根据j的不同而不同，增加了阈值的实用性，从而提高了效率。

上面结合附图对本发明的实施例进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨和权利要求所保护的范围情况下，还可做出很多形式，这些均属于本发明的保护之内。

Claims (6)

1.一种含噪IMF分量的去噪方法，包括如下步骤：

S10、对第一个含噪的IMF分量进行小波分解，分解N层，N为自然数；

S20、前K个层次的小波系数利用双阈值方法处理，K为自然数且K＜N；

S30、第K+1到N层中的小波系数利用单阈值方法处理；

S40、阈值处理后进行小波重构，得到第一个去噪的IMF分量；

S50、重复上述S1～S4四个步骤，将所有含噪的IMF分量都进行阈值去噪处理。

2.根据权利要求1所述的含噪IMF分量的去噪方法，其特征在于，所述步骤S20包括如下步骤：

S21、小波分解后得到第i层的小波系数di，i≤K且i为自然数；

S22、求取第i层小波系数di的平均值g，计算公式如下：

其中，n为信号长度，di为小波系数；

S23、根据该平均值计算出第i层小波系数中的上阈值Ht、下阈值Lt，计算公式如下：

Ht＝g+(Max-g)·θ

Lt＝g-(g-Min)·θ，其中，Max和Min分别为第i层小波系数中的最大值和最小值，θ为上下阈值选取参数，其取值范围为(0,1)；

S24、通过上阈值Ht、下阈值Lt对di的小波系数进行处理，保留上下阈值之间的小波系数，达到去噪的目的；

S25、重复上述S21～S24四个步骤，直到前K个尺度都进行双阈值去噪处理。

3.根据权利要求1所述的含噪IMF分量的去噪方法，其特征在于，所述步骤S30包括如下步骤：

S31、小波分解后得到第i层的小波系数di，K+1≤i≤N且i为自然数；

S32、使用如下阈值函数处理第i层小波系数di：

其中，sgn()为符号函数，d_i为小波系数，α和β为调整阈值函数结构的参数；λ为改进阈值，计算公式为：

其中，N为信号长度，j为分解尺度，σ为高斯白噪声的标准差，

S33、重复上述S21和S21两个步骤，直到第K+1～N个尺度都进行单阈值去噪处理。

4.一种心电信号去噪方法，其特征在于，包括如下步骤：

S100、将含噪心电信号进行EEMD分解；

S200、根据噪声统计特性判定含噪IMF和信号IMF的分界点序号；

S300、使用如权利要求1～3中任意一种方法处理所述含噪IMF分量，得到小波阈值处理后的IMF分量；

S400、将所述小波阈值处理后的IMF分量与所述信号IMF分量重构，得到去噪的心电信号。

5.根据权利要求1所述的心电信号去噪方法，其特征在于，所述步骤S100包括如下步骤：

S101、将信号f加入均值为零、方差为常数的高斯白噪声s(t)，得到一个新的信号X(t)，即X(t)＝f(t)+s(t)；

S102、将信号X(t)进行EMD分解为有限个IMF分量和趋势项，即其中j表示第j个IMF分量，j≤N且j为自然数；

S103、每次加入不同的白噪声，重复上述S101和S102两个步骤，得到其中i表示第i次加入白噪声，i≤N且i为自然数；

S104、将上述结果进行总体平均运算，用来消除多次加入白噪声对真实IMF的影响，得到IMF分量，即其中s表示IMF分量的个数。

6.根据权利要求1所述的心电信号去噪方法，其特征在于，所述步骤S200包括如下步骤：

S201、计算各阶IMF与原始心电信号x(t_i)的相关系数R(x(t_i)，IMF_j(t_i))；它由下式确定：

其中，cov()代表协方差；

S202、含噪IMF和信号IMF的分界点序号k由下式确定：

其中，firstlocal表示第一个局部极小值点；

S203、根据分界点k，将其前面视为含噪IMF分量，其后面视为信号IMF分量。