一種濾除高采樣心電工頻干擾的改進算法

王吉鳴 呂穎瑩 董 晗 吳 韜 包 濤 馮靖祎*

1(浙江大學醫(yī)學院附屬第一醫(yī)院醫(yī)學工程部， 杭州 310003)2(浙江省醫(yī)藥經(jīng)濟發(fā)展中心， 杭州 310012)

一種濾除高采樣心電工頻干擾的改進算法

王吉鳴1呂穎瑩1董 晗2吳 韜1包 濤1馮靖祎1*

1(浙江大學醫(yī)學院附屬第一醫(yī)院醫(yī)學工程部， 杭州 310003)2(浙江省醫(yī)藥經(jīng)濟發(fā)展中心， 杭州 310012)

為了解決高采樣頻率下50 Hz工頻及其倍頻對心電波形的影響并保留有用的心電高頻信息，提出一種濾除高采樣心電工頻干擾的改進算法。首先，通過修改零極點得到改進的IIR梳狀陷波器，對帶工頻干擾的ECG信號進行濾波處理；然后，將濾波后的ECG信號進行線性段判別，并對線性段進行平滑濾波處理以消除振鈴現(xiàn)象。通過R波幅值、信噪比(SNR)和均方誤差(RMSE)對ECG信號工頻干擾的消除結(jié)果進行評價，與其他4種濾波方法(平滑濾波、Notch濾波、自適應濾波和Levkov濾波)相比，IIR梳狀陷波器改進方法得到的ECG信號R波幅值與原始未加工頻的ECG信號最為接近，恢復率達到98%以上，且SNR提高約11倍、RMSE降低約30倍。所提出方法不僅能濾除高采樣下的工頻及其倍頻的干擾，且能最大程度保留有用信號。

高采樣頻率；工頻干擾；IIR梳狀陷波器；線性段；平滑濾波

引言

在心電圖中，大于100 Hz的高頻成分對許多心臟疾病(如冠心病、心肌病、心肌炎等)都有重要的診斷價值[1-3]，并已逐漸引起心血管臨床工作者的廣泛關(guān)注。目前，相關(guān)的心電高頻分析技術(shù)如高頻心電圖、心室晚電位等已在臨床上得到應用，且正處于逐步深入與擴大的研究中。為了獲取心電的高頻信息，需要對心電信號進行高頻率采樣，但在高采樣頻率下，會引入更多的噪聲，如肌電、工頻干擾、隨機噪聲等，從而對ECG信號的噪聲濾除提出了更高的要求。特別地，在工頻干擾噪聲的濾除中，由于其倍頻信號的存在，在濾波器設(shè)計中必須做針對性處理。同時，由于高采樣下的高尖型QRS波(高幅值、短寬度)信號的頻譜與工頻干擾的頻譜相混疊，在濾波過程中易產(chǎn)生振鈴現(xiàn)象[4]，從而對心電圖的診斷造成影響。

為解決工頻干擾的影響，前期研究者已提出多種濾波算法，如平滑濾波法[5]、簡單整系數(shù)帶阻濾波(Notch濾波)法[6]、自適應濾波法[7]、自適應相干模板法[8]、小波變換法[9]、Levkov法[10]及其改進算法[11-12]等，這些算法均有較好的工頻濾除效果，但在高采樣下進行心電濾波時都存在一定的局限性：一是平滑濾波法算法簡單，處理速度快，濾波效果較好，但對心電波形的R波有較大削峰；二是整系數(shù)Notch濾波法具有嚴格的線性相位且只需較少的計算量，但延時較大且在工頻干擾波動時濾波效果下降明顯；三是自適應濾波法能夠跟隨干擾頻率和幅度的變化，取得較高的信噪比，但需要一個與噪聲有關(guān)而與信號無關(guān)的參考信號；四是自適應模板法抑制工頻干擾的效果較好，但在低頻處與窄帶干擾信號處的抑制帶寬必須相等，導致截止頻率高于0.05 Hz的限制；五是小波變換計算量比較大，不適用于實時處理，且小波基和閾值的選擇都需要通過大量的實驗進行研究論證；六是Levkov法計算量相對小，易于實時處理，但由于其涉及線性段和非線性段的選取，要依據(jù)信噪比來決定，較難適應復雜情況下的心電信號處理[13]。

為解決高采樣情況下的心電工頻干擾濾除方法中存在的問題，本研究基于MATLAB濾波器設(shè)計工具箱(FDATool)提出了一種IIR梳狀陷波器改進方法，在最大程度保留有效心電高頻成分的前提下，能夠有效濾除50 Hz工頻及其倍頻干擾，并通過對ECG波形進行線性段判別和平滑處理，有效地消除了振鈴現(xiàn)象的影響。

1 方法

首先對FDATool工具箱設(shè)計的IIR梳狀陷波器進行了改進，讓其僅對50 Hz及其倍頻窄帶位置陷波，而對0 Hz窄帶位置處不進行陷波。然后，針對陷波器在處理具有高尖型QRS波的ECG信號時容易引入振鈴現(xiàn)象的問題，提出了一種ECG信號線性段判別方法，對經(jīng)過IIR梳狀陷波器處理的ECG信號的線性段和非線性段分別進行相應處理，進而完成工頻干擾的完整濾除?？傮w流程如圖1所示。

圖1 去工頻流程Fig.1 The process of elimination of power-line interference

1.1 IIR梳狀陷波器原理及其改進

由于數(shù)字濾波器的頻率響應就是其單位沖激響應在單位圓上的Z變換，在單位圓上相應于所需濾波器阻帶位置的頻率處設(shè)置零點可以使濾波器的頻率響應在所需阻帶頻率處為零。同時為了得到非常陡峭的過渡帶和常數(shù)幅度的通帶特性，必須在Z平面上為每一個零點再配置一個相應的極點。

基于FDATool工具設(shè)計IIR梳狀陷波器的傳遞函數(shù)，因為每一個極點都有一個零點與之相對應，所以能夠用于抑制某一頻率成分及其諧波，有

(1)

式中，N為濾波器階數(shù)，由采樣頻率和工頻決定，參數(shù)a、b的計算方法為

(2)

(3)

由于IIR梳狀陷波器對0 Hz窄帶成分也有抑制作用，且抑制帶寬與干擾信號抑制帶寬相等，而在ECG信號中這部分信號成分需要單獨處理，不能被濾除，因此需要對該IIR陷波器進行改進。

對式(1)進行分解，可以得到式(6)，z=1處的零點對應的是0 Hz阻帶，z=a1/N為相應極點。去掉0 Hz對應的零極點，可以消除陷波器對0 Hz處的影響，最后改進的IIR梳狀陷波器的傳遞函數(shù)如下：

(6)

(7)

1.2 ECG線性段判別和處理

實際的信號處理過程中振鈴現(xiàn)象的產(chǎn)生具有兩個條件：從濾波器的角度來講是因為其具有單位脈沖響應；從信號的角度來說是因為輸入為高尖信號[14]，其表現(xiàn)為持續(xù)一段時間并逐漸衰減的紋波。由于心電波形的QRS波又高又尖，類似于脈沖響應，用陷波器對其處理必然引起振鈴效應，從而影響濾波結(jié)果。

為消除振鈴現(xiàn)象，本研究提出一種ECG線性段判別方法，并對得到的信號線性段進行平滑濾波處理。具體過程如圖2所示，其中閾值T由前一個R波幅值決定，取其1/4。當檢測到一個線性段后，判斷前面數(shù)據(jù)段是否為線性段，若是則對連續(xù)線性段進行N點平滑濾波，否則僅對該線性段進行N點平滑濾波。當檢測到一個非線性段后，則不作任何處理，保持原來的濾波結(jié)果。

圖2 ECG波形線性段判別標準Fig.2 The linear criterion of ECG waveform

1.3 工頻干擾的去噪效果評價

為更加準確地評價工頻干擾的去噪效果，本研究采用R波幅值、信噪比(signal to noise ratio，SNR)及標準均方誤差(root mean square error，RMSE)等指標進行定量分析。QRS波幅值作為心電圖診斷的一項指標，當R波被削峰太多，有些疾病(如心肌梗死)將無法診斷，濾波過程中保留原始R波幅值程度越大，保證相關(guān)分析的準確性就越高。

SNR和RMSE是兩個最基本的信號處理統(tǒng)計評價指標，其計算公式[15]為

(8)

(9)

式中，xi表示原始含工頻干擾的ECG信號，yi表示去工頻干擾后的ECG信號，N表示信號長度。

一般地，SNR越高，RMSE越小，表明去噪后的信號就越接近原始信號，去噪的效果也就越好。

2 結(jié)果

測試數(shù)據(jù)均來源于FlukeMPS450標準心電采集設(shè)備，其參數(shù)設(shè)置為采樣頻率1 000Hz、幅值4mV、心率80。測試數(shù)據(jù)分兩段，第1段是未受工頻干擾的原始ECG信號，第2段是在原始信號上疊加50Hz工頻干擾的ECG信號。

2.1 IIR梳狀陷波器的改進結(jié)果

對于采樣頻率1 000 Hz、工頻50 Hz的信號濾波，IIR梳狀陷波器的幅頻響應和零極點分布如圖3所示；改進的IIR梳狀陷波器，其幅頻響應和零極點分布如圖4所示；可以看出，0 Hz窄帶阻滯消失。

圖3 IIR梳狀陷波器。(a)幅頻響應；(b)零極點分布Fig.3 IIR notch filter. (a)Amplitude-frequency response； (b)Pole-zero diagram

圖4 改進IIR梳狀陷波器。(a)幅頻響應；(b)零極點分布Fig.4 Improved IIR notch filter. (a)Amplitude-frequency response; (b)Pole-zero diagram

圖5 心電波形濾波前后比較。(a) 原始ECG；(b) 工頻ECG；(c)改進的IIR梳狀陷波濾波和ECG；(d)改進的IIR梳狀陷波器處理后再經(jīng)過線性段平滑后的ECGFig.5 ECG waveform comparison before and after filtering.(a)original ECG；(b) ECG with power-line interference；(c) ECG afterimproved IIR notch filter；(d) ECG after improved IIR notch filter and smoothing filter

2.2 濾波結(jié)果

圖5顯示了本濾波方法的濾波結(jié)果，其中(a)為原始ECG波形，(b)為疊加50 Hz工頻干擾的ECG波形，(c)為經(jīng)過改進的IIR梳狀陷波器后的ECG波形，(d)為經(jīng)過改進的IIR梳狀陷波器處理后再經(jīng)過線性段平滑后的ECG波形。從(c)可看到，經(jīng)過改進的IIR梳狀陷波器處理后的50 Hz工頻干擾得到了很好的抑制，但是在QRS波后的心電波形出現(xiàn)了振鈴現(xiàn)象。從(d)可看出，經(jīng)過線性段判別和處理后，振鈴現(xiàn)象基本得到抑制。

2.3 與其他濾波算法結(jié)果比較

圖6 不同方法處理后的ECG波形。(a)原始ECG；(b)工頻ECG；(c)平滑濾波后ECG；(d)Notch濾波后ECG；(e)自適應濾波后ECG；(f)Levkov濾波后ECG；(g)本方法濾波后ECGFig.6 ECG waveform after different methods.(a)Original ECG；(b) ECG with power-line Interference；(c) ECG after smoothing filter；(d) ECG after notch filter；(e) ECG after adaptive filter；(f) ECG after Levkov filter；(g)ECG after the method in this article

圖6、7分別為本方法與平滑濾波、Notch濾波、自適應濾波及Levkov濾波對一段和一個周期的ECG信號進行處理的結(jié)果比較。從中可看出，平滑濾波和Levkov濾波極大地削弱了R波幅值信息，而Notch濾波和自適應濾波雖然保留了大部分的R波幅值信息，但仍存在一部分噪聲干擾。相比其他方法，本濾波方法能夠在獲得較好的濾波效果的同時，最大程度地保留R波幅值信息。

圖7 不同方法處理后的ECG波形(一個周期)。(a)原始ECG；(b)工頻ECG；(c)平滑濾波后ECG；(d)Notch濾波后ECG；(e)自適應濾波后ECG；(f)Levkov濾波后ECG；(g)本方法濾波后ECGFig.7 ECG waveform after filtering by different methods(one period). (a)Original ECG；(b) ECG with power-line interference；(c) ECG after smoothing filter；(d) ECG after notch filter；(e) ECG after adaptive filter；(f) ECG after Levkov filter；(g)ECG after the method in this article

2.4 去噪結(jié)果評價

表1給出了不同方法處理后的R波幅值、SNR和RMSE等參數(shù)對比?？梢钥闯?，本方法、Notch濾波法及自適應濾波方法均能最大化保持R波幅值信息，R波幅值恢復率達到98%以上，但是本方法在SNR和RMSE上表現(xiàn)更佳，SNR提高了約11倍而RMSE降低了約30倍，說明有更好的信號去噪能力。

表1 不同方法處理后的ECG信號的R波幅值、SNR和RMSE

Tab.1 Amplitude of R wave SNR and MSE of ECG signals after different treatments

不同處理方法R波幅值/mVSNR/dBRMSE加工頻后4262150651平滑濾波后2428780164Notch濾波后39420170051自適應濾波后39020230051Levkov濾波后23272240191本方法濾波后395232140024

3 討論和結(jié)論

在心電分析中，為了提取高頻成分，需要更高的采樣頻率。但采樣頻率提高后，需要考慮工頻干擾的倍頻成分對心電波形的影響，使用傳統(tǒng)的濾波方法，相對低頻的工頻陷波的計算量會大大增加。本研究根據(jù)濾除50 Hz工頻及其倍頻并保留心電高頻信息的要求，采用了IIR梳狀陷波器，并通過改變零極點來改變陷波器結(jié)構(gòu)以消除對0 Hz窄帶的影響，從而達到濾波效果，且實現(xiàn)簡單。與傳統(tǒng)的很多工頻干擾濾波方法相比，本方法僅僅濾除了心電信號中的工頻及其倍頻窄帶而保留了其他所有頻率信息，以供后續(xù)分析。局限是IIR梳狀陷波器的階數(shù)是由采樣頻率和工頻決定的，所以要求采樣頻率必須為工頻的整倍數(shù)。

本研究對濾波后的振鈴現(xiàn)象處理方法參考了Levkov法對線性段進行平均濾波的思想，不同的是先對采集的心電信號進行了窄帶寬50 Hz工頻陷波，以便提高信號的信噪比并使R波削峰程度最小，再根據(jù)前面提到的ECG波形線性段判別標準來確定線性段并進行平滑濾波。這樣的處理方式在最大化保持R波幅值信息的同時比Levkov法直接根據(jù)信噪比來判定線性段和非線性段的方法穩(wěn)定很多，進而去振鈴的方法也比較穩(wěn)定。

另外，在心電波形的實際采集過程中工頻頻率并不是固定不變的，工頻波動會影響濾波結(jié)果。本研究提出的線性段平滑濾波方法不僅能消除高尖ECG信號濾波后的振鈴現(xiàn)象，相同原理對由于工頻干擾波動引起的濾波后出現(xiàn)的振鈴現(xiàn)象也有很好的抑制作用。如圖8所示，A、B、C所指分別為工頻波動處、陷波后振蕩及線性段平滑濾波對振蕩的抑制，可以看出，C相對于B，振蕩小了很多。

圖8 工頻干擾波動(A—工頻波動；B—陷波后振蕩；C—線性段平滑濾波對振蕩的抑制)Fig.8 Power-line interference fluctuation (A: Power frequency fluctuation, B: Oscillation after notch filter, C: Inhibition of oscillation by smoothing the linear segment)

心電高頻成分能夠提供對疾病有診斷價值的信息，但在對心電高采樣的同時會引入更多的噪聲，從而對ECG信號的噪聲濾除方法有更高的要求。作為其中的一個重要干擾因素，工頻干擾在高采樣頻率下的濾除仍存在問題。本研究提出一種改進的工頻干擾濾除算法，通過結(jié)合改進的IIR梳狀陷波器以及線性段判別處理，可在最大化保持R波幅值信息的前提下，有效地濾除ECG信號中的工頻及其倍頻成分，并可消除由于陷波器造成的振鈴現(xiàn)象，從而為后續(xù)心電的準確分析提供了保證。

[1] 張開滋，肖傳實，郭繼紅，等.中國心電信息學圖解集成：下冊[M]. 長沙：湖南科技出版社，2010：1656.

[2] 盛海民．心電高頻信息的臨床應用[J]. 醫(yī)學綜述，1996，2(6)：269-270.

[3] 劉恒亮，單新中．心電高頻信息研究的現(xiàn)狀與展望[J]. 心血管病學進展，2000，21(1)：52-55.

[4] 蘭瑞芬，胡廣書．高采樣率下簡單整系數(shù)工頻陷波器的設(shè)計[J]. 航天醫(yī)學與醫(yī)學工程，2008(2)：152-156.

[5] 王鎮(zhèn)，蔡萍．用于去除心電信號中工頻干擾數(shù)字濾波技術(shù)[J]. 電子測量技術(shù)，2000(2)：43-47.

[6] 周靜．心電信號中工頻干擾的消除[J]. 生物醫(yī)學工程研究，2003，22(4)：61-64.

[7] 林紹杰，賴麗娟，吳效明．基于阻抗檢測的自適應消除心電運動偽跡方法[J]. 生物醫(yī)學工程學雜志，2010，27(3)：529-532.

[8] 于超，李剛，林凌．基于簡化FFT可跟蹤工頻干擾的自適應相干模板法的實現(xiàn)[J]. 生物醫(yī)學工程學雜志，2007，24(4)：780-784.

[9] 張雪，王海燕，李保軍．臨床心電信號工頻干擾小波去噪方法對比分析[J]. 計算機測量與控制，2010，18(4)：902-905.

[10] Levkov C, Michov G, Ivanov R, et al. Subtraction of 50Hz interference from the electrocardiogram [J]. Med Bio Eng Compute, 1984, 22(4):371-373.

[11] 孫京霞，白延強，楊玉星. 一種抑制心電信號50 Hz 工頻干擾的改進Levkov方法[J].航天醫(yī)學與醫(yī)學工程.2000，13(3)：196-199.

[12] Levkov C, Michov G, Ivanov R, et al. Removal of power-line interference from the ECG: a review of the subtraction produre[J]. BioMedical Engineering OnLine 2005, 4(1): 18-36.

[13] 劉勝洋，張根選，曹陽．心電信號中的工頻干擾濾除方法[J].生物醫(yī)學工程學雜志，2014(3)：577-582.

[14] 羅旻暉．心電圖機自動分析算法的研究與實現(xiàn)[D]. 哈爾濱：哈爾濱工業(yè)大學，2012.

[15] 張喜紅，王玉香．心電信號工頻干擾陷波器的設(shè)計與實現(xiàn)[J]. 湖南工業(yè)大學學報，2014(2)：72-76.

An Improved Algorithm for Removing the Power-Line Interference from ECG Signals in High Sampling Rate

Wang Jiming1Lv Yingying1Dong Han2Wu Tao1Bao Tao1Feng Jingyi1*

1(ClinicalEngineeringDepartment，TheFirstAffiliatedHospital，SchoolofMedicine，ZhejiangUniversity，Hangzhou310003，China)2(ZhejiangPharmaceuticalEconomicDevelopmentCenter，Hangzhou310012，China)

high sampling rate; power-line interference; IIR comb notch filter; linear segment; smoothing filter

10.3969/j.issn.0258-8021. 2016. 06.015

2015-12-20， 錄用日期:2016-08-21

R318

D

0258-8021(2016) 06-0744-05

# 中國生物醫(yī)學工程學會高級會員(Senior member, Chinese Society of Biomedical Engineering)

*通信作者(Corresponding author)， E-mail: casper_feng@163.com