Abstract
This paper presents an efficient electrocardiogram (ECG) signals compression technique based on QRS detection, estimation, and 2D DWT coefficients thresholding. Firstly, the original ECG signal is preprocessed by detecting QRS complex, then the difference between the preprocessed ECG signal and the estimated QRS-complex waveform is estimated. 2D approaches utilize the fact that ECG signals generally show redundancy between adjacent beats and between adjacent samples. The error signal is cut and aligned to form a 2-D matrix, then the 2-D matrix is wavelet transformed and the resulting wavelet coefficients are segmented into groups and thresholded. There are two grouping techniques proposed to segment the DWT coefficients. The threshold level of each group of coefficients is calculated based on entropy of coefficients. The resulted thresholded DWT coefficients are coded using the coding technique given in the work by (Abo-Zahhad and Rajoub, 2002). The compression algorithm is tested for 24 different records selected from the MIT-BIH Arrhythmia Database (MIT-BIH Arrhythmia Database). The experimental results show that the proposed method achieves high compression ratio with relatively low distortion and low computational complexity in comparison with other methods.
1. Introduction
The electrocardiogram (ECG) is an invaluable tool for diagnosis of heart diseases. ECG signals are usually sampled at 200–500 samples/s with 8–12 bits resolution. Considering long monitoring periods, compression is required to handle such vast amount of data. It can increase the capacity of databases where hundreds of thousands of ECG signals are stored for subsequent monitoring and evaluation. Other applications of ECG compression include transmission via telephone or mobile radio to an ECG center for further processing.
In recent years, many schemes have been suggested for compression of ECG data [1–11]. These compression techniques can be broadly classified into three groups: direct data handling techniques [1–3], transformation-based techniques [3–11], and parameterized model-based techniques [9]. Our approach belongs to the second group. In general, transform techniques involve expanding a signal into a weighted sum of basis functions. The coefficients of this sum are properly encoded and stored or transmitted instead of the original data. The best transform is the one which requires the least number of basis functions to represent the input signal for a given mean-square error (MSE). Transform techniques include several wavelet-based compression methods. The good time-frequency localization properties of wavelets make them especially suitable for data compression applications. Wavelets have been used in many data compression applications recently and have produced good results [4–8, 10].
Many of the previous works has employed the QRS-complex detection in compressing the ECG waveform. The detection of the QRS complex specifically, the detection of the peak of the QRS complex, or wave in an electrocardiogram (ECG) signal is a difficult problem since it has a time-varying morphology and is subject to physiological variations due to the patient and to corruption due to noise. For a tutorial on ECG signals, readers are referred to [12]. As noted in [13], most of the current QRS detectors can be divided into two stages: a preprocessor stage to emphasize the QRS complex and a decision stage to threshold the QRS-enhanced signal. Typically, the preprocessor stage consists of both linear and nonlinear filtering of the ECG. The ECG signal is first bandpass filtered to reduce noise and differentiated to emphasize the large slope of the wave; it is then squared to further exploit the high-frequency content of the QRS complex.
This paper is organized as follows. Section 2 presents the QRS significant features extraction and QRS-complex estimation. Section 3 is concerned with the preprocessing of the original ECG signal. Section 4 is a review of the two-dimensional discrete wavelet transform. Section 5 is an overview of coefficients grouping and thresholding and entropy principle. Section 6 shows the coding technique. Section 7 displays the experimental results of the compression algorithm on the selected records and a comparison with other coders in the literature. Finally, the paper is concluded in Section 8.
2. QRS-Complex Detection and Estimation
The ECG signal is composed from three main components, namely, the QRS complexes, -wave: and -wave as shown in Figure 1. Since the QRS complexes have a time-varying morphology, they are not always the strongest signal component in an ECG signal. Therefore, -waves or -waves with characteristics similar to that of the QRS complex, as well as spikes from high-frequency pacemakers, can compromise the detection of the QRS complex. In addition, there are many sources of noise in a clinical environment that can degrade the ECG signal. These include power line interference, muscle contraction noise, poor electrode contact, patient movement, and baseline wandering due to respiration [14]. Therefore, QRS detectors must be invariant to different noise sources and should be able to detect QRS complexes even when the morphology of the ECG signal is varying with respect to time. The significant features of the ECG signal are the , QRS, and -waves. In this work only the QRS waves are detected. Moreover, for best quality of QRS estimation and the duration between the - and - are detected and transmitted with the code stream. The QRS-complex estimation produce the typical QRS-complex waveform using the parameters extracted from the original ECG signal. The estimation algorithm is a MATLAB based estimator and is able to produce normal QRS waveform using a shifted and scaled versions of a triangular and sinusoidal wave forms.
3. Preprocessing of the Original ECG Signal
ECG signal preprocessing consists of three steps: (i) curettage, (ii) normalization and mean removal, (iii) and conversion the error matrix into 2D matrix.
3.1. The ECG Signal Curettage
The curettage process means to subtract the estimated QRS signal from the original ECG waveform upon the following equation: where is the error signal, , and is the length of the original ECG signal. The benefit of this step is to preserve the QRS signal which contains the most energy of the ECG signal.
3.2. Normalization and Mean Removal
Normalization of the error signal is done by dividing the error signal on the maximum value of the error signal , which is calculated using (2) as follows: Mean removal is performed by subtracting the mean value from each sample of the signal method. The mean of the ECG signal is calculated using the following equation: where is the number of samples in the error signal and is the samples of the error signal. Normalization and mean removal is done by the following equation: where is the normalized mean removed error signal. The main benefit of this step is to guarantee that the absolute values of all DWT coefficients are less than one.
3.3. Conversion into 2D Matrix
The dependencies in ECG signals can be broadly classified into two types: the dependencies among a single ECG cycle and the dependencies across ECG cycles. These dependencies are sometimes referred to as intrabeat and interbeat dependencies, respectively. An efficient compression scheme needs to exploit both dependencies to achieve maximum data compression. To compress the ECG data using the proposed coding scheme, the one-dimensional ECG signal has to be converted into a two-dimensional matrix. Thus, the first step in the presented algorithm is to separate each “period” of the ECG. This step is done with assistant of QRS complex detection which is described in Section 2. Each period is then stored as row in the 2D Matrix. It can be seen that the intrabeat dependencies are in the horizontal direction of the matrix and the interbeat dependencies are in the vertical direction. The created image is shown in Figure 2. Since each ECG period can have different lengths, using the constructed 2D ECG matrix using will have a different number of data points in each row. In the literature, there are many approaches tried to overcome this problem. In [16], the compression system using JPEG2000 normalizes each ECG period to the same length. Some approaches [20] tend to add anumber of zeros at the end of each heartbeat data sequence. Some drawbacks appear in approaches like the bits needed to be sent “or stored” to identify the length of each period [16, 20] which reduce the compression ratio (CR). Another drawback is the discontinuity in the matrix resulting from padding the end of each heartbeat with appropriate number of zeros [20]. Here, a novel technique is proposed for converting the 1D ECG signal into 2D array and warding the drawbacks of the previous approaches. Firstly, 32 periods are aligned with reference to the signal and all periods are cut at a certain length . The length is chosen at the minimum period length downward to a multiple of 32 “where 32 is the length of the constructed 2D ECG array.” For example, if the minimum period length is 267, then the length . The residual of all periods is assembled into one row then segmented into 32 by 32 matrices. After that, the matrices are added beside the first array “32 by matrix.” This approach exploits the benefit of the interbeat dependencies without additional bits.
4. Two-Dimensional Discrete Wavelet Transform (2D-DWT)
The continuous wavelet transform (CWT) is provided by (5), where is the signal to be analyzed. is the mother wavelet or the basis function. All the wavelet functions used in the transformation are derived from the mother wavelet through translation (shifting) and scaling (dilation or compression). Consider the following: The mother wavelet is used to generate all the basis functions, the translation parameter relates to the location of the wavelet function as it is shifted through the signal, and the scale parameter is defined as |1/frequency| and corresponds to frequency information. Notice that the wavelet transform merely performs the convolution operation of the signal and the basic function. The wavelet series is just a sampled version of CWT and its computation may consume significant amount of time and resources, depending on the resolution required.
The discrete wavelet transform (DWT), which is based on subband coding, is found to yield a fast computation of wavelet transform. It is easy to implement and reduces the computation time and required resources. In DWT, a time-scale representation of the digital signal is obtained using digital filtering techniques. The signal to be analyzed is passed through filters with different cut-off frequencies at different scales. Wavelets can be realized by iteration of filters with rescaling. The resolution of the signal, which is a measure of the amount of detail information in the signal, is determined by the filtering operations, and the scale is determined by upsampling and downsampling (subsampling) operations.
The DWT is computed by successive lowpass and high-pass filtering of the discrete time-domain signal , where is an integer, as shown in Figure 3. The low-pass filter is denoted by while the high-pass filter is denoted by . At each level, the high-pass filter produces detail information , while the low-pass filter associated with scaling function produces coarse approximations . Figure 4 shows the reconstruction of the original signal from the wavelet coefficients. Basically, the reconstruction is the reverse process of decomposition. The approximation and detail coefficients at every level are upsampled by two, passed through the low-pass and high-pass synthesis filters, and then added. This process is continued through the same number of levels as in the decomposition process to obtain the original signal. In most wavelet transform applications, it is required that the original signal be synthesized from the wavelet coefficients. So the analysis and synthesis filters have to be selected carefully to achieve perfect reconstruction. The 2D-DWT is simply the application of the 1D-DWT repeatedly to first horizontal data of the image, then the vertical data of the image as shown in Figure 5 and the inverse 2D-DWT is shown in Figure 6.
5. Coefficients Grouping and Thresholding
As a consequence of DWT decomposition, the resulted coefficients are classified into a few coefficients that contain the most energy of the signal “important coefficients” and a huge number of coefficients which contain less amount of energy of the signal “less important coefficients.” For efficient compression ratio and signal quality, this property have to be exploited correctly by dedicating more bits to represent the important coefficients and less bits to represent less important coefficients. Thresholding is a process that converts a range of coefficients values to a certain level based on threshold value . Normally, hard thresholding and soft thresholding techniques are used for such denoising process. Hard and soft thresholding with threshold are defined as follows.
The hard thresholding operator is defined as
The soft thresholding operator on the other hand is defined as Hard thresholding is “keep or kill” procedure and is more intuitively appealing and also it introduces artifacts in the recovered images. But soft thresholding is more efficient and it is used for the entire algorithm. In [21], thresholding methods are categorized into the following five groups based on the information the algorithm manipulates. (1)Histogram shape-based methods, where, for example, the peaks, valleys, and curvatures of the smoothed histogram are analyzed.(2)Entropy-based methods result in algorithms that use the entropy of the original and the transformed signals and the cross-entropy between the original and binarized signals.(3)Attribute-based methods search for a measure of similarity between the adjacent samples and periods.(4)The spatial methods use higher-order probability distribution and/or correlation between samples.(5)Local methods adapt the threshold value on each subband signal to the signal energy.
In this paper, entropy-based thresholding technique is adopted as a thresholding method. This class of algorithms exploits the entropy of the distribution of the coefficients levels. For illustration, consider a source that generates random symbols . For example, may be a digital image, and represents one of possible pixel levels. If denotes the probability of occurrence of symbol , then where is defined as the self-information for the symbol , that is, the information we get from receiving . If the base of the logarithm is two, then self-information is measured in bits. According to Shannon, the average information “Entropy” of a source is defined as For information theory, if the symbols are distinct, then the average number of bits needed to encode them is always bounded by their entropy. In practice, the entropy of a source in general is unknown, and estimates for the entropy depend on the probability model adapted for the structure of the symbols. For example, in a probability model, assume a source generates six distinct symbols (, , , , , and ) and that each symbol is generated with equal probability, that is, , for . The self-information of each symbol is . From the definition of the entropy [7], for this probability model . In other words, using this model a coding scheme for this sequence cannot be found that can code better than 2.585 bit per sample.
In a different model, assume that the probability of each symbol is and . Using this model, , , and . As a consequence of this discussion, it can be concluded that coefficients which occurred a few times yield to more information “more bits in coding” and those coefficients which occurred many times yield to less amount of information “less bits in coding.”
Hard thresholding of a certain group of subbands coefficients is done by eliminating all coefficients that are smaller than a certain threshold level . This process introduces distortion in a certain aspect in the reconstructed signal. To decrease the effect of thresholding, a soft thresholding technique is applied in this work. The coefficients which value with high self-information are thresholded to the nearest value that has low . The thresholding value in this proposal is not a constant value as normally used in the literature [22], the thresholding value is a variable depending on the self-information of the coefficients values.
To explore the effect of proposed thresholding technique, many thresholding rules “coefficients grouping” for the DWT coefficients had been suggested. In [23], the authors propose a list of six rules to group the wavelet coefficients. In this paper, two rules for “coefficients grouping” are listed to group the 2D wavelet coefficients. After wavelet transformation of the error signal, the wavelet coefficients are grouped to be thresholded upon a certain threshold value . The results coefficients of wavelet transform can be presented as in the following equation: where is the row vectors of the approximation coefficients, is the row vectors of the horizontal detail coefficients, is row vectors of the vertical detail coefficients, and is the row vectors of the diagonal detail coefficients. The first grouping scheme arranges the wavelet coefficients into three groups as in the following equations: The second grouping scheme divides the wavelet coefficients into five groups as in the following equations: Figures 7 and 8 show the first and second grouping schemes of 3-level 2D DWT decomposition.
The following steps describe the entropy thresholding algorithm.(1)Calculate the histogram and pdf distribution function of levels in each group of coefficients and the self-information of all levels as in (8) and the entropy.(2)Arrange the levels from the highest to the lowest based on the self-information of levels.(3)Threshold all coefficients which have levels with the highest self-information.(4)Repeat step 1, 2, and 3 until reaching the desired entropy.
6. The Coding Technique System
Figure 9 illustrates the QRS-complex estimation compression algorithm. The coding process is manipulated by dividing the coded stream into two parts: the header part and the thresholded coefficients part. The header part consists of two sections. The first section has 3 bits dedicated for pointing out the dimension of the 2D ECG signal (e.g., the 32 × 32 matrix is coded as “101”), 9 bits dedicated for storing the length of each period, 12 bits dedicated for storing the maximum value in the original signal, and 12 bits dedicated for storing the mean of the normalized signal. The second section has 36 bits to represent the , , and values and 12 bits to represent and duration. Figures 10(a) and 10(b) illustrate the coding stream that represents the header part.
(a)
(b)
(c)
(d)
The significant and insignificant coefficients are coded separately. The runs of significant coefficients are coded as follows.(i)One bit of value “1” identifies the run of significant coefficients.(ii)A sign bit to encode the sign of the significant coefficient.(iii)Eight bits to encode the value of the significant coefficient.
Figure 10(c) illustrates the coding stream that represents runs of significant coefficients.
The runs of insignificant coefficients are coded with variable-length coding VLC based on run length encoding as follows.(i)One bit of value “0” identifies the run of insignificant coefficients.(ii)Four bits to represent the number of bits needed to code the run length.(iii)Variable in length (from 1 to 16 bits) to encode the run length.
Figure 10(d) illustrates the coding stream that represents runs of insignificant coefficients.
The compression ratio CR, the percent RMS difference PRD, and peak signal-to-noise ratio PSNR will be used as a performance measure. The CR and PSNR are calculated respectively by (13) as follows: where and represent the original and the reconstructed signal, respectively. Three different PRD formulas can be found in the open literature. The first formula is described by (14) as follows: The second formula is expressed by (15) as follows: Finally, the last formula is expressed by The RMS of a set of values is the square root of the arithmetic mean (average) of the squares of the original values as
7. Experimental Results
The MIT-BIH Arrhythmia Database [24] has been used to evaluate the performance of the proposed compression algorithm. The ECG signals used here were sampled at 360 Hz and each sample is represented by 11 bits/sample. We used two datasets formed by taking certain records from the MIT-BIH Arrhythmia Database. These datasets were chosen because they were used in earlier studies and allow us to compare the performance of the proposed method with other coders in the literature. The first dataset consists of 10 min of data (each) from record numbers 100, 101, 102, 103, 107, 109, 111, 115, 117, 118, and 119. The second dataset consists of 1 min of data (each) from record numbers 104, 107, 111, 112, 115, 116, 117, 118, 119, 201, 207, 208, 209, 212, 213, 214, 228, 231, and 232.
7.1. Experiment 1
This experiment conclude testing the proposed algorithm on records 100 and 220 in order to explore and reveal the effectiveness of the proposed method on the clinical diagnosis information of the reconstructed ECG signal. Figure 11 shows the compression performance of the proposed method of the selected two records of the MIT-BIH Arrhythmia Database. The evaluation of the figure shows that the reconstructed ECG signal had preserved the most of the clinical diagnosis information of the original ECG signal.
(a)
(b)
7.2. Experiment 2
In this experiment the proposed algorithm was applied to 10 ECG records selected randomly from the MIT-BIH arrhythmia database. These records are 100, 101, 102, 103, 107, 109, 111, 115, 117, and 119. Figure 12 explore the results of this experiment. The results demonstrate that the compression results versus the percentage RMS difference for the tested records are converging to each other, which mean that the proposed compression method is opportune for all ECG signals.
7.3. Experiment 3
This experiment is a comparison of the proposed compression algorithm with previous compression methods based on 2D ECG signal in the literature as follows. In [16], Bilgin et al. applied the image coding standard—JPEG to ECG signal compression. In [20], Tai et al. adopted a modified SPIHT method to compress ECG signals. In [18], Lee and Buckley proposed an ECG compression method based on the 2D DCT. In [25], Xingyuan and Juan have applied wavelet-based hybrid ECG compression technique on ECG signals.
8. Conclusion
In this paper, new wavelet based ECG compression technique is proposed, associated with acceptable PRD and PSNR values. It is based on converting 1D ECG signal with irregular periods into 2D matrix of size . The length is chosen at the minimum period length downward to a multiple of 32. The residual of all periods is assembled into one row then segmented into 32 by 32 matrices. This approach exploits the benefit of the interbeat dependencies without additional bits. Numerical results show that the 2D ECG compression is better than 1D ECG compression methods as shown in Table 1. Table 1 illustrates also the detailed comparison between the proposed algorithm and the last listed methods for 117, 119, 205, and 220 ECG records. It is clear that the proposed algorithm has the highest CR with acceptable PRD. Also, it is clear that the coding and decoding stages in the proposed algorithm are fast and easy to implement. Moreover, the proposed algorithm has a good performance for different ECG signals considered from clinical point of view and all the clinical information is preserved after reconstruction. Finally, it should be noted that a further improvement in results may be achieved with sophisticated implementation of compression system in any simulation program as LABVIEW.