CN106889984B - 一种心电信号自动降噪方法 - Google Patents

Abstract

本发明涉及一种心电信号的自动降噪方法，其实现过程是：a）检测心电信号的R波波峰位置，并进行初步滤波；b）从初步滤波后的心电信号中选择若干连续心拍的心拍段做均值，得平均模板；c）用平均模板替换所获取的心电信号相对应的心拍段，得到指导信号；d）在指导信号上依次设定若干窗口，然后在每个窗口中，利用指导滤波模型将指导信号经过线性变换得到滤波输出；e）对所有窗口计算完线性系数后，代入指导滤波模型中求得q_i的值作为最终结果。本发明利用平均模板的方法可以在一定的程度上保留住心电信号大部分的幅值特性，然后通过替代心电信号主要部分来构造出指导信号，这样得到的指导信号保留住了心电信号大部分的幅值特征。

Description

一种心电信号自动降噪方法

技术领域

本发明涉及心电信号的自动检测与分析技术，具体的说是一种心电信号自动降噪方法。

背景技术

近年来，心血管疾病已成为威胁人类生命和健康的头号杀手。心血管疾病具有突发性，发病就在一瞬间。往往因为症状短暂，患者赶到医院进行心电采集时，症状已经消失，这时候的心电图往往难以捕捉到有效的诊断依据，不能反映患者的实际心脏健康状态，因为缺少疾病发病期间的心电数据导致患者无法进行相应的治疗。因此，针对心血管疾病最好的办法就是心血管患者住院进行24小时的动态心电监测。在现代化的远程通信技术突飞猛进的社会背景下，远程心电监护仪无疑是最好的解决办法。远程心电监护就需要患者实时穿戴检测设备，在互联网的支持下，数据可以实时传输到定点医院。这其中最大的问题是，穿戴的设备与医院的心电图机不同，数据是在患者进行各种日程生活中采集的，与在医院静态采集的心电数据相比存在很大的干扰信号，例如最常见的电极接触噪声(EM)、肌肉颤动引起的肌电干扰(MA)、人体呼吸移动等引起的基线漂移(BW)。这些在很大程度上影响了心电信号的正常形态，而已有算法在降噪过程中容易丢失心电信号的某些重要特征，对医生的诊断造成很大的偏差，误诊率提高。因此恢复心电信号的原有形态对疾病诊断有十分重要的意义。

发明内容

本发明的目的是提供一种心电信号自动降噪的方法，以解决现有算法在降噪过程中丢失心电信号重要特征的问题。

本发明的目的是这样实现的：

一种心电信号的自动降噪方法，其包括以下步骤：

a)获取一段人体的心电信号，先检测所获取的心电信号的R波波峰位置，然后进行初步滤波以滤除低频信号，去除基线漂移；

b)从初步滤波后的心电信号中选择若干个连续的心拍的心拍段做均值，即得平均模板；从每个心拍中所选择的心拍段的长度包含270个采样点，所述270个采样点由R波波峰之前的89个采样点、R波波峰处的1个采样点及R波波峰之后的180个采样点构成；

c)用步骤b)所得到的平均模板替换步骤a)中所获取的人体心电信号中每一个心拍相对应的心拍段，得到指导信号；

进行替换时，如果一个RR间期的长度小于所述平均模板的长度，会有部分采样点重叠，此时将重叠的采样点进行均值处理；

d)在所述指导信号上依次设定若干窗口ω，然后在每一个窗口中，利用指导滤波模型将所述指导信号经过线性变换得到滤波输出：

其中1：对于在所述指导信号上依次设定的各窗口，(a,b)是在一个窗口中保持不变的一系列线性系数，其值为一个窗口中各采样点的滤波输出与对应的待滤波信号差值最小时所对应的线性系数值；

其中2：I_i表示所述指导信号在所在窗口的第i个采样点的值；q_i表示所述指导信号在所在窗口的第i个采样点的滤波输出值；

所述待滤波信号即为步骤a)中所获取的人体心电信号；

e)对所有窗口ω计算完线性系数(a,b)后，将各窗口内a、I_i和b分别代入指导滤波模型中求得q_i的值作为最终结果，如果该过程中出现不同的窗口涉及同一个采样点的情况时，该采样点的滤波输出值根据其在不同窗口中的计算结果的平均值确定即可。

所述的心电信号的自动降噪方法，步骤d)中，采用最小二乘法计算(a,b)的值。

本方法最大的优点就是可以在滤波的过程中保持住心电信号的微弱特征，心电信号幅值特性十分微弱，这就决定了在去噪的过程中很容易丢失微弱的特征，指导滤波就很好的解决了这个问题，我们可以在指导信号的指导下找回在滤波过程中丢失的特征，我们利用平均模板的方法可以在一定的程度上保留住心电信号大部分的幅值特性，然后通过替代心电信号主要部分来构造出指导信号，这样得到的指导信号保留住了心电信号大部分的幅值特征，同时在进行指导滤波的过程中可以根据原始的待滤波信号找回丢失的那些能量很小的微弱特征，例如P波和T波等。

附图说明

图1是实施例1的步骤(1)所获取的人体心电信号。

图2是实施例1的步骤(5)所得到的指导信号。

图3是实施例1的最终滤波输出结果。

图4是对比例1的滤波输出结果。

图5是对心电信号进行巴特沃斯高通滤波器处理后选取到的极大值点。

图6是优化极大值点流程图。

具体实施方式

下面通过实施例和对比例对本发明的技术方案做进一步解释和说明。

实施例1

(1)获取人体心电信号：采集设备为北京蓬阳丰业的MedSun18导联Holter长时间采集人体的心电信号，其采样输出频率为360Hz，采集心电数据以TXT的形式存储。所采集的心电信号可以很容易的读取到Matlab环境中进行显示，本实施例中截取的人体心电信号包含大概10个心拍，共3000个采样点，如图1所示。

(2)对所采集的心电原始信号数据进行滤波处理：利用巴特沃斯高通滤波器滤除基线漂移噪声，阻带的截止频率为1Hz，阻带衰减最小为30dB，通带的衰减最大为15dB。

(3)对经过巴特沃斯高通滤波器的心电信号进行能量窗变换，并选取极大值点：

(3-1)能量窗变换：按下式，将经过巴特沃斯高通滤波器的心电信号p由时间域分析变换到能量域分析，得到心电信号能量曲线：

其中，E_n表示第n个采样点的能量值；N为所选的窗口长度，取值26；M为总的采样点数；p_n表示所述利用巴特沃斯高通滤波器去噪后的心电信号p的第n个采样点；

(3-2)选取极大值点：将所得到的心电信号能量曲线进行硬阈值化处理，即：

其中，T_h为所选取的阈值，取T_h＝0.3*median(E_n)，

然后选择经过硬阈值化处理后的心电信号能量曲线的波峰位置作为极大值点，如图5；

(3-3)优化极大值点：按现有算法(可参见CN103156599A)，如图6所给出的流程图，设定2个时间阈值t₁和t₂，且t₁＜t₂当任意两个极大值点的时间间隔小于t₁时，就去掉这两个极大值点之间幅值较小的那个；当任意两个极大值点的时间间隔大于t₂时，就在这两个极大值点之间寻找另一未被识别的极值点；如两个极大值点的时间间隔既大于t₁，又小于t₂，则该两个极大值点均保留，如此最终得到的经优化的每个极大值点都对应一个QRS波群。

图6中，E_t表示步骤(4-2)所得到的所有极大值点的时间间隔的平均值，t₁＝0.5×E_t，t₂＝1.5×E_t。

(3-4)根据步骤(3-3)中每个极大值点所在的时间点，在步骤(1)所获取的人体心电信号和步骤(2)中滤波后的心电信号上相应的时间点左右各7个采样点的范围内搜寻信号幅值最大的点，作为检测到的R波波峰。

(4)在得到R波波峰位置后构建平均模板，该平均模板由同一个人的初步滤波后的心电信号中1000个连续心拍的心拍段取均值得到；

所截取的心拍段包含270个采样点，该270个采样点是根据检测到的R波波峰，由R波波峰所在采样点，R波波峰左边(即其之前)的89个采样点和R波波峰右边(即其之后)的180个采样点构成。

(5)将上述得到的平均模板替换步骤(1)中所获取的人体心电信号中每一个心拍中相应的心拍段的位置，得到指导信号，如图2所示；

由于心电信号是一种周期性信号，所以将上述得到的平均模板替换待滤波信号(即步骤(1)中所获取的人体心电信号原始波形)的每个心电信号周期内的相应位置上的心拍段时，会出现两种情况：

①当一个RR间期(即一个心电信号周期)的长度小于上述得到的平均模板的长度：由于平均模板的长度大于一个心电信号的周期，这样就会造成有部分采样点会重叠，此时需将重叠的采样点再进行均值处理，从而得到与原始信号长度相同的重构信号。

②当一个RR间期(即一个心电信号周期)的长度大于或等于上述得到的平均模板的长度：由于平均模板的长度小于或等于一个心电信号的周期，所以可以正常的完成替换。

(6)利用指导信号对待滤波信号中的噪声进行滤除：

指导滤波的模型是由指导信号I到滤波输出q的一个局部线性变换模型，本例中，设定窗口长度为|ω|＝401，所以对于这个指导信号I而言，第201个采样点为第一个窗口ω₁的中心点，第202个采样点为第二个窗口ω₂的中心点，依此类推，第j(j＝201,202，…，3000)个采样点为第(j-200)个窗口ω_k(k＝j-200)的中心点，在设定的窗口ω_k(k＝1，2，3,…,j-200)中，q_i由I_i经过线性变换得到：

公式1：

公式1中，(a_k,b_k)是在窗口ω_k中保持不变的线性系数；I_i表示指导信号在所在窗口的第i个采样点的幅值；q_i表示指导信号在所在窗口的第i个采样点的滤波输出幅值。

为确定各窗口的线性系数(a_k,b_k)的值，需要求出公式1的一组使得q_i与输入信号p_i(即待滤波信号)相差最少的解，为此我们在窗口ω_k中利用如下最小二乘法进行计算：

公式2：

这里p_k表示窗口ω_k中的原始待滤波信号在相应采样点的幅值，ε是正则化系数，防止a_k变得过大，取值ε＝0.2²，通过线性回归的方法得出公式2的解：

其中：μ_k与分别为在窗口ω_k中的指导信号I的各采样点的均值与方差，|ω|为ω_k中的采样点个数，为在窗口ω_k中的输入信号p_k的各采样点的均值。

(7)在对指导信号中所有窗口ω_k计算完(a_k，b_k)以后，对整个心电信号的所有局部窗口应用公式1的线性模型，但是出现的问题是所有包含采样点i的窗口ω_k都要牵涉到采样点i，导致公式1中的q_i在不同窗口中的计算结果不同，此时所采取的策略是把所有计算得到的q_i取平均值作为最终结果，即：

公式3：

此处是在求q_i的平均值，有多个窗口会包含采样点i，不同的窗口里面的i对应不同的系数(a_k，b_k)，因此，我们求多个窗口中(a_k，b_k)的平均值，因为这些窗口ω_k里面都含有采样点i和系数(a_k，b_k)，因此有k|i∈ω_k，

其中

这里k∈ω_i中，k表示所有包含采样点i的窗口的个数，其值正好等于窗口长度401，ω_i表示包含同一个采样点i的窗口；表示所有包含同一个采样点的窗口的系数的平均值，根据计算，包含采样点i的窗口正好与设定窗口的长度相等。

根据以上处理步骤得滤波结果，如图3所示。

对比例1

同样取实施例中所选取得心电信号原始波段，首先选择db6作为基础小波基，对该心电原始信号进行8层小波分解，得到各个尺度上的小波系数 这里i为层标号。

为了保证去除噪声的同时尽可能多的保留有用信号，本研究采用具有尺度自适应性的阈值法：

其中，Tⁱ为设定的阈值，Nⁱ为第i层系数的个数，e是自然常数，

采用软阈值方法处理各个尺度上的小波系数：

对阈值化处理后的小波系数进行重构，得到滤波结果，如图4所示。

对比例1详细处理步骤可参照Reddy,G.U.,Muralidhar,M.,&Varadarajan,S.(2009).Ecg denoising using improved thresholding based on wavelettransforms.International Journal of Computer Science&Network Security,253-257.。

观察原始心电信号可以看出，在心电信号发生T波倒置后，原始采集的心电图会被噪声严重干扰，导致T波倒置发生严重偏差。对比图3与图4可看出，本发明方法的输出结果(图3)很好的保持住了心电信号的细微特征，很明显可以看出，T波倒置被完全的识别，同时去除大部分噪声，恢复了心电信号原有的形态特征，本发明的方法根据指导信号提供信号的大体形态，综合输入噪声信号的细节特征，结合二者的特征，完美的恢复出T波和P波这类能量极其微弱的信号的形态特征，找回了极具医学参考价值的信息，为医生的诊断提供更好的帮助。采用现有普通滤波方法得输出结果(图4)可以看出，利用小波去除心电信号的噪声，导致信号偏离基线的基线漂移噪声基本滤除，使得信号回到基线附近，但是仍有大量的噪声干扰波形，心电信号典型的P波、QRS波群、以及T波都被噪声严重影响导致无法看出原有正常的形态，这也就是说大量的有用心电信息被丢失，或者说是没有被恢复，这样的滤波结果没有达到很理想的去除噪声的目的。

Claims (2)

1.一种心电信号的自动降噪方法，其特征是，包括以下步骤：

a)获取一段人体的心电信号，先检测所获取的心电信号的R波波峰位置，然后进行初步滤波以滤除低频信号，去除基线漂移；

b)从同一个人的初步滤波后的心电信号中选择1000个连续的心拍的心拍段做均值，即得平均模板；所选择的心拍段的长度包含270个采样点，所述270个采样点由R波波峰之前的89个采样点、R波波峰处的1个采样点及R波波峰之后的180个采样点构成；

c)用步骤b)所得到的平均模板替换步骤a)中所获取的人体心电信号中每一个心拍相对应的心拍段，得到指导信号；

进行替换时，如果一个RR间期的长度小于所述平均模板的长度，会有部分采样点重叠，此时将重叠的采样点进行均值处理；

d)在所述指导信号上依次设定若干窗口ω，然后在每一个窗口中，利用指导滤波模型将所述指导信号经过线性变换得到滤波输出：

其中1：对于在所述指导信号上依次设定的各窗口，(a,b)是在一个窗口中保持不变的一系列线性系数，其值为一个窗口中各采样点的滤波输出与对应的待滤波信号差值最小时所对应的线性系数值；

其中2：I_i表示所述指导信号在所在窗口的第i个采样点的值；q_i表示所述指导信号在所在窗口的第i个采样点的滤波输出值；

所述待滤波信号即为步骤a)中所获取的人体心电信号；

e)对所有窗口ω计算完线性系数(a,b)后，将各窗口内a、I_i和b分别代入指导滤波模型中求得q_i的值作为最终结果，如果本步骤计算过程中出现不同的窗口涉及同一个采样点的情况时，该采样点的滤波输出值根据其在不同窗口中的计算结果的平均值确定即可。

2.根据权利要求1所述的心电信号的自动降噪方法，其特征是，步骤d)中，采用最小二乘法计算(a,b)的值。