一种滤除心电信号50Hz工频干扰的算法

孙九菊，郭峰林，杨 茜

(武汉工业学院数学与计算机学院，湖北武汉430023)

心电信号是由一系列的波形组成的，通过电极对心电信号进行提取，可以画出心电信号的电压幅度随着时间变化的图形，一个典型的完整的正常心电波形如图1所示。由图1可以看出：心电波形主要是由P波、QRS复合波、T波和U波组成，各波段都有各自的生理意义。信号的幅度很小，一般为10 μV—4 mV，典型值是1 mV［1］。

图1 正常心电图波形

一般正常心电图(Electro Cardio Gram，ECG)信号在0.05—100 Hz频率范围内，而90%的心电信号频谱能量主要集中在0.25 Hz到35 Hz之间［2］，其中频率最高的QRS综合波群能量主要集中在2 Hz至20 Hz之间，在12 Hz左右最大。

从人体采集到的心电图常常会受到各种干扰的影响，其中，以50 Hz工频干扰对信号的影响最大，因此，如何抑制50 Hz工频干扰是最重要的问题之一。目前比较常用的方法有：自适应滤波［3］，平滑滤波［4］，小波变换［5］，levkov 滤波等等。平滑滤波虽然算法简单，易于快速实时处理，但是伴有严重的削峰现象同时对信号的衰减很大，所以这种方法无法满足诊断要求;自适应滤波效果虽然好，并且能够跟随干扰频率的变化，但是自适应滤波计算量比较大，而且需要一定的学习时间。所以，该文在Levkov滤波的基础上，结合基于小波变换的QRS检测函数，对心电信号中工频干扰进行滤除。此算法既能比较好地实现滤波，又可以保持波形不失真，同时也克服了前几种滤除工频干扰算法的不足。

1 算法原理

1984年levkov首先对ECG信号的线性段和非线性段采用不同的数字滤波法［6］。1988年Christov改进该算法［7］，引入了ECG信号的线性段判断M来加快滤波的速度。在同样的硬件环境下，不同人的ECG噪声水平不完全相同。对于噪声较大的ECG信号，该方法滤波效果不是很好，这与M的选择有关，M较小，滤波效果较差;M较大，滤波效果较好，但是QRS削峰较为严重，所以在此引入了基于小波变换的QRS检测函数，同时通过实验，该方法取得了很好的效果。

1.1 Levkov滤波原理

Levkov滤波原理：设原始信号为Si，干扰信号为 Ni，包含噪声的信号为：Xi=Si+Ni，i∈0，1，2，….设Si为线性段信号，则Si+1-Si=d，d为常数。

对于采样频率为500 Hz的ECG信号，相邻间距为2 ms，故50 Hz干扰经历的一个周期内有10个采样数据，且有于是：

因N10为50 Hz干扰的下一周期的起点，故：

X10-X0=(S10+N10)-(S0+N0)=S10-S0=10d.即：

由(1)得：

在S4处的干扰信号为：N4=X4-S4.

对于线性段可以按照上述模式处理。当出现非线性段时，考虑到工频干扰的相位在短时间之内不会有很大的变化，可以用最近计算出来的并且与当前同相位处的噪声作为当前噪声的预测值，则当前的信号值为：S4=X4-N4;实验中发现对于噪声较大的ECG波形，该方法滤波效果和M取值有关。当M大于一定的值时，工频干扰滤除的较干净，但是出现QRS波严重削峰的现象，所以在此基础上该文引入了基于小波变换的QRS检测函数。

在心电判别中，该文使用了二进制离散小波变换［8-9］：

通过分析ECG信号功率谱的特点以及小波变换的尺度和信号频率间的关系，这里对ECG信号进行5级小波分解。可以发现：QRS波群的能量大部分集中在23尺度上，在大于24尺度上就大大减小，因此在对心电信号进行处理时，通常取21到24尺度上的小波变换结果来进行分析。

根据小波变换理论，信号的奇异点对应于其小波变换的一个正模极大值和一个负模极大值，其位置对应于正模极大值和负模极大值的过零点。具体的算法是基于小波变换的模最大值max|Wf(2j，τ)|，当计算的信号模大于一定的阈值时，就判定为QRS波，通过检测过零点，就可以判定具体的R波位置，同时用R峰的值去更新阈值。当检测出一个QRS波后，求出其最大值，然后再调整阈值，并置一标志，把其作为levkov法中单独处理QRS波的依据。阈值变化曲线可以用以下的分段函数来表示：

其中：a1≥f(θ+D1+1)＞...＞ f(θ+D2)≥a2.

θ代表最近检测到的QRS波的位置，D1和D2分别代表阈值变化的两个转折点。

1.2 算法步骤

设原始信号为Si，干扰信号为Ni，包含噪声的信号为：Xi=Si+Ni，i∈ 0，1，2，….

(1)通过QRS检测函数判断是不是QRS波。基于小波变换的模最大值是max|Wf(2j，τ)|，首先在24尺度上去寻找模最大值的位置，然后在此位置附近去寻找23尺度上的模最大值，如果不是该尺度上模最大值，该段就不是QRS波;这样循环下去，一直找到21尺度上，只有在21到24尺度上都是模最大值的才是QRS波。

(2)如果是QRS波，则对该点信号实行levkov法单独处理，单独处理过程为：

在S4处的干扰信号为：

如果不是QRS波：判断|D1-D2|是否小于M?如果|D1-D2|＜M，则对该点波形进行线性处理。如果|D1-D2|不小于M，则对该点波形进行非线性处理。

(3)依次重复(1)和(2)的过程，直至波形结束。

2 实验结果及算法分析

2.1 实验结果

图2和图3以及图4分别是滤波前后的比较。

图2 原始心电图

图3 levkov滤波产生的削峰现象

图4 新的levkov滤波

由图2—图4可见，新的levkov滤波对于抑制噪声效果很好，并且对信号几乎没有削峰。这种算法在滤除50 Hz工频干扰时兼有各种其他(如平滑滤波、带阻滤波等)算法的优点。

2.2 算法分析

2.2.1 对于一些不规则波形，较难判断波峰的位置，而该算法可以很好的判断出QRS波的位置，从而能在有效抑制工频干扰的情况下，提高整体的运行速度。

2.2.2 该算法简单并且易于实现，同时该算法对波形的改变基本上不会超过2 μv，即它能够较好的保持原始波形的形状。

3 结论

通过实验，该算法在Levkov滤波的基础上，结合基于小波变换的QRS检测函数，不仅对心电信号中工频干扰进行了有效的滤除，同时也确保信号的不失真。

［1］ Willis J，Tom Pkins.Biomedical Digital Signal Processing-C language Examples and Laboratory Experiments for the IBM PC［M］.London：Prentice Hall，1993，11-20.

［2］ 张正国，程小明，林金森.正常人体表心电图形卷的逐拍变化［J］.生物物理学报，2000，16(1)：105-113.

［3］ Wu G B，Zhu L Y.A LMS adaptive filtering algorithm with variable step size［J］.Acta Electronica Sinica，1994，22(1)：55-60.

［4］ Hu G S.Digital signal processing-theory，algorithm and implement［M］.Beijing：Tsinghua University，1997：330-333.

［5］ Donoho D L.Denoiaing by soft-threah holding［J］.IEEE Trans on Information Theory，1995(3)：613-627.

［6］ Xue Q，Hu Y H，Tompkins W.Neural-network-based adaptive matched filtering for detection［J］.IEEE Trans on BME，1992，39：317-318.

［7］ Barbaro V，Bartolinia P，Fierli M.New algorithm for the detection of the ECG fiducial point in the averaging technique［J］.IEEE MBEC，1991，21：129-131.

［8］ Li C，Zheng C，Tai C.Detection of ECG Characteristic Points Using Wavelet Transforms［J］.IEEE Trans Biomed Eng，1995，42(1)：21-28.

［9］ Kadambe S，Murray R.Wavelet Transform-Based QRS Complex Detector［J］.IEEE Trans Biomed Eng，1999，46(7)：838-848.