針對心電信號(ECG)在采集和轉換的過程中容易受到工頻干擾(PLI)的問題, 提出了一種基于形態分量分析(MCA)和集合經驗模態分解(EEMD)的PLI消除新算法。首先根據ECG特征波形的形態差異性, 利用MCA將其分解為突變成分、平滑成分和殘余白噪聲成分, 然后對含PLI的平滑成分進行EEMD, 再濾除PLI的本征模態函數(IMF), 最后重構ECG信號。文中采用噪聲抑制率(NSR)和信號失真率(SDR)來評價算法的降噪效果。通過實驗發現, 該算法不僅能夠有效地濾除工頻干擾, 而且SDR值較小, 濾波效果優于改進的Levkov算法。
引用本文： 趙偉, 肖世校, 張保燦, 黃曉菁, 游榮義. 基于形態分量分析和集合經驗模態分解的心電信號工頻干擾消除法. 生物醫學工程學雜志, 2015, 32(6): 1179-1184. doi: 10.7507/1001-5515.20150209 復制
0 引言
心電信號(electrocardiogram, ECG)在診斷心血管疾病、評估各種治療方法的有效性等臨床應用中有重要意義[1-3]。體表采集到的ECG是一種非線性非平穩的微弱信號，易受到來自外界的各種干擾，主要有工頻干擾(power line interference, PLI)、肌電干擾和基線漂移等，降低了正常診斷率。其主要干擾之一為醫療設備電源PLI，有時其幅度之大甚至可能淹沒需檢測的ECG，且落在ECG有效頻帶內，將直接影響到系統的正常工作以及醫生對疾病的準確診斷。因此，如何有效地抑制PLI噪聲是心電智能診斷的一個重要內容。
目前濾除ECG中PLI的主要方法包括：濾波器濾波法[2-4](如陷波濾波法、自適應濾波法、Levkov濾波法等)和小波變換法[1, 5]。濾波器濾波法通常利用頻帶差異性進行濾波，由于PLI和ECG信號在頻域有交疊，因此濾波的同時會造成ECG特征波形的損失, 從而導致錯誤的診斷結果；小波變換法則存在小波基的選擇和閾值的設定問題。
形態分量分析(morphological component analysis, MCA)根據各組成成分的形態差異性進行分離。參考文獻[6]根據ECG特征波形的形態差異，利用MCA先將ECG分解為突變成分、平滑成分和殘余成分，通過濾除殘余成分有效地濾除白噪聲的干擾，取得了較好的效果，但它無法濾除PLI。由于PLI主要為變化較為平滑的信號，為了避免消除PLI的同時對特征波形(如QRS波群)造成損傷，本文提出了一種基于MCA和集合經驗模態分解(ensemble empirical mode decomposition, EEMD)的消除ECG中PLI的新算法。在文獻[6]的基礎上，先采用MCA提取出突變成分(如QRS波群)，再利用EEMD對平滑成分進行濾波，濾除PLI的本征模態函數(intrinsic mode functions, IMF)，最后重構ECG信號。
1 MCA基本原理
MCA是Starck等[7]近來基于信號的形態多樣性和稀疏表示提出的一種基于多種基函數的稀疏信號分解方法。該方法的主要思想是利用信號組成成分的形態差異性，通過構建不同形態的稀疏表示字典來實現信號中各形態成分的分離，已廣泛應用于圖像處理、醫學信號處理、機械故障診斷等領域[8-10]。
對于任意的觀測信號s∈R，假設s是K個不同形態分量sk的線性組合，即。MCA的目標是從混合信號中恢復出不同形態的各個源信號。這是一個難解的逆問題。MCA假定字典是由一系列基本變換組成，其中源信號sk在變換基Φk上的表示系數αk是稀疏的，而在其他變換基上不能稀疏表示，即字典Φk相互之間不相關。因此字典Φk可以用來區分信號類型，從混合信號中提取出源信號sk。信號s的稀疏分解可歸結為變換系數{α1opt, α2opt, …, αkopt}的優化問題：
| $\left\{ \alpha _{1}^{opt},\alpha _{2}^{opt},\cdots ,\alpha _{\text{K}}^{opt} \right\}=\text{Arg}\underset{\left\{ {{\alpha }_{1}},{{\alpha }_{2}},\cdots ,{{\alpha }_{K}} \right\}k=1}{\mathop{\text{min}}}\,{{\left\| {{\alpha }_{k}} \right\|}_{1}}+\lambda \left\| s-\sum\limits_{k=1}^{K}{{{\Phi }_{k}}{{\alpha }_{k}}} \right\|_{2}^{2}$ |
式中‖αk‖0表示l0范數。放寬式(1)中的約束條件，當信號稀疏分解后的剩余成分類似高斯白噪聲時，可以選擇l2作為誤差范數，則式(1)轉化為：
| $\left\{ {\alpha _1^{opt}, \alpha _2^{opt}, \cdots, \alpha _k^{opt}} \right\}={\text{Arg}}\mathop {{\text{min}}}\limits_{\left\{ {{\alpha _1}, {\alpha _2}, \cdots, {\alpha _K}} \right\}} \sum\limits_{k=1}^K {{{\left\| {{\alpha _k}} \right\|}_1}} + \lambda \left\| {s-\sum\limits_{k=1}^K {{\Phi _k}{\alpha _k}} } \right\|_2^2$ |
式中λ為給定的閾值。可見，MCA方法不僅可以將信號分解成不同的形態分量，還能濾除高斯白噪聲。Starck等[7]在BCR(Block-Coordinate Relaxation)[11]方法的基礎上給出了MCA的數值實現。
2 EEMD基本原理
EEMD是Wu等為了解決模式混疊問題，提出的一種噪聲輔助數據分析方法[12-13]。它完全依賴待檢測信號本身的特征進行模態分解，非常適用于非線性、非平穩信號的處理。
EEMD方法的核心是經驗模態分解(empirical mode decomposition, EMD)。EMD是由Huang等提出的一種自適應的、高效的信號分解方法。它根據信號自身特征將其從高頻到低頻分解成有限個具有物理意義的IMF，這些模態函數能夠內在地描述信號的特征。EMD的基本理論與分解步驟詳見文獻[12]。當待分析信號中某個頻段的分量不連續時，就會導致一個IMF中出現多個不同尺度的信號成分，或相似尺度的成分出現在不同的IMF中，即出現模式混疊現象。為了避免模式混疊，EEMD利用了白噪聲信號具有在各個頻段能量一致的特點，對待分析信號添加白噪聲后再進行EMD分解以保證每個IMF時域的連續性。在待分析信號中添加多次白噪聲后，利用白噪聲均值為零的特性，對各次EMD分解得到的同階IMF分量求取均值，以消除附加白噪聲的影響。
EEMD具體分解步驟如下：
(1)在原始信號s(t)中加入多次幅值系數為k的白噪聲，即：
| ${s_i}\left(t \right)=s\left(t \right) + k \times {n_i}\left(t \right)\; \;\; \left({i=1, 2, 3, \cdots, N} \right)$ |
式中，N為加入白噪聲的總次數，si(t)為第i次加入白噪聲后的信號，ni(t)為第i次加入的白噪聲，N和k通常根據經驗來選取。
(2)依次對si(t)進行EMD分解，分別得到各階的IMFij。IMFij是si(t)分解的第j個IMF，即：
| ${s_i}\left(t \right)=\sum\limits_{j=1}^J {{\text{IM}}{{\text{F}}_{ij}}} $ |
(3)利用不相關隨機序列的統計值為零的原理，將各階相應的IMF進行總體求平均運算便可獲得EEMD的各個IMF分量，即EEMD分解獲得的第j個IMF為IMFj(t)，公式如下：
| ${\text{IM}}{{\text{F}}_j}\left(t \right)=\frac{1}{N}\sum\limits_{i=1}^J {{\text{IM}}{{\text{F}}_{ij}}\left(t \right)} $ |
3 基于MCA和EEMD的工頻干擾消除法
ECG信號的主要特征波形包含突變成分(如QRS波群)和平滑成分，PLI為變化較為平滑的50 Hz正弦波信號，因此消除PLI只需對平滑成分進行濾波即可，從而有效地避免對QRS波群的損傷。對含PLI的平滑成分進行EEMD，平滑成分被自適應地從高頻到低頻依次進行分解。多次實驗結果表明，PLI絕大部分的能量集中在EEMD的第一個本征模態函數IMF1上, 而且IMF1分量上幾乎全為PLI，因此濾除IMF1和MCA分解的白噪聲成分后重構ECG，即可達到消除PLI的目的，同時濾除了白噪聲。
假設無污染的ECG信號為s(t)，PLI的強度為β，受PLI的ECG信號為s′(t)，消噪后的ECG信號為。根據ECG的特征波形，采用非抽樣小波字典來稀疏表示s′(t)的突變成分mm(t)，采用離散余弦變換字典來稀疏表示s′(t)的平滑成分ms(t)，MCA分解的殘余噪聲成分用r(t)表示，則添加噪聲、MCA、EEMD及濾除PLI模型分別為式(6)、式(7)、式(8)、式(9)。
| $s'\left(t \right)=s\left(t \right) + \beta \times \sin \left({2\pi \times 50t} \right)$ |
| $s'\left(t \right)={m_{\text{m}}}\left(t \right) + {m_{\text{s}}}\left(t \right) + r\left(t \right)$ |
| ${m_s}\left(t \right)=\sum\limits_{i=1}^J {{\text{IM}}{{\text{F}}_i}\left(t \right)} $ |
| $\hat s\left(t \right)={m_{\text{m}}}\left(t \right) + \sum\limits_{i=2}^J {{\text{IM}}{{\text{F}}_i}\left(t \right)} $ |
3.1 算法評價指標
為驗證算法的有效性，采用信噪比(signal noise ratio, SNR)、噪聲抑制率(noise suppression ratio, NSR)和信號失真率(signal distortion ratio, SDR)作為評估指標，其定義如下：
| ${\text{SNR=10}} \times {\log _{10}}\left({\sum\limits_{t=1}^N {{s^2}\left(t \right)} /\sum\limits_{t=1}^N {{n^2}\left(t \right)} } \right)$ |
| ${\text{NSR=1-}}{\left[{\sum\limits_{t=1}^N {{{\left| {\hat s\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}{\text{/}}{\left[{\sum\limits_{t=1}^N {{{\left| {s'\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}$ |
| ${\text{SDR=}}{\left[{\sum\limits_{t=1}^N {{{\left| {\hat s\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}{\text{/}}{\left[{\sum\limits_{t=1}^N {{{\left| {s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}$ |
評價指標中的信噪比SNR反映信號受噪聲干擾的程度，其值越小表示受噪聲干擾越強，反之亦然。消噪算法的噪聲抑制率NSR越大，信號失真率SDR越小，說明算法的性能越優。
3.2 基于MCA的ECG分解
實驗數據選用MIT-BIH心率失常數據庫(采樣頻率360 Hz，10 mV/11 bit量化)中相對無污染的ECG片段(編號103，0~2.8 s)。實驗通過添加不同強度的PLI，來驗證算法的性能。原始ECG信號s(t)和受PLI的信號s′(t)及其頻譜如圖 1所示，其中圖 1(a)為s(t)，圖 1(b)為s(t)的頻譜圖，圖 1(c)為受50 Hz PLI的s′(t)(SNR=11.52 dB)，圖 1(d)為s′(t)的頻譜圖。可見，ECG受到了較強的PLI。
圖1
原始ECG和工頻干擾ECG及其頻譜圖
(a)原始ECG信號；(b)原始ECG信號的頻譜圖；(c)受工頻干擾的ECG信號；(d)受干擾ECG的頻譜圖
Figure1. Original and disturbed ECG with its spectrum chart(a) original ECG; (b) spectrum chart of original ECG; (c) disturbed ECG by PLI; (d) spectrum chart of disturbed ECG
設置λ=4對s′(t)進行MCA分解，結果如圖 2所示，圖 2(a)為由突變成分mm(t)，主要包含QRS波群能量；圖 2(b)為平滑成分ms(t)，主要為P波及PLI；圖 2(c)為MCA分解的殘余項，主要包含白噪聲及部分PLI。由圖 2可知，MCA分解后ECG的PLI成分幾乎全部集中于平滑成分，因此要消除PLI只需對平滑成分進行濾波即可。
圖2
受干擾ECG的各形態成分
(a)突變成分；(b)平滑成分；(c)殘余成分
Figure2. Each morphological components of ECG disturbed by PLI(a) mutated component; (b) smooth component; (c) residual component
3.3 基于EEMD的平滑成分分解
EEMD的參數設置如下，白噪聲添加的總次數N為200，白噪聲的幅值系數k為0.01，對MCA分解的平滑成分ms(t)進行EEMD，得到的各階IMF如圖 3所示，各階IMF對應的頻譜圖如圖 4所示。由圖 3、圖 4可知，EEMD將平滑成分ms(t)從高頻到低頻依次進行自適應分解，PLI成分全集中在IMF1分量上, 而且IMF1基本不含其它有用成分。因此，直接濾除其IMF1分量不僅可以方便地濾除PLI，而且不會損傷其它有用的信息。
圖3
平滑成分的各階IMF
Figure3.
Each IMF of smooth component
圖4
各IMF對應的頻譜圖
Figure4.
Spectrum chart of each IMF
3.4 實驗結果
按照公式(9)對圖 1(c)(含PLI的s′(t))進行濾波得到降噪后的ECG信號。與原始ECG(s(t))的對比圖如圖 5所示，其中黑色曲線為s(t)，紅色曲線為。改進的Levkov算法[4]對圖 1(c)的降噪圖與s(t)對比圖如圖 6所示，其中黑色曲線為s(t)，紅色曲線為Levkov的降噪圖。由圖 5、圖 6可知兩種算法都可以濾除PLI。本文算法先用MCA將QRS突變成分分離出來，因此圖 5比圖 6更好地保留了QRS波的細節信息。
圖5
原始ECG與本文算法的降噪對比圖
Figure5.
Comparison between original ECG and ECG de-noised with the proposed algorithm
圖6
原始ECG與改進Levkov算法的降噪對比圖
Figure6.
Comparison between the original ECG and the de-noised signal by improving Levkov algorithm
采用兩種算法對受不同干擾強度β下的ECG進行降噪，實驗仿真具體結果如表 1所示。高信比SNR情況下，由于噪聲較小，即噪聲抑制率NSR相對較低，同時信號失真率SDR也較小，隨著信噪比SNR的減小，噪聲抑制率NSR相對較大，同時信號失真度有所增加，但仍然很小。對比可知，基于本文算法降噪和改進的Levkov算法都有較高的噪聲抑制率NSR和較低的信號失真率SDR，并且本文降噪方法的噪聲抑制率NSR和信號失真率SDR結果都優于改進的Levkov算法。
表1
不同強度工頻干擾ECG的降噪結果
Table1.
De-noising results of ECG with different intensity PLI
4 結束語
本文提出了一種基于MCA和EEMD的ECG PLI消除新算法，充分利用了MCA基于形態差異性分解的特性，將QRS特征波群先提取出來，再利用EEMD的自適應分解特性分離出PLI成分。實驗結果表明，提出的算法可以很有效地濾除ECG信號采集過程中的PLI，同時信號失真率較小，可以較好地保留ECG的特征波形，為ECG信號的進一步處理奠定了良好的基礎。
0 引言
心電信號(electrocardiogram, ECG)在診斷心血管疾病、評估各種治療方法的有效性等臨床應用中有重要意義[1-3]。體表采集到的ECG是一種非線性非平穩的微弱信號，易受到來自外界的各種干擾，主要有工頻干擾(power line interference, PLI)、肌電干擾和基線漂移等，降低了正常診斷率。其主要干擾之一為醫療設備電源PLI，有時其幅度之大甚至可能淹沒需檢測的ECG，且落在ECG有效頻帶內，將直接影響到系統的正常工作以及醫生對疾病的準確診斷。因此，如何有效地抑制PLI噪聲是心電智能診斷的一個重要內容。
目前濾除ECG中PLI的主要方法包括：濾波器濾波法[2-4](如陷波濾波法、自適應濾波法、Levkov濾波法等)和小波變換法[1, 5]。濾波器濾波法通常利用頻帶差異性進行濾波，由于PLI和ECG信號在頻域有交疊，因此濾波的同時會造成ECG特征波形的損失, 從而導致錯誤的診斷結果；小波變換法則存在小波基的選擇和閾值的設定問題。
形態分量分析(morphological component analysis, MCA)根據各組成成分的形態差異性進行分離。參考文獻[6]根據ECG特征波形的形態差異，利用MCA先將ECG分解為突變成分、平滑成分和殘余成分，通過濾除殘余成分有效地濾除白噪聲的干擾，取得了較好的效果，但它無法濾除PLI。由于PLI主要為變化較為平滑的信號，為了避免消除PLI的同時對特征波形(如QRS波群)造成損傷，本文提出了一種基于MCA和集合經驗模態分解(ensemble empirical mode decomposition, EEMD)的消除ECG中PLI的新算法。在文獻[6]的基礎上，先采用MCA提取出突變成分(如QRS波群)，再利用EEMD對平滑成分進行濾波，濾除PLI的本征模態函數(intrinsic mode functions, IMF)，最后重構ECG信號。
1 MCA基本原理
MCA是Starck等[7]近來基于信號的形態多樣性和稀疏表示提出的一種基于多種基函數的稀疏信號分解方法。該方法的主要思想是利用信號組成成分的形態差異性，通過構建不同形態的稀疏表示字典來實現信號中各形態成分的分離，已廣泛應用于圖像處理、醫學信號處理、機械故障診斷等領域[8-10]。
對于任意的觀測信號s∈R，假設s是K個不同形態分量sk的線性組合，即。MCA的目標是從混合信號中恢復出不同形態的各個源信號。這是一個難解的逆問題。MCA假定字典是由一系列基本變換組成，其中源信號sk在變換基Φk上的表示系數αk是稀疏的，而在其他變換基上不能稀疏表示，即字典Φk相互之間不相關。因此字典Φk可以用來區分信號類型，從混合信號中提取出源信號sk。信號s的稀疏分解可歸結為變換系數{α1opt, α2opt, …, αkopt}的優化問題：
| $\left\{ \alpha _{1}^{opt},\alpha _{2}^{opt},\cdots ,\alpha _{\text{K}}^{opt} \right\}=\text{Arg}\underset{\left\{ {{\alpha }_{1}},{{\alpha }_{2}},\cdots ,{{\alpha }_{K}} \right\}k=1}{\mathop{\text{min}}}\,{{\left\| {{\alpha }_{k}} \right\|}_{1}}+\lambda \left\| s-\sum\limits_{k=1}^{K}{{{\Phi }_{k}}{{\alpha }_{k}}} \right\|_{2}^{2}$ |
式中‖αk‖0表示l0范數。放寬式(1)中的約束條件，當信號稀疏分解后的剩余成分類似高斯白噪聲時，可以選擇l2作為誤差范數，則式(1)轉化為：
| $\left\{ {\alpha _1^{opt}, \alpha _2^{opt}, \cdots, \alpha _k^{opt}} \right\}={\text{Arg}}\mathop {{\text{min}}}\limits_{\left\{ {{\alpha _1}, {\alpha _2}, \cdots, {\alpha _K}} \right\}} \sum\limits_{k=1}^K {{{\left\| {{\alpha _k}} \right\|}_1}} + \lambda \left\| {s-\sum\limits_{k=1}^K {{\Phi _k}{\alpha _k}} } \right\|_2^2$ |
式中λ為給定的閾值。可見，MCA方法不僅可以將信號分解成不同的形態分量，還能濾除高斯白噪聲。Starck等[7]在BCR(Block-Coordinate Relaxation)[11]方法的基礎上給出了MCA的數值實現。
2 EEMD基本原理
EEMD是Wu等為了解決模式混疊問題，提出的一種噪聲輔助數據分析方法[12-13]。它完全依賴待檢測信號本身的特征進行模態分解，非常適用于非線性、非平穩信號的處理。
EEMD方法的核心是經驗模態分解(empirical mode decomposition, EMD)。EMD是由Huang等提出的一種自適應的、高效的信號分解方法。它根據信號自身特征將其從高頻到低頻分解成有限個具有物理意義的IMF，這些模態函數能夠內在地描述信號的特征。EMD的基本理論與分解步驟詳見文獻[12]。當待分析信號中某個頻段的分量不連續時，就會導致一個IMF中出現多個不同尺度的信號成分，或相似尺度的成分出現在不同的IMF中，即出現模式混疊現象。為了避免模式混疊，EEMD利用了白噪聲信號具有在各個頻段能量一致的特點，對待分析信號添加白噪聲后再進行EMD分解以保證每個IMF時域的連續性。在待分析信號中添加多次白噪聲后，利用白噪聲均值為零的特性，對各次EMD分解得到的同階IMF分量求取均值，以消除附加白噪聲的影響。
EEMD具體分解步驟如下：
(1)在原始信號s(t)中加入多次幅值系數為k的白噪聲，即：
| ${s_i}\left(t \right)=s\left(t \right) + k \times {n_i}\left(t \right)\; \;\; \left({i=1, 2, 3, \cdots, N} \right)$ |
式中，N為加入白噪聲的總次數，si(t)為第i次加入白噪聲后的信號，ni(t)為第i次加入的白噪聲，N和k通常根據經驗來選取。
(2)依次對si(t)進行EMD分解，分別得到各階的IMFij。IMFij是si(t)分解的第j個IMF，即：
| ${s_i}\left(t \right)=\sum\limits_{j=1}^J {{\text{IM}}{{\text{F}}_{ij}}} $ |
(3)利用不相關隨機序列的統計值為零的原理，將各階相應的IMF進行總體求平均運算便可獲得EEMD的各個IMF分量，即EEMD分解獲得的第j個IMF為IMFj(t)，公式如下：
| ${\text{IM}}{{\text{F}}_j}\left(t \right)=\frac{1}{N}\sum\limits_{i=1}^J {{\text{IM}}{{\text{F}}_{ij}}\left(t \right)} $ |
3 基于MCA和EEMD的工頻干擾消除法
ECG信號的主要特征波形包含突變成分(如QRS波群)和平滑成分，PLI為變化較為平滑的50 Hz正弦波信號，因此消除PLI只需對平滑成分進行濾波即可，從而有效地避免對QRS波群的損傷。對含PLI的平滑成分進行EEMD，平滑成分被自適應地從高頻到低頻依次進行分解。多次實驗結果表明，PLI絕大部分的能量集中在EEMD的第一個本征模態函數IMF1上, 而且IMF1分量上幾乎全為PLI，因此濾除IMF1和MCA分解的白噪聲成分后重構ECG，即可達到消除PLI的目的，同時濾除了白噪聲。
假設無污染的ECG信號為s(t)，PLI的強度為β，受PLI的ECG信號為s′(t)，消噪后的ECG信號為。根據ECG的特征波形，采用非抽樣小波字典來稀疏表示s′(t)的突變成分mm(t)，采用離散余弦變換字典來稀疏表示s′(t)的平滑成分ms(t)，MCA分解的殘余噪聲成分用r(t)表示，則添加噪聲、MCA、EEMD及濾除PLI模型分別為式(6)、式(7)、式(8)、式(9)。
| $s'\left(t \right)=s\left(t \right) + \beta \times \sin \left({2\pi \times 50t} \right)$ |
| $s'\left(t \right)={m_{\text{m}}}\left(t \right) + {m_{\text{s}}}\left(t \right) + r\left(t \right)$ |
| ${m_s}\left(t \right)=\sum\limits_{i=1}^J {{\text{IM}}{{\text{F}}_i}\left(t \right)} $ |
| $\hat s\left(t \right)={m_{\text{m}}}\left(t \right) + \sum\limits_{i=2}^J {{\text{IM}}{{\text{F}}_i}\left(t \right)} $ |
3.1 算法評價指標
為驗證算法的有效性，采用信噪比(signal noise ratio, SNR)、噪聲抑制率(noise suppression ratio, NSR)和信號失真率(signal distortion ratio, SDR)作為評估指標，其定義如下：
| ${\text{SNR=10}} \times {\log _{10}}\left({\sum\limits_{t=1}^N {{s^2}\left(t \right)} /\sum\limits_{t=1}^N {{n^2}\left(t \right)} } \right)$ |
| ${\text{NSR=1-}}{\left[{\sum\limits_{t=1}^N {{{\left| {\hat s\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}{\text{/}}{\left[{\sum\limits_{t=1}^N {{{\left| {s'\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}$ |
| ${\text{SDR=}}{\left[{\sum\limits_{t=1}^N {{{\left| {\hat s\left(t \right)-s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}{\text{/}}{\left[{\sum\limits_{t=1}^N {{{\left| {s\left(t \right)} \right|}^2}} } \right]^{{\text{0}}{\text{.5}}}}$ |
評價指標中的信噪比SNR反映信號受噪聲干擾的程度，其值越小表示受噪聲干擾越強，反之亦然。消噪算法的噪聲抑制率NSR越大，信號失真率SDR越小，說明算法的性能越優。
3.2 基于MCA的ECG分解
實驗數據選用MIT-BIH心率失常數據庫(采樣頻率360 Hz，10 mV/11 bit量化)中相對無污染的ECG片段(編號103，0~2.8 s)。實驗通過添加不同強度的PLI，來驗證算法的性能。原始ECG信號s(t)和受PLI的信號s′(t)及其頻譜如圖 1所示，其中圖 1(a)為s(t)，圖 1(b)為s(t)的頻譜圖，圖 1(c)為受50 Hz PLI的s′(t)(SNR=11.52 dB)，圖 1(d)為s′(t)的頻譜圖。可見，ECG受到了較強的PLI。
圖1
原始ECG和工頻干擾ECG及其頻譜圖
(a)原始ECG信號；(b)原始ECG信號的頻譜圖；(c)受工頻干擾的ECG信號；(d)受干擾ECG的頻譜圖
Figure1. Original and disturbed ECG with its spectrum chart(a) original ECG; (b) spectrum chart of original ECG; (c) disturbed ECG by PLI; (d) spectrum chart of disturbed ECG
設置λ=4對s′(t)進行MCA分解，結果如圖 2所示，圖 2(a)為由突變成分mm(t)，主要包含QRS波群能量；圖 2(b)為平滑成分ms(t)，主要為P波及PLI；圖 2(c)為MCA分解的殘余項，主要包含白噪聲及部分PLI。由圖 2可知，MCA分解后ECG的PLI成分幾乎全部集中于平滑成分，因此要消除PLI只需對平滑成分進行濾波即可。
圖2
受干擾ECG的各形態成分
(a)突變成分；(b)平滑成分；(c)殘余成分
Figure2. Each morphological components of ECG disturbed by PLI(a) mutated component; (b) smooth component; (c) residual component
3.3 基于EEMD的平滑成分分解
EEMD的參數設置如下，白噪聲添加的總次數N為200，白噪聲的幅值系數k為0.01，對MCA分解的平滑成分ms(t)進行EEMD，得到的各階IMF如圖 3所示，各階IMF對應的頻譜圖如圖 4所示。由圖 3、圖 4可知，EEMD將平滑成分ms(t)從高頻到低頻依次進行自適應分解，PLI成分全集中在IMF1分量上, 而且IMF1基本不含其它有用成分。因此，直接濾除其IMF1分量不僅可以方便地濾除PLI，而且不會損傷其它有用的信息。
圖3
平滑成分的各階IMF
Figure3.
Each IMF of smooth component
圖4
各IMF對應的頻譜圖
Figure4.
Spectrum chart of each IMF
3.4 實驗結果
按照公式(9)對圖 1(c)(含PLI的s′(t))進行濾波得到降噪后的ECG信號。與原始ECG(s(t))的對比圖如圖 5所示，其中黑色曲線為s(t)，紅色曲線為。改進的Levkov算法[4]對圖 1(c)的降噪圖與s(t)對比圖如圖 6所示，其中黑色曲線為s(t)，紅色曲線為Levkov的降噪圖。由圖 5、圖 6可知兩種算法都可以濾除PLI。本文算法先用MCA將QRS突變成分分離出來，因此圖 5比圖 6更好地保留了QRS波的細節信息。
圖5
原始ECG與本文算法的降噪對比圖
Figure5.
Comparison between original ECG and ECG de-noised with the proposed algorithm
圖6
原始ECG與改進Levkov算法的降噪對比圖
Figure6.
Comparison between the original ECG and the de-noised signal by improving Levkov algorithm
采用兩種算法對受不同干擾強度β下的ECG進行降噪，實驗仿真具體結果如表 1所示。高信比SNR情況下，由于噪聲較小，即噪聲抑制率NSR相對較低，同時信號失真率SDR也較小，隨著信噪比SNR的減小，噪聲抑制率NSR相對較大，同時信號失真度有所增加，但仍然很小。對比可知，基于本文算法降噪和改進的Levkov算法都有較高的噪聲抑制率NSR和較低的信號失真率SDR，并且本文降噪方法的噪聲抑制率NSR和信號失真率SDR結果都優于改進的Levkov算法。
表1
不同強度工頻干擾ECG的降噪結果
Table1.
De-noising results of ECG with different intensity PLI
4 結束語
本文提出了一種基于MCA和EEMD的ECG PLI消除新算法，充分利用了MCA基于形態差異性分解的特性，將QRS特征波群先提取出來，再利用EEMD的自適應分解特性分離出PLI成分。實驗結果表明，提出的算法可以很有效地濾除ECG信號采集過程中的PLI，同時信號失真率較小，可以較好地保留ECG的特征波形，為ECG信號的進一步處理奠定了良好的基礎。
