第三章 傳統濾波方法與研究理論
3.2 經驗模態分解法
3.2.3 停止準則
(5) 頻率相近之分量無法被解析[39]。
(a)本質模態函數數目為 8 個
(b)本質模態函數數目為 9 個
圖3.5 不同停止準則所得到之不同本質模態函數
Mixer
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -2
0 2
EMD
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -0.5
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 5.5 6
Mixer
濾波器分別由帶通濾波器(Bandpass Filter)、微分器與積分器(Integrator)所 構成。
ECG original
data
Bandpass filter
Differentiator
Squaring process
Moving window integration
End
圖3.7 QRS 綜合波檢測流程
(1) 帶通濾波器
QRS 檢測方法第一個步驟即是使用帶通濾波器,目的主要是能減低 肌肉波、60Hz市電干擾、基準線漂移及T 波所造成的影響與誤判。而理 想中,帶通濾波器參數設定上下為5~15Hz。
(2) 微分器
利用微分器的目的主要是希望能得到QRS 綜合波的斜率資訊。
(3) 平方處理
經過此程序,使得所有資料點皆變為正且非線性的放大,目的主要突 顯訊號中頻率較高的成分。
)]2
( [ )
(nT x nT
y (3.22) (4) 移動的窗型積分
可獲取R 波斜率之外的波形特徵。
經由上述數位濾波器的一系列結合與處理方法,再配合閥值的選定,
最終會輸出R 波脈衝波形,也就是可得到 QRS 波的所在位置資訊。
圖3.8 移動窗型積分器與 QRS 綜合波的關係
第四章 研究步驟與結果
4.1 MIT/BIH 心律不整資料庫簡介[43]
心律不整的發生主要是由於心電訊號觸發點及電訊號傳導路徑的異 常所產生,其定義為連續 RR 間距變化大於 0.16 秒即為心律不整。心律 不整為常見之心血管疾病,目前最有效偵測方式還是以心電圖記錄最為 精準正確。
而 MIT/BIH 心律不整資料庫為一全世界普遍使用之標準心電訊號測 試資料庫,其可做為心律不整檢測器效能評估之用。MIT/BIH 心律不整 資料庫目前共有48 筆資料,每筆資料含有兩個導程之記錄。目前網路上 所能取得之完整心律資料(30 分鐘)有 25 筆,而不完整心律資料(10 分鐘) 則有23 筆。而本研究將從 MIT/BIH 心律不整資料庫中,擷取當中三筆心 電訊號實例記錄做為演算法測試對象。
4.2 心電訊號雜訊濾除 4.2.1 心電訊號之分解
首先,本研究從MIT/BIH 心律不整資料庫中,選取編號 103 第二通 道、104 第二通道與 109 第一通道之心電訊號實例做為演算法測試對象。
接下來第一步驟即為心電訊號的解析,而由前章介紹得知經驗模態分解 法對於非線性與非穩態訊號有良好的解析度,相當適合做為心電訊號之 分解方法,因此本研究將使用前章所述之經驗模態分解法來拆解心電訊 號。本研究將所選取之心律不整心電資料經由經驗模態分解法拆解後,
可得到如圖所示之結果。
Channel 2
0 1 2 3 4 5 6 7 8
time ( sec ) 0
200
EMD
0 1 2 3 4 5 6 7 8
time ( sec ) -100
0 100 200 300 400 500 600 700 800 900 1000
圖4.1 MIT/BIH 資料庫編號 103 第二通道心電圖資料
Channel 2
0 1 2 3 4 5 6 7 8
time ( sec ) -200
0 200
EMD
0 1 2 3 4 5 6 7 8
time ( sec ) -100
0 100 200 300 400 500 600 700 800 900 1000 1100
圖4.2 MIT/BIH 資料庫編號 104 第二通道心電圖資料
Channel 1
0 1 2 3 4 5 6 7 8
time ( sec ) -200
0 200
EMD
0 1 2 3 4 5 6 7 8
time ( sec ) -200
0 200 400 600 800 1000 1200
圖4.3 MIT/BIH 資料庫編號 109 第一通道心電圖資料
4.2.2 本質模態函數之選擇
對於所拆解出之本質模態函數,首先針對各個本質模態函數做傅立葉 轉換,目的是將訊號由時域轉至頻域觀察。再由前面章節所提至基本雜 訊號所存在之頻段,如高頻的市電干擾為60Hz與低頻的基準線漂移約為
3Hz .
0 等,即可得知經驗模態分解法所拆解出之分量何者為雜訊,又何者
為重建之分量。
Channel 1-FFT
Channel 2-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
2 4
Channel 3-FFT
0 20 40 60 80 100 120 140 160 180
Channel 4-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
5 10
Channel 5-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
8
Channel 6-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
5 10
Channel 7-FFT
Channel 8-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
2 4
Channel 9-FFT
0 20 40 60 80 100 120 140 160 180
Channel 1-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
1 2
Channel 3-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
8
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
5 10
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
15
Channel 6-FFT
0 10 20 30 40 50 60 70
Channel 7-FFT
0 10 20 30 40 50 60 70
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 9-FFT
Channel 1-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
0.5 1
Channel 2-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
6
Channel 3-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
5 10
Channel 4-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
10 20
Channel 5-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
10 20
Channel 6-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
5 10
Channel 7-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
5 10
Channel 8-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
10
Channel 9-FFT
0 20 40 60 80 100 120 140 160 180
frequency ( Hz ) 0
20
圖4.6 MIT/BIH 資料庫編號 109 第一通道心電圖資料各個分量在頻域的分 布情形
4.3 心電訊號之重建 利用MSE(mean square error)來比較最小平方法與部分重建法的優劣。首 先,我們假設訊號源為 態函數。由傅立葉轉換(Fourier transform, FT)可觀察各個本質模態函數所包 含的頻率成份,發現IMF3~IMF4皆為重建的分量。
Mixer
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -10
0 10
EMD
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) 0
10 20 30 40 50 60
圖4.7 Case1(k=0.05)
Channel 1-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.06
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.02 0.04
Channel 3-FFT
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.6
Channel 6-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2
Channel 7-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.8
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
Channel 9-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.15
圖4.8 Case1各個分量在頻域上之分佈情形
Case2:k=0.1
圖4.9 為 k=0.1 時,訊號源經由經驗模態分解法所拆解出的各個本質 模態函數。由傅立葉轉換可觀察各個本質模態函數所包含的頻率成份,
發現IMF3~IMF4 皆為重建的分量。
Mixer
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -10
0 10
EMD
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) 0
10 20 30 40 50 60 70
Channel 1-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.06
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.15
Channel 3-FFT
0 10 20 30 40 50 60 70
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 6-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.6
Channel 7-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
Channel 9-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
圖4.10 Case2 各個分量在頻域上之分佈情形
Case3:k=0.2
圖 4.11 為 k=0.2 時,訊號源經由經驗模態分解法所拆解出的各個本 質模態函數。由傅立葉轉換可觀察各個本質模態函數所包含的頻率成 份,發現IMF3~IMF4 皆為重建的分量。
Mixer
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -10
0 10
EMD
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) 0
10 20 30 40 50 60 70
圖4.11 Case3(k=0.2)
Channel 1-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.05
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1
Channel 3-FFT
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.6
Channel 6-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.3
Channel 7-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
Channel 9-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1 0.2
圖4.12 Case3 各個分量在頻域上之分佈情形
Case4:k=0.5
圖 4.13 為 k=0.5 時,訊號源經由經驗模態分解法所拆解出的各個本 質模態函數。由傅立葉轉換可觀察各個本質模態函數所包含的頻率成 份,發現IMF3~IMF4 皆為重建的分量。
Mixer
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) -10
0 10
EMD
0 1 2 3 4 5 6 7 8 9 10
time ( sec ) 0
10 20 30 40 50 60 70
圖4.13 Case4(k=0.5)
Channel 1-FFT
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1 0.2
Channel 3-FFT
0 10 20 30 40 50 60 70
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 6-FFT
0 0.2 0.4
Channel 7-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 9-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
圖4.14 Case4 各個分量在頻域上之分佈情形
Case5:k=2
圖4.15 為 k=2 時,訊號源經由經驗模態分解法所拆解出的各個本質 模態函數。由傅立葉轉換可觀察各個本質模態函數所包含的頻率成份,
發現IMF2~IMF4 皆為重建的分量。
Mixer
0 1 2 3 4 5 6 7 8 9 10
tim e ( sec ) -10
0 10
EMD
0 1 2 3 4 5 6 7 8 9 10
tim e ( sec ) 0
10 20 30 40 50 60 70
圖4.15 Case5(k=2)
Channel 1-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
1 2
Channel 2-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
Channel 3-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
1.5
Channel 4-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
2 4
Channel 5-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.2 0.4
Channel 6-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1 0.2
Channel 7-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1 0.2
Channel 8-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.1 0.2
Channel 9-FFT
0 10 20 30 40 50 60 70
frequency ( Hz ) 0
0.5 1
圖4.16 Case5 各個分量在頻域上之分佈情形
在最小平方法與部分重建法的比較上(表 4.1),由於最小平方法可針 對每個本質模態函數乘上不同的權重係數,因此重建後的訊號效果較 佳。另外,我們發現在雜訊干擾小的情況下,所計算出的均方誤差是使 用最小平方法可以達到較佳的重建效果;當雜訊干擾較為強烈時(k=2),
使得來源訊號較為接近雜訊,而直接將雜訊分量刪除即可獲得接近原始 訊號的重建訊號,因此採用較簡單的部分重建法即可獲得與較複雜的最 小平方法同樣甚至較佳的效果。
表4.1 最小平方法與部分重建法比較表
4.3.2 重建無干擾訊號
由上述比較得知,在重建步驟上使用最小平方法可以獲得較佳的效 果,因此本研究即採用了最小平方法做為原始心電訊號重建的基礎理 論。在MIT/BIH 心律不整資料庫編號 103 之資料經由經驗模態分解法之
分量所存在的頻段為何,發現IMF2~IMF5 皆為重建之分量,再利用最小 平方法求出每個分量所需乘上之權重係數,最後將成上權重係數之各個 分量相加總合在一起,結果如圖所示。而MIT/BIH 心律不整資料庫標號 104 第二通道,IMF1~IMF4 為其重建分量,圖 4.12 為重建波形;MIT/BIH 心律不整資料庫標號109 第一通道中,IMF2~IMF5 為其重建分量,圖為 重建波形。圖中顯示,經由上述濾波之步驟,發現能有效濾除高頻市電 干擾雜訊與低頻基準線漂移雜訊。
圖4.17 MIT/BIH 編號 103 第二通道重建後波形(藍色為原始訊號,紅色為重 建訊號)
圖4.18 MIT/BIH 編號 104 第二通道重建後波形(藍色為原始訊號,紅色為 重建訊號)
圖4.19 MIT/BIH 編號 109 第一通道重建後波形(藍色為原始訊號,紅色為 重建訊號)
4.4 QRS 檢測結果
在雜訊濾除之後,接著即為計算RR 間距並藉以判斷是否存在心律不 整的病症。本研究利用 3.4 小節所提至之 QRS 檢測方法並自行設計 MATLAB 圖形使用者介面(GUI),針對雜訊濾除後的心電訊號進行診斷,
結果皆顯示為心律不整,此結果與本研究擷取MIT/BIH 心律不整的原始 心電訊號相互吻合。
圖4.20 MIT/BIH 編號 103 第二通道 QRS 檢測結果為心律不整
圖4.21 MIT/BIH 編號 104 第二通道 QRS 檢測結果為心律不整
圖4.22 MIT/BIH 編號 109 第一通道 QRS 檢測結果為心律不整
第五章 結論與未來展望
本文應用經驗模態分解法針對非穩態與非線性訊號之良好解析度,以 及最小平方法在數學領域上之最佳化技術成功有效濾除心電訊號中所存 在之雜訊成分。在訊號分解上,經驗模態分解法能將訊號中所隱含的成 分由高頻至低頻拆解出來,其中本質模態函數屬於高頻部份,殘餘值則 屬於低頻部份,因此整個訊號拆解的過程為一濾波的動作。而重建的部 份則利用最小平方法來進行,先求出欲重建的各個本質模態函數需要相 乘的權重係數,再乘上相對應的本質模態函數並相加之後,即可得到重 建訊號。本研究亦在最小平方法與部分重建法之間利用均方誤差針對訊 號重建效果做了相關比較,最後實驗驗證雜訊干擾小的情況下,使用最 小平方法可以達到較佳的重建效果;當雜訊干擾較為強烈時,訊號源會 近似雜訊,則採用過程較簡單的部分重建法即可獲得與過程較複雜的最 小平方法同樣甚至更佳效果。
本文所提出之雜訊濾除方法雖然已能有效的濾除心電訊號中雜訊成 分,但經驗模態分解法還存在著許多限制與缺點,期望未來能有更多研 究者針對此領域繼續深入探討及改善演算法,為全民福祉做出更大的貢 獻。
參考文獻
[1] 行政院衛生署
http://www.doh.gov.tw/cht2006/index_populace.aspx [2] 基本心電圖圖例簡介
http://www.dcbiomed.com/material/ECG3CH.pdf
[3] J.A.V. Alste, T.S. Schilder, “Removal of Base Line Wander and Power Line Interference from the ECG by an Efficient FIR Filter with a Reduced Number of Taps”, IEEE transactions on biomedical engineering, Vol. Bme-32. No. 12, Dec. 1985, pp.1052-1060.
[4] A. Gutierrez, P.R. Hernandez, M. Lara, S. Perez, “A QRS Detection Algorithm Based on Haar Wavelet”, Computers in Cardiology 1998 Sep. 1998, pp.353-356.
[5] A. Haar, “Zur theorie der orthogonalen funktionen-systeme”, Math. Ann., vol. 69, 1910, pp.331-373.
[6] I. Daubechies, Ten Lectures on Wavelets, SIAM, 1992.
[7] H. Xing, M. Huang, “A New QRS Detection Algorithm Based on Empirical Mode Decomposition”, Bioinformatics and biomedical engineering, May 2008, pp.693-696.
[8] D. Davis (1991), 基本心電圖判讀 (黃天守、陳清輝編譯),台北:眾 文圖書股份有限公司。(原著出版於 1985)
[9] http://www.adam.com/democontent/hie/images/en/1135.jpg
[10] http://static.howstuffworks.com/gif/adam/images/en/ecg-electrode-place ment-picture.jpg
[11] G.M. Friesen, T.C. Jannett, M.A. Jadallah, S.L. Yates, S.R. Quint, H.T.
Algorithm”, IEEE transactions on biomedical engineering, Vol. 37. No.
1. January, 1990, pp.85-98.
[12] G.D. Clifford, F. Azuaje, P.E. mcSharry, Advanced method and tools for ECG data analysis, Artech House, INC. 685 Canton Street Norwood.
[13] S.K. Mitra, Digital signal processing, McGraw.Hill international 3rd edtion 2006.
[14] A.G. Kleber, M.J. Janse, F.J. van Capelle, D. Durrer, “Mechanism and time course of S-T and T-Q segment changes during acute regional ischemia in the pig determined by extracellular and intracellular recordings”, Circ. Res., Vol. 42, 1978, pp.603-613.
[15] R. Coronel, J.R. DE Groot, M.J. Janse, “Laplacian electrograms and the interpretation of complex ventricular activation patterns during ventricular fibrillation”, J Cardiovasc Electrophysiol, Vol. 11, 2000, pp.1119-1128.
[16] N. E. Huang, Z. Shen, S. R. Long, M.C. Wu, H.H. Shih, Q. Zheng, N.C.
Yen, C.C. Tung, H.H. Liu, “The Empirical Mode Decomposition and
Yen, C.C. Tung, H.H. Liu, “The Empirical Mode Decomposition and