Abstract: QRS complex detecting algorithm was core of ECG auto-diagnosis method, heart rate variability analysis and deeply influences cardiac cycle division for signal compression. However, ECG signals collected by noninvasive surface electrodes were confused by several kinds of noise and its waveform variation was the main reasons for the hard realization of 100% detection accruracy. QRS complex detecting algorithms based on mixed methods were studied. This study proposed a QRS complex detecting algorithm based on wavelet transform and multi-resolution mathematical morphological decomposition (WMR algorithm). This algorithm possessed superiorities in R peak detection of the two methods. Moreover, a pre-processing method based on lifting scheme constructing multi-resolution morphological decomposition was adopted to reduce noise affection. And an efficient R peak search-back algorithm was employed to reduce the False Positives (FP) and False Negatives (FN). According to simulation results in MIT-BIH Arrhythmia Database, QRS detection accuracy was over 99.8%.
INTRODUCTION
QRS complex detection was the most important stage in ECG analyzing process. Accurate detection of QRS complex was the precondition of detecting other characteristic waves. This was the basis of ECG automatic classification and diagnosis (Hendel et al., 2010; El-Asir and Mamlook, 2002). Meanwhile, it ensured dividing cardiac cycles correctly which was usually the first step before ECG compression and had deep impact on utilization of correlations between cardiac cycles or different leads. The QRS complex detection consisted of R peak detection and onset/offset decision and R peak detection was focused in this study.
There were many algorithms designed for QRS detection purpose. Before wavelet techniques were adopted, the derivative and its related information of ECG were mainly used for QRS detection, as it could access the deep slops of the R waves (Pan and Tompkins, 1985; Bereksi-Reguig and Slimane, 2000). Such kind of methods had the characteristics of less calculation and easy implementation but with lower positive rate and noise immunity. The QRS detection algorithms based on wavelet techniques had become popular since last decade. The main idea was that QRS complex’s analyzing components were different between levels and by choosing the level with local maximum, it was easier to distinguish modulus maximum which were related to QRS complexes (Gautam and Sharma, 2010; Chouakri et al., 2011, Zidelmal et al., 2012). However, for QRS complex with steep slope, methods based on wavelet techniques were likely to mix it with EMG or other noise and resulted in leak detection. Mathematical Morphology (MM) methods had been introduced into ECG signal processing field in 1980s (Chu and Delp, 1989). Mostly, algorithms based on MM methods were adopted in the pre-procession for its robustness and self-adaptibility in extracting morphological information (Sun et al., 2002). And then QRS detection based on MM method was mainly implemented by triangle Structure Elements (SE) (Trahanias, 1993; Zhang et al., 2007). Goutsias and Heijmans (2000) proposed Nonlinear Multi-resolution Signal Decomposition Scheme based on MM which could achieve good performance for ECG signals procession, especially in ECG pre-procession (Zhang et al., 2013). The QRS complex detection could also be realized by Multi-resolution Decomposition Scheme based on MM (MMMD) (Hu and Bao, 2010). And QRS complex detection methods based on MM, distinguished steep R peaks more easily than gentle ones.
This study proposed a QRS complex detection method based on lifting scheme constructing MMMD and wavelet decomposition which was implemented after ECG signal pre-processing by LMW algorithm described by Zhang et al. (2013). Furthermore, an effective R peak search-back scheme was adopted to reduce the False Positives (FP) and False Negatives (FN).
METHODOLOGY
This proposed QRS complex detection algorithm (WMR algorithm) was based on wavelet transform and lifting scheme constructing MMMD which was also adopted in previous pre-processing method. The basic idea was described as follows.
Lifting scheme constructing MMMD: Define ECG signal as f(n), n = 0, 1…, N-1, symmetrical SE as B(m), m = 0, 1, 2,…., M-1, erosion and dilation of ECG signal could be denoted as following:
| (1) |
| (2) |
Opening operator was
According to MMMD theory, ECG signals at level j could be decomposed as following:
| (3) |
| (4) |
where, MFj(f) was defined as multi-resolution morphological filter at level j:
| (5) |
where, Bj was a linear SE with length equaled to j+1 while xj was the approximate signal and yj was the detail signal at level j.
More improvements in specified performance of algorithms such as noise suppression and QRS detection accuracy and so on, were able to be achieved by adopting lifting scheme constructing method with appropriate operators. Consider lifting scheme constructing new detail and approximate component at level j as:
| (6) |
| (7) |
where, π: Vj→Wj was the predicting operator and λ: Wj→Vj was the updating operator. Reconstructions at level j and j-1 were defined as:
| (8) |
| (9) |
WMR method for QRS complex detection: WMR method for QRS complex detection could be summarized by application of the following steps and the flow chart was illustrated in Fig. 1:
| Step 1: | ECG signal pre-processing by LMW algorithm, noise suppression in ECG signals could be divided into two steps: Remove the baseline wandering caused by electrodes movement and respiration, suppress high frequency noise caused by EMG and other interference. Therefore, LMW algorithm first did baseline wandering correction by using the following MM filter (Chu and Delp, 1989): |
| Fig. 1: | Work flow of the proposed QRS complex detection method |
| (10) |
where, f0 was the original ECG signal, fb was the detected baseline wandering, B0 and Bc were two linear SEs which depended on the duration of characteristic waveforms in ECG signals Tω and sampling rate Fs. Commonly, Tw<0.2 sec, therefore, let B0 = 0.2 Fs and Bc = 1.5 B0 = 0.3 Fs. The ECG signals after baseline wandering correction could be decomposed at level j as following:
| (11) |
and let:
x1 = x0
Then high frequency noise suppression was conducted by lifting scheme constructing MMMD. Block effect was the main cause of waveform distortion in morphological filter based ECG pre-processing method. To address this issue, an updating and predicting operator based on cubic spline interpolation was proposed as:
| (12) |
| (13) |
where, the predicting operator was the difference between adjacent samples and updating operator was the interpretation according to the forward and backward sample. Because of the smoothing characteristic of cubic spline interpolation, LMW algorithm could reduce the block effect and improve the waveform distortion. The good performance in noise suppression is also proved by previous study.
Take sample No. 222 in MIT-BIH Arrhythmia Database with visible baseline wandering and high frequency noise as example, as was shown in Fig. 2, baseline wandering was corrected firstly, then ECG signals were decomposed to level 2 by lifting scheme constructing MMMD, finally approximation component at level 2 was chosen as the pre-processing results and original data for QRS detection:
| Step 2: | Segment ECG data every 10 sec as a result of compromise between data size and the number of calculation |
| Step 3: | Detect QRS complex by db6 wavelet decomposition. Apply wavelet decomposition up to level 4 by db6 wavelet. Specific to ECG signals in MIT-BIH Arrhythmia Database with sampling rate of 360 Hz, frequency component of QRS complex mainly focused on detail signals on level 3 and 4 (D3 and D4). Therefore, modulus maximum of D3 and D4 and the corresponding zero crossing points were searched and fed back to original signals as R peaks |
| Step 4: | Detect QRS complex by lifting scheme constructing MMMD with maximization operator, using results from ECG pre-procession. Let predicting operator was the difference between adjacent samples and updating operator was the maximum of adjacent three samples, as describe in Eq. 14 and 15: |
| (14) |
| (15) |
| and ECG signals could be decomposed through 4-10. Decompose the signals to level 4, (i.e., do twice lifting scheme constructing with maximization operator after pre-processing). Locate modulus maxima of detail component at level 3 (D3) and detail component at level 4 (D4). Also take sample No. 222 for example, D3 and D4 were shown in Fig. 3. Then feed them back to original signals as R peaks | |
| Step 5: | Do “OR” operation between the R-peak positions from step 3 and 4, if the result equaled to 1, reserve it as new R-peak positions |
| Step 6: | R peak search back scheme was adopted to reduce false and leak detection and the study flow was illustrated in Fig. 4. Set initial value of threshold as one-tenth the maximum of samples in the segment. Define R was RR interval and aveR was the average RR intervals of two preceding segments. It was critical for both medical staff and researchers to restrict R and aveR in reasonable ranges. Jezior et al. (2005), presented that R<250 msec is an important symbol of sudden death and by Bjerregaard (1983), it was proved that only 5% of patients have R>1750 msec occasionally. Therefore, reasonable range of 250 msec<R<1750 msec was considered. Because the sampling rate of MIT-BIH Arrhythmia Database was 360 Hz, the limits of R were set as 263 and 1675 msec which meant if R<263 msec, there was a False Peak (FP) and if R>1675 msec, there was a leak detection (FN). Define Ra = R/aveR, if Ra>1.66, there was a FN. When there was a FP, the peak was deleted; when there was a FN, then the threshold was reduced by half and the search back algorithm was conducted one more time, until the values of R and Ra were reasonable and R peaks would be decided |
| Fig. 2(a-f): | Sample No. 222 in MIT-BIH Arrhythmia Database and pre-processed by LWM algorithm. Take
approximation component at level 2 as preprocessing results (a) Sample No. 222, (b) ECG signal after
baseline wander correction, (c) Approximation component at level 1, (d) Detail component at level 1,
(e) Approximation component at level 2 and (f) Detail component at level 2 |
| Fig. 3(a-b): | Detail components at (a) Level 3 and (b) Level 4 (D3 and D4) of sample No. 222 |
| Fig. 4: | Work flow of R peak search back scheme |
| Step 7: | Locate onset/offset of QRS complex. According to step 3, choose extremum nearest the zero crossing points as the onset/offset of QRS complex |
| Step 8: | Update threshold. Set initial value of threshold as one-tenth the maximum of samples in the segment for step 4. For step 3, let θ3 be D3 threshold, define it as Root of Mean Square (RMS) of previous 3 sec samples including the sample; while θ3 be D4 threshold, define it as half the RMS of previous 3 sec samples including the sample; as were described in Eq. 15 and 16. Where W3x[m] and W4x[m] were ECG components at D3 and D4, k was ECG sample number per second. For MIT-BIH Arrhythmia Database, k = 360: |
| (16) |
| (17) |
RESULTS
The proposed QRS complex detection algorithm (WMR algorithm) was applied to the first channel of 48 ECG signals in MIT-BIH Arrhythmia Database. Each signal is of length 30 min and sampled at 360 Hz. Three evaluation index have been chosen: FP, FN and sensitivity that was define as follows:
where, RR Numbers denoted overall RR Numbers, i.e., heart beats in the sample. The algorithm was compared with Method 1 (Bahoura et al., 1997) based on wavelet decomposition, Method 2 (Adnane et al., 2009) based on positive and negative slope of QRS complex, Method 3 (Ning and Selesnick, 2013) based on sparse derivatives, Method 4 (Arzeno et al., 2008) based on first-derivative and Method 5 (Lee et al., 1996) based on topological mapping.
The performances of above algorithms were shown in Table 1.
| Table 1: | R peak detection results of proposed method and comparison with 5 other algorithms applied to MIT-BIH arrhythmia database |
Although, the overall heart beats for simulation were slightly different, it was obvious to conclude that, compared with other algorithms; the proposed algorithm achieved better performance on FP, FN and sensitivity. Not only the improving QRS detection method but also the good pre-processing performance contributed to this priority.
DISCUSSION
Samples with abnormal waveforms and serious noise which would cause FP and FN, were also focused. Sample No. 104, 108, 200, 203, 208 and 228 for instance were chosen as examples to inspect the advantages of the proposed method and failure analysis, as illustrated in Fig. 5.
Figure 6a showed segment of sample No. 104 with Premature Beats (PB). After pre-processing, the periodicity became stronger, the waveform pattern of R peaks were distinguish with other characteristic waves. Therefore, FP and FN were hard to find.
Figure 6b showed segment of sample No. 108 with abnormal P waves.
| Fig. 5(a-f): | Abnormal waveforms information in (a) Sample No. 104 with PB, (b) Sample No. 108 with abnormal
P waves, (c) Sample No. 200 with QRS MC, (d) Sample No. 203 with iRR and QRS MC,
(e) Sample No. 208 with fVN and (f) Sample No. 228 with bcR |
| Fig. 6(a-f): | Segment of R peak detection results in sample No. (a) 104, (b) 108, (c) 200, (d) 203, (e) 208 and (f) 228 |
Generally, high P waves did not cause false detection (FP), but when RR intervals became longer suddenly and high P waves existed, it was possible to cause FP, as samples between 400 and 500. On the other hand, for samples around 2000, two adjacent R peaks were very close, RR = 230 msec, Ra ≈ 1.724. Therefore, if there was no search back process, a leak detection (FN) would definitely occur.
Figure 6c showed segment of sample No. 200 with Ventricular Premature Contractions (VPC) in the 2nd, 4th and 6th RR intervals. The slopes of R peaks were big enough and RR intervals were regular. Therefore, WMR algorithm could ensure the right detection.
Figure 6d showed segment of sample No. 203 with irregular RR intervals (iRR) and QRS Morphological Changes (QRS MC). Obviously, irregular RR intervals and noise interference were the reason that causes the first FN and amplitude changes rapidly caused the second FN.
Figure 6e showed segment in V1 lead of sample No. 208 with fusion ventricular heartbeats (fVN) which was more typical for analysis. The adjacent R peak amplitudes varied. Meantime the RR intervals were irregular but the sum of adjacent RR intervals was no more than 1675 msec. Therefore, when R peak amplitudes were high and led to higher threshold, FN would occur.
Figure 6f showed segment of sample No. 228 with broad changes of R peaks (bcR). Although, there was one extremely high R peak around sample 2000, it had no critical affection to threshold for its discontinuity. And QRS waveforms had good patterns which ensured good detection accuracy.
In conclusion, the proposed WMR algorithm had robustness to abnormal waveforms and irregular RR intervals and decreased FP and FN to some extent.
CONCLUSION
A novel QRS complex detection method based on lifting scheme constructing MMMD and wavelet transformation was proposed. Cooperating with good performance pre-procession and effective R peak search back strategy, this method achieved sensitivity of 99.8% as applied to MIT-BIH Arrhythmia Database and possessed certain robustness to abnormal ECG signals. Furthermore, it didn’t add much more calculation complexity than methods based on wavelet or other multi-resolution decomposition methods. Therefore, it could be adopted by mobile medical devices using in Telemedicine system as part of auto-diagnosis or compressing and transmission design in the future.
ACKNOWLEDGMENT
This project is supported by Communication Engineering Research Center of Shenzhen Graduate School of Harbin Institute of Technology and Konaka Laboratory of The University of Tokushima.