A singular spectrum analysis-based model-free electrocardiogram denoising technique

Journal Pre-proof

A singular spectrum analysis-based model-free electrocardiogram denoising technique Sourav Kumar Mukhopadhyay , Sridhar Krishnan PII: DOI: Reference:

S0169-2607(19)30980-0 https://doi.org/10.1016/j.cmpb.2019.105304 COMM 105304

To appear in:

Computer Methods and Programs in Biomedicine

Received date: Revised date: Accepted date:

20 June 2019 19 December 2019 25 December 2019

Please cite this article as: Sourav Kumar Mukhopadhyay , Sridhar Krishnan , A singular spectrum analysis-based model-free electrocardiogram denoising technique, Computer Methods and Programs in Biomedicine (2020), doi: https://doi.org/10.1016/j.cmpb.2019.105304

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier B.V.

Highlights    

It is a data-driven, model-independent and high performance ECG denoising technique. The proposed ECG denoising technique outperforms state-of-the-art ECG denoising methods. The mean opinion score of the denoised ECG signal falls under the category „very good‟ as per the gold standard subjective measure. The proposed technique could be adapted for denoising other biomedical signals exhibiting periodic or quasi-periodic nature.

1

A singular spectrum analysis-based model-free electrocardiogram denoising technique

Sourav Kumar Mukhopadhyay Department of Electrical, Computer, and Biomedical Engineering, Ryerson University, 350 Victoria Street, Toronto, Ontario, Canada, M5B 2K3 Email: [email protected], Phone: +1 416-979-5000, Extn.: 2048

Sridhar Krishnan Department of Electrical, Computer, and Biomedical Engineering, Ryerson University, 350 Victoria Street, Toronto, Ontario, Canada, M5B 2K3 Email: [email protected], Phone: +1 416-979-5000, Extn.: 4931

2

Abstract Background and objective: An efficient and robust electrocardiogram (ECG) denoising technique caters three-fold benefits in the subsequent processing steps: first, it helps improving the accuracy of extracted features. Second, the improved accuracy in the extracted features enhances the performance as well as the reliability of computerised cardiovascular-disease diagnosis systems, and third, it also makes the interpretation task easier for the clinicians. Albeit a number of ECG denoising techniques are proposed in the literature, but most of these techniques suffer from one or more of the following drawbacks: i) model or function dependency, ii) sampling-rate dependency, or iii) high time-complexity. Methods: This paper presents a singular spectrum analysis (SSA)-based ECG denoising technique addressing most of these afore-mentioned shortcomings. First, a trajectory matrix of dimension K × L is formed using the original one-dimensional ECG signal of length N. In SSA operation the parameter L, which is denoted as the window-length, plays a very important role and is related to the sampling frequency of the signal. In this research the value of L is calculated dynamically based on the morphological property of the ECG signal. Then, the matrix is decomposed using singular value decomposition technique, and the principal components (PC) of the original ECG signal are computed. Next, the reconstructed components (RC) are calculated from the PCs, and all the RCs are filtered through Butterworth bandpass and notch filters. An optimum number of filtered RCs are retained based on their significance. Finally, these retained RCs are summed up to obtain the denoised ECG signal. Results: Evaluation result shows that the proposed technique outperforms state-of-the-art ECG denoising methods; in particular, the mean opinion score of the denoised signal falls under the category „very good‟ as per the gold standard subjective measure. Conclusions: Both the quantitative and qualitative distortion measure metrics show that the proposed ECG denoising technique is robust enough to filter various noises present in the signal without jeopardizing the clinical content. The proposed technique could be adapted for denoising other

biomedical

signals

exhibiting

periodic

or

quasi-periodic

nature

such

as

photoplethysmogram and esophageal pressure signal.

Keywords: ECG distortion, dynamic window-length, ECG denoising, mean opinion score, singular spectrum analysis, data-driven denoising. 3

1. Introduction Electrocardiogram (ECG), which is considered as the prime human physiological signal [1], is the graphical representation of heart‟s electrical activity recorded by means of electrodes from different standardized locations on the body surface. However, ECG signals get routinely contaminated during acquisition by the influence of several high and low frequency noises, which not only exacerbates the task of interpretation for doctors, but also degrades the overall performance of computerized ECG processing systems. Respiration, patients‟ movement, powerline interference, electromyogram (EMG), electro-surgical and instrumentation noises are the common types of noises, which contaminate the ECG signal. Researchers around the globe are continuing to develop robust and high performance noise elimination techniques from ECG signal for over 30 years [2]. In recent times, wavelet transform (WT) [3]-[6], [37], empirical mode decomposition (EMD) and ensemble empirical mode decomposition (EEMD) [7]-[11]-based ECG denoising techniques have drawn significant attention. In [3], WT is used in combination with β-hill climbing metaheuristic algorithm for denoising ECG signal, where the original ECG signal is corrupted by mixing with synthetic white Gaussian noise (WGN), and then the optimum filtering parameters such as mother wavelet function, number of decomposition levels, thresholding method are selected for denoising. In [3] the authors have claimed that even for the same ECG signal at a constant sampling rate: (i) the selection of mother wavelet function varies with varying WGN characteristics, and also (ii) the number of decomposition level varies with varying WGN characteristics. But, the reason behind such an interesting fact is left unclarified. Jain et al. have proposed an EMD-based ECG denoising method in [7]. In [7], first, the signal is decomposed into intrinsic mode functions (IMF), and a few high-order IMFs are empirically considered to be corrupted with low frequency noises. Then, these considered IMFs are scaled, and denoised using a Savitzky Golay filter for the removal of low frequency noises. Rest of the IMFs are filtered using a Riegmann Liouvelle fractional integral filter for high frequency noise removal. However, the scaling factors are chosen based on the performance of the algorithm only on one ECG segment of 800 samples (2.22 seconds of duration). In [8], Rakshit et al. have proposed an EMD, WT and adaptive switching mean filter (ASMF)-based ECG denoising technique, where the noisy signal is decomposed into six IMFs and one residue signal. Then, each of the first three IMFs are again decomposed into two levels using WT, and the wavelet coefficients are thresholded to discard 4

the high frequency noises. Next, these three filtered IMFs are added with other three unfiltered ones to get a partially-filtered ECG signal. Then, the R-peaks are detected from this partially filtered ECG. Next, the partially-filtered signal is again processed using ASMF. Since the ASMF method highly distorts the QRS complexes regions of the ECG, the previously detected R-peak information are used to replace the distorted QRS complexes with their corresponding undistorted ones obtained before applying ASMF. Use of extended Kalman filter (EKF) is proposed in [13] and [14] for denoising ECG. In [14], first, the noisy ECG signal is preprocessed using variational mode decomposition technique along with WT-based method. Then, the Rpeaks are detected, and these detected ECG-beats are classified into normal and abnormal with the help of k-nearest neighborhood algorithm. Finally, the signal is denoised using marginalized particle-extended Kalman filter. However, the performance of the technique strongly depends on the accuracy of the classification algorithm. Sparse decomposition and dictionary learning-based ECG denoising techniques are proposed in [15] and [16]. In [15], a moving average filter is used to extract the baseline wander noise, and in total eight different threshold values are used. However, no proper explanation is given for using such threshold values. Besides these above mentioned techniques, a number of efficient ECG denoising methods based on spectro-temporal filtering [12], adaptive Fourier decomposition [17], fractional zero-phase filtering [18, 33], non-local means [19], Bayesian framework [20] and adaptive filtering [21] are also proposed in the literature in recent times. Singular spectrum analysis (SSA) is a model-free and data-driven time-series-decomposition method, which decomposes a time series into three components: trend, seasonal components and noise [22]. These components are generated based on the covariance property of the original time series data, and they are independent and additive in nature [34, 35]. SSA consists of two complementary stages: decomposition and reconstruction. The decomposition stage is divided into two sub-stages: embedding and singular value decomposition (SVD), and the reconstruction stage is also sub-divided into Eigen-triplet grouping and diagonal averaging. Use of SSA-based methods are proven to be very efficient in a wide range of applications including climatology [23], electricity consumption forecasting [24], biomedical image processing [25], [26], gait parameter estimation [27]. However, the potential of denoising biosignals; in particular, ECG signal, using SSA-based method is not explored to a good extent to date. Yang et al. have proposed an SSA-based ECG denoising technique in [28], where a guideline, however empirical, 5

is provided for selecting the window-length L, which is an important parameter to set in SSA operation. In [28], first, synthetic ECG signals are taken and mixed with power-line interference noise. Then, the noisy signals are processed using SSA at different values of L, and the denoised ECGs are shown. In [28], the SSA-based denoising technique is tested on ECG signals recorded at 250 Hz and 1250 Hz sampling rates with four different values of . The values of , which are used at 250 Hz sampling rate are 10, 13, 15 and 18, and the values of , which are used at 1250 Hz sampling rate are 50, 65, 75 and 90. The assumption, based on which the values of L are chosen in [28], is that

. However, neither the value of

is optimized for

different sampling rates, nor the noise removal efficiency at different values of

and (or) at

different noise levels are studied. Motivation behind the proposed research work is to develop a high performance and model-free ECG denoising technique. The proposed technique also addresses most of these above mentioned hurdles and shortcomings of other ECG denoising techniques. The main challenges, which are addressed in this research work are (1) denoising WGN and random noisecontaminated ECG signals; which is the hardest problem faced by the researchers as the bandwidth of the ECG signal overlaps with the bandwidths of these two types of noises, (2) efficient enough to expel out other types of noises such as baseline-wander noise, instrumentation noise, power-line interference noise, (3) designing an ECG denoising technique that is adaptive to the sampling rate, and (4) designing an ECG denoising technique that is datadriven, and whose performance does not depend on any pre-assumed mathematical function/model unlike the wavelet transform-based ones. In this proposed research work, first, the original ECG signal of length N is arranged in form of a trajectory matrix of dimension K × L, where the window-length L is calculated dynamically. Next, the trajectory matrix is decomposed using SVD method. Then, the principal components (

) of the original ECG

signal are computed by a linear combination of the trajectory matrix and the matrix of eigenvectors. Next, reconstructed components ( the

) are calculated from these

. All

are filtered through Butterworth bandpass and notch filters. Thereafter, the dynamic

weights of the filtered

are computed, and an optimum number of filtered

based on their significance. Finally, these retained

are retained

are summed up to obtain the denoised

ECG signal.

6

The remainder of this paper is organized as follows. The proposed SSA-based ECG denoising technique is presented in Section II. Denoising performance of the proposed technique is presented in Section III. Section IV compares the denoising performances of the proposed technique with other state-of-art techniques. Finally, discussions and conclusions are drawn in Section V.

2. SSA-based ECG denoising technique The original ECG signal (denoted as ) of length N is mapped into L-lagged vectors of size K, and a trajectory matrix, which is denoted as rows and columns of

, is formed with these lagged vectors. Both the

are the subseries of . The first column of the matrix

the entire original-ECG-signal

(i.e., all the QRS-complexes, P-waves, T-waves, ST-segments

and PR-segments); the second column of the matrix

contains the original ECG signal

shifted by 1 samples; the third column contains the original ECG signal and likewise, the

contains

column contains the original ECG signal

shifted by 2 samples,

shifted by L-1 number of

samples.

[

where,

]

.

From Equation 1 it is evident that the elements of all the anti-diagonals of the trajectory matrix are the same, i.e.,

. Therefore, the trajectory matrix

is a Henkel matrix.

Selection of the proper value of the window-length L is of great importance in SSA operation. The efficiency of SSA in separating different signal-components from a composite-signal depends on the value of L. An SSA-based method is able to separate two signals from a composite signal if their frequency-difference is ≥

. On the other hand, SSA fails to

separate two signals from a composite signal if (1) the amplitudes of both the signals are equal, or (2) their frequency-difference is less than

. A detailed mathematical-analysis and

theoretical-explanation of SSA method can be found in [40]. However, The window-length L should be chosen in such a way that all the ECG-beats (i.e., the QRS-complexes) which are

7

present in , are also to be there in all the lagged vectors, i.e., in all the columns of

.A

small value of L may not help denoising the signal at all, and the use of a large L value (i) is redundant, (ii) increases the computational burden, and most importantly (iii) may not guarantee that all the ECG-beats are present in all the columns of

, which leads to a significant loss of

clinical information. In this research work the window-length L is calculated based on the inherent periodicity and morphological property of the ECG. Since, the maximum possible heart rate is 240 beats per minute [29], the optimum value of the window-length L can be calculated as below.

Equation 2 guarantees that (i) all the ECG-beats which are present in , are also present in all the columns of

, and (ii) the ECG-beats present in the adjacent columns of the

matrix

would never overlap with their preceding QRS complexes. Let us take an example. It is well known that the clinical bandwidth of an ECG signal lies between 0.5 to 100 Hz [1]. Now, if an ECG signal is recorded at a sampling rate of 1 kHz, then, following Equation 2, the value of L = 250. Therefore, (1) any high frequency noise above (100 Hz + Hz) 100.004 kHz, and (2) any low frequency noise below (0.5 Hz - Hz) 0.496 Hz could be expelled out from the ECG using an SSA-based technique. Structures of the L-lagged vectors containing the subseries of Equation 2 it can be noted that columns of the matrix number of rows of the matrix the dimension of the matrix

are shown in Figure 1. From

of the ECG signal. Therefore, the number of

depends on the sampling rate of the original ECG signal, and the depends on the length of the original ECG signal. Hence, is dynamic.

8

Figure 1: Structures of the L-lagged vectors containing the subseries of . In the next step,

is decomposed using SVD to represent it as a sum of rank-one

biorthogonal elementary matrices. The SVD technique factorizes a matrix into the product of another three matrices: an orthogonal matrix ( ), a diagonal matrix ( ) and the transpose of an orthogonal matrix (V). The eigenvalues of the matrix S

are denoted as

, which are in the decreasing order of magnitude, i.e., corresponding eigenvectors of (

) is call the

expressed as 40]. If

. The

-eigentriple of the SVD. The energy-contribution of the

∑

collection of -eigentriple is

, which is also known the singular-spectrum of the time series data [39, √

, then it is possible to write the trajectory matrix as:

∑

where

are denoted as

, and the

∑√

,

√

right eigenvector [27]. In SVD operation

is the ith left eigenvector, and

is the ith

is denoted as the right singular

matrix, which contains the information about the most important axis of the data. The first right singular vector

reveals the direction having the most variance, and

reveals the direction

having the least variance. The variance is precisely defined by the singular values. The strong inter-sample and inter-heartbeat correlation of the ECG signal enhances the covariance among the eigenvector of the matrix . Since the SVD operation arranges the singular values in decreasing order of magnitude, the small singular values barely carry any signal-information, and 9

these small singular values could be discarded [36]. The original signal

can be considered as

the composition of the ECG signal and the noise. If the indices

of the

eigenvalues, which hold the most of the ECG-signal-information are known, then the matrix corresponding to the ECG-signal-only can be written as follows [28, 38]. ̅̅̅̅̅ ̅̅̅̅̅ ̅̅̅̅̅ ̅̅̅̅̅

∑

[̅̅̅̅̅

̅̅̅̅̅ ̅̅̅̅̅̅

̅̅̅̅̅

̅̅̅̅̅̅]

However, the indices, which hold the most of the ECG-signal-information are unknown. Next, the

are computed by a linear combination of the trajectory matrix

and the matrix of

eigenvectors

The columns of the matrix dimension

are the principal components of . Then, a matrix of

is created by an inverse projection of the

Finally, the reconstructed components (

and the matrix of eigenvectors .

s) of the original ECG signal are computed by

averaging along the anti-diagonals of this matrix [28]. ∑ ̅̅̅̅̅ where,

,

. The original ECG signal

is the number of elements in

can be obtained by summing up all the

.

∑ Now, all the

s are filtered using a zero-phase 4th order Butterworth bandpass filter having

lower and upper cut-off frequencies 0.5 Hz and 100 Hz [1], respectively, and a notch filter to remove the power-line interference noise. Filtered example of first and last two

s and

s are denoted as

s. Figure 2A shows an

, and Figure 2B shows their corresponding frequency

spectrums.

10

Figure 2A: An example of first and last two

s and

.

11

Figure 2B: (A) Frequency spectrum of (C) Frequency spectrum of and

and , (B) Frequency spectrum of , and (D) Frequency spectrum of

From Figure 2A it can be noted that the amplitudes of the high-order

and and

s

, . get

noticeably suppressed compared to that of the low-order counterparts afterward the bandpass and notch filtering operations, which strongly suggests that the high-order

s are noise dominant.

Figure 2B corroborates this notion. Now, to discard the less significant high-order

s, a

dynamic weight calculation-based method is used. Dynamic weights of all the

s are

calculated using the following equation. ∑ ∑ where

∑

is the dynamic weight of the denoted the

coefficients of the

,

denotes the total number of

. Since the amplitude-bands of

difference between the maximum and minimum amplitude of each

, and i.e., the

decreases progressively,

so does their dynamic weights. An example of the distribution of dynamic weights among are shown in Figure 3.

Figure 3. Distribution of dynamic weights among Next, starting from the first

(

), an optimum number of

. are summed up until

their total weight accumulates a predefined threshold value (γ). This summed up considered as the denoised ECG, which is denoted as

. The

is

-summing algorithm is

given below.

12

where

is the optimum number

.

3. Experimental results Performance of the proposed SSA-based ECG denoising technique is tested through both quantitative and qualitative distortion measures. Digitized ECG signals are collected from five different databases [30] and are used as the evaluation testbeds of the implemented ECG denoising technique. Details of these collected ECG data are shown in Table 1. Table 1: Details of collected ECG data. Database

Sampling rate in Hz

European ST-T Database (edb)

250

MIT-BIH arrhythmia database (mitdb) MIT-BIH noise stress test database (nstdb) Motion Artifact Contaminated ECG Database (macecgdb)

360

PTB diagnosis ECG database (ptbdb)

1k

360 500

No. of lead taken into consideration 2 (ML III and V4) 1 (ML II) 2 (ML II and V1) 4 (ECG-1, 2, 3, 4) 12 (lead I-III, aVR, aVL, aVF, V1-V6)

No. of samples in each lead 10 105 21.6

103

No. of ECG files taken 89 48 12

4 60

103 103

27 519

In total, the technique is tested on 695 ECG data-files. Out of 695 files, 519 files are taken from PTBDB. Among these 519 files, 17 files contain bundle branch block ECG, 16 files contain cardiomyopathy ECG, 15 files contain dysrhythmia ECG, 3 files contain heart failure ECG, 7

13

files contain hypertrophy ECG, 367 files contain myocardial infraction ECG, 4 files contain myocarditis ECG, 80 files contain normal ECG, 1 files contains palpitation ECG, 2 files contain stable angina ECG, 1 files contains unstable angina ECG, and 6 files contain valvular heart disease ECG. All the ECG recordings of nstdb are highly corrupted with baseline wander, electromyogram and electrode motion artifact noises, which are recorded from physically active volunteers. ECG recordings of macecgdb are highly corrupted with motion artifact noise. ECG signals of these above mentioned databases are recorded using wired ECG acquisition modules. No such database is available online to date, where wireless ECG recordings are freely available. Therefore, the performance of the technique is tested only on wired ECG recordings. All these digitized ECG signals collected from different databases, already contain low and high frequency noises, and hence, a little preprocessing is required. Therefore in the testing phase of the proposed ECG denoising technique, first, all these collected ECG signals are preprocessed through a zero-phase Butterworth bandpass filter to remove the low frequency baseline wander as well as high frequency noises above 100 Hz. It is important to note that the preprocessing is only performed at the testing phase of the technique, and is not required when denoising realtime ECG signals. The preprocessed signal is denoted as synthetic noises are added with

. Next, four different types of

. Details of these synthetic noises are given below.

(1) 50 Hz powerline interference noise

.

where where,

,

(2) White Gaussian noise

of power

at different standard deviation values (σ=0.1

to 5.0) and zero mean (μ = 0). (3) Random noise

in the frequency-band 20 Hz - 500 Hz with different standard deviation

values (σ=0.1 to 5.0) and zero mean (μ = 0). (4) 0.2 Hz baseline wander noise

The noisy ECG is denoted as

. Now, the

signal is denoised using the proposed SSA-based denoising technique, which is described in Section 2. The denoised ECG signal is denoted as

.

14

In the performance evaluation phase, the metrics: input signal-to-noise ratio ( signal-to-noise ratio (

), improvement in signal-to-noise ratio (

root-mean-square difference

), output ), and percent-

are used as quantitative non-diagnostic distortion measures,

and wavelet energy-based diagnostic distortion ( diagnostic distortion measure.

,

,

) [31] is used as the quantitative , and

are calculated using the

following equations. ∑

(8)

∑

where

, and ∑ ∑

(9)

where ∑

(10)

∑

∑ √

∑

The higher the value of performance. However, except

, and the lower the value of and

and

, the better the

the other metrics are not standardized yet. As

per the criteria given in [31], a denoised/reconstructed ECG signal is considered as: (1) „excellent‟, if the WEDD value is <4.517%, (2) „very good‟, if it is within the range 4.517%– 6.914%, (3) „good‟ if it falls in the range 6.914%-11.125%, (4) „not bad‟ if it falls in the range 11.125%-13.56%, and (5) „bad‟ beyond that. Also, according to the globally considered standard [32], quality of a denoised/reconstructed ECG signal is considered to be (1) „very good‟ if the PRD value lies in between 0 to 2%, (2) „good‟ if the PRD value lies in between 2% to 9%, and (3) „not good‟ if the PRD value is >9% . The ratio of denoted as

. The closer the value of

to

is calculated and

to unity, the better the performance.

All the 48 ECG signals (lead ML II) of MIT-BIH arrhythmia database are mixed with PLI noise (generated using Equation 6), BWN (generated using Equation 7), and both 30dB WGN and RN at σ =1.5 and μ=0. Then, these noisy signals are denoised using the proposed technique, and the performance is shown in Table 2.

15

Figures 4 to 7 show the performance of the proposed ECG denoising technique at different sampling rates and at

. Figure 4 shows the

different types of synthetic noises which are added to

,

,

as well as four

. Synthetic PLI noise and BWN,

which are shown in Figure 4 are generated using Equations (6) and (7), respectively, and 30dB WGN and RN are generated at σ=1.5 and μ=0. Exactly the same amount of noises are also added with other signals, which are shown in Figures 5 to 7, and hence these synthetic noises are not shown further in these figures. Table 2: Performance of the proposed ECG denoising technique on the first lead (ML II) of MITBIH arrhythmia database at γ=99.999%. PLI noise is generated using Equation 6, BWN is generated using Equation 7, and both WGN (30dB) and RN have σ=1.5 and μ=0. File ID δ 100 101 102 103 104 105 106 107 108 109 111 112 113 114 115 116 117 118 119 121 122 123 124 200 201 202 203 205 207 208 209 210 212 213

-21.93 -19.70 -18.48 -16.75 -17.56 -19.14 -16.29 -12.20 -21.35 -15.97 -16.69 -14.17 -12.16 -24.12 -15.55 -16.86 -14.67 -15.10 -13.47 -16.03 -11.81 -16.78 -16.20 -16.10 -15.25 -15.29 -16.88 -19.52 -16.64 -14.31 -16.81 -17.10 -15.12 -12.54

-21.95 -19.73 -18.51 -16.75 -17.6 -19.22 -16.35 -12.21 -21.42 -15.97 -16.74 -14.17 -12.18 -24.22 -15.55 -16.89 -14.71 -15.10 -13.54 -16.08 -11.81 -16.79 -16.21 -16.13 -15.28 -15.31 -16.92 -19.52 -16.68 -14.36 -16.83 -17.12 -15.14 -12.55

0.999089 0.998479 0.998379 1 0.997727 0.995838 0.99633 0.999181 0.996732 1 0.997013 1 0.998358 0.995871 1 0.998224 0.997281 1 0.99483 0.996891 1 0.999404 0.999383 0.99814 0.998037 0.998694 0.997636 1 0.997602 0.996518 0.998812 0.998832 0.998679 0.999203

45.38 43.29 39.60 44.56 39.14 42.89 39.83 44.44 42.03 48.44 38.71 40.80 39.55 42.30 45.89 43.23 42.10 45.80 38.57 39.76 42.72 47.00 48.38 43.48 39.01 41.16 41.95 47.37 41.07 38.56 41.02 42.02 38.98 41.57

6.72 6.62 8.80 4.07 8.34 6.50 6.65 2.44 9.25 2.38 7.93 4.66 4.27 12.34 3.04 4.80 4.25 2.92 5.56 6.51 2.85 3.08 2.46 4.28 6.48 5.09 5.58 4.05 6.00 6.13 6.16 5.67 6.41 3.54

3.65 3.17 3.74 2.01 3.43 2.84 2.72 0.62 4.38 0.77 3.33 2.51 1.93 6.91 1.95 1.71 2.88 1.36 2.25 3.54 1.24 2.03 1.35 1.77 3.35 2.10 1.82 2.86 2.54 2.14 3.16 2.26 2.64 1.12 16

214 215 217 219 220 221 222 223 228 230 231 232 233 234 Average

-15.35 -16.84 -13.27 -13.63 -16.57 -15.73 -18.28 -17.23 -17.52 -17.69 -16.23 -17.51 -13.77 -15.01 -16.32

-15.35 -16.85 -13.30 -13.64 -16.57 -15.75 -18.43 -17.23 -17.65 -17.69 -16.26 -17.52 -13.78 -15.01 -16.35

1 0.999407 0.997744 0.999267 1 0.998730 0.991861 1 0.992635 1 0.998155 0.999429 0.999274 1 0.998285

45.14 42.09 41.17 44.03 45.37 41.64 35.45 49.26 38.90 48.41 43.67 42.07 46.33 45.45 42.70

3.24 5.47 4.03 3.02 3.63 5.06 13.86 2.50 8.54 2.91 4.24 5.92 2.35 3.00 5.28

1.06 2.29 1.33 1.28 2.07 2.15 7.78 1.31 4.01 1.64 2.74 3.85 1.03 1.71 2.51

From Table 2, it is noted that the PRD and WEDD values are maximum for the ECG signal contained in the file #222. Figure 8 shows the original ECG signal contained in file #222 as well as its denoised version. From this figure it can be seen that the original ECG signal itself is very noisy, and the denoising technique extract the added synthetic noise and also the other noises which were present in the signal beforehand. This is the reason for getting such high values of and value of

. Moreover, the value of

is minimum for this ECG signal; suggesting a high

, which indicates that the power of the denoised signal is higher than that of the

power of the extracted noise. Moreover, it can also be seen from Table 2, that in most of the cases, the value of

is found to be very close to unity, which suggests that the denoising

technique is robust enough to extract most of the synthetically added noises.

17

Figure 4: (a) : MIT-BIH 100, sampling rate 360 Hz, (b) PLI noise, (c) 30 dB WGN, μ=0, σ=1.5, (d) RN, μ=0, σ=1.5, (e) BWN, (f) , (g) and on top of each other: , , , , and .

Figure 5: (a) : edb-e0103, lead-V4, sampling rate 250 Hz, (b) on top of each other: , and .

, (c) ,

and ,

18

Figure 6: (a) and

: macecgdb test01_00s, signal ECG-3, sampling rate 500 Hz, (b) on top of each other: , , , and .

,(c) ,

Figure 7: (a) (b) , (c) ,

: PTBDB s0431_re (myocardial infraction), lead-I, sampling rate 1 kHz, and on top of each other: , , and .

19

Figure 8: Original ECG signal (file #222, MIT-BIH arrhythmia database) and its denoised version plotted on top of each other. Now, the variation in average , , and values at different levels of input noises (PLI noise is generated using Equation 6, BWN noise is generated using Equation 7, 30dB WGN and RN with μ=0 and σ=0.1 to 5.0) on MIT-BIH arrhythmia database are shown in Figure 9. Form this figure it can be noted that the variation patterns of the average values of

, and

with σ take the shape of polynomials of order 5.

Figure 9: (a) variation in the average values of with σ ( =0.0033σ50.05σ4+0.31σ3+0.87σ2-0.36σ+44), (b) variation of with σ ( =-0.0017σ5 4 3 2 +0.023σ -0.14σ +0.48σ +0.036σ+4.5), (c) variation of with σ ( = -0.0003σ5+0.0043σ4-0.026σ3+0.085σ2+0.83σ+1.1). The σ of WGN and RN is varied together. 20

Now, the variation in the denoising performance of the proposed technique on mitdb for different values of γ at a fixed noise level is shown in Table 3. Table 3: Variation in the denoising performance of the proposed technique for different values of γ at a fixed noise level on MIT-BIH arrhythmia database. PLI noise is generated using Equation 6, BWN noise is generated using Equation 7, 30dB WGN and RN with μ=0 and σ = 2.5. γ in % 90.00 26.95 29.57 16.24 92.00 27.90 26.49 14.03 94.00 29.15 22.94 11.60 96.00 30.71 19.17 9.26 98.00 33.57 13.85 6.36 99.00 35.78 10.80 4.98 99.99 40.84 6.44 3.50 99.995 40.91 6.40 3.50 99.999 41.24 6.19 3.48 100.00 41.24 6.19 3.48

Figure 10 shows the denoising performance of the proposed technique for different values of γ on a particular ECG signal (file #200) taken from mitdb. From Table 3 it can be noted that both the values of

and

values of

,

decreases with increasing the value of , and vice-versa. Also, the changes rapidly up to a γ value of 99.99%, and then get

and

saturated at a γ value of 99.999%. As per the

criteria, the denoised ECG signal falls under

the category „not good‟ below γ=99.999%, and as per the

criteria the denoised ECG

signal falls under the category „not bad‟ below γ=94%. It important to note that both the values of

and

are identical at γ=99.999% and 100%, and the reason is as follows. It is seen

in Figure 2A that the amplitudes of the high-order low (of the order of

s

are significantly

) compared to that of the low-order

s

. It

suggests that most of the signal-information is contained in the low-order addition of those high-order

s

does not help enhancing the quality of

the reconstructed ECG signal farther, and therefore, the values of change noticeably above a γ value of 99.999%. However, the value of number of

s. Therefore,

and

do not

, i.e., the optimum

required for γ=99.999% is less than that of for γ=100%, i.e., . In our research it is observed that

γ at a fixed value of σ, and

sampling frequency of the signal when both σ and γ are constant. Figure 11 shows the variation of

with σ for different values of γ on two different ECG signals: (i) file#100 taken from 21

mitdb (sampling rate 360 Hz) and (ii) file#s0273lre (sampling rate 1 kHz) taken from ptbdb. From Figure 11a, it can be noted that the value of double than that of the value of last-43

) for γ=100% is almost

(

) for γ=99.999% at σ=0.1. It suggests that the

(

s

) do not contain any such signal-information, which

might help improving the quality of the reconstructed ECG signal. Similarly, Figure 11b shows that the value of (

) at γ=100% is almost five times than that of the value of

( )

for

γ=99.999%

s

at

σ=0.1.

It

suggests

that

the

last-197

) do not contain any such signal-information, which could

help enhancing the quality of the reconstructed ECG signal. It is also to be noted from both Figure 11a and 11b, that (i) the difference between the values of

at γ=99.999% and 100%,

respectively, reduces with increasing the value of σ, and (ii) the values of

do not change

significantly below a γ value of 99%. Therefore, in this research, 99.999% is found as the optimum ( value of

value

of

γ

for

two

reasons:

(1)

the

performance

measure

metrics

) provide identical results both at γ = 99.999% and 100%, and (2) the is less at γ = 99.999% compared to that of at 100%; particularly at low value σ.

Henceforth, a γ value of 99.999% is used in preparing all the tables and figures in the rest of this paper. One can also chose a γ value of 100%, but it would not help improving the quality of the reconstructed signal, and also it would increase the processing time. Now, Table 4 shows the performance of the proposed denoising technique on ECG signals collected from different databases at different noise levels.

22

Figure 10: (a) (PLI noise is generated using Equation 6, BWN noise is generated using Equation 7, 30dB WGN and RN with μ=0 and σ = 2.5), (b) at γ=90%, (c) at γ=92%, (d) at γ=94%, (e) at γ=96%, (f) at γ=98%, (g) at γ=99%, (h) at γ=99.999%, (i) at γ=100%. From (b) to (i) the and signals are plotted on top of each other.

Figure 11: (a) variation of with σ (0.1 to 5) at different values of γ on the ECG signal #100 of mitdb (sampling rate 360 Hz, L=90), (b) variation of with σ (0.1 to 5) for different values of γ on the ECG signal #s0273lre of ptbdb (sampling rate 1000 Hz, L=250). 23

Striking information, which can be extracted from Table 4 are: 1) in most of the cases the denoising performance of the technique enhances with increasing the sampling rate, 2) in most of the cases both

and

values decrease with increasing sampling rate at a fixed value of

σ, and 3) in most of the cases the rate of change of both

and

values with σ

decreases with increasing sampling rate. Variation of

and

with σ of WGN

,

and RN at different sampling rates are shown in Figure 12. Now, the performance of the proposed ECG denoising technique on real-time ECG signals, which are corrupted with real noises (i.e., ECG signals, which are not contaminated with synthetic noises) is shown in Figure 13. This figure demonstrates that the proposed technique is not only efficient in removing synthetic noises, it also works proficiently in real-time scenarios. Table 4: Performance of the proposed ECG denoising technique on different ECG databases at different noise levels (γ=99.999%).

nstdb, sampling rate 360 Hz

mitdb, sampling rate 360 Hz

edb, sampling rate 250 Hz

PLI noise is generated using Equation 6, BWN noise is generated using Equation 7, 30dB WGN and RN with μ=0 and σ=0.1 to 5.0. σ of WGN and RN 0.1 37.68 3.15 0.87 0.5 36.78 4.65 1.99 1.0 35.04 5.52 2.87 1.5 33.33 6.64 3.78 2.0 31.89 7.88 4.75 2.5 30.49 9.25 5.74 3.0 29.30 10.63 6.72 3.5 28.18 12.01 7.73 4.0 27.20 13.58 8.72 4.5 26.32 14.99 9.72 5.0 25.56 16.37 10.70 0.1 44.33 4.52 1.22 0.5 44.05 4.63 1.57 1.0 43.39 4.92 2.03 1.5 42.70 5.28 2.51 2.0 41.96 5.71 2.99 2.5 41.24 6.19 3.48 3.0 40.50 6.71 3.98 3.5 39.86 7.18 4.47 4.0 39.24 7.72 4.96 4.5 38.68 8.23 5.46 5.0 38.21 8.67 5.94 0.1 28.13 4.22 1.16 0.5 27.98 4.27 1.41 1.0 27.50 4.46 1.75 1.5 26.80 4.77 2.10 2.0 26.03 5.15 2.45 2.5 25.23 5.59 2.81 24

macecgdb, sampling rate 500Hz ptbdb, sampling rate 1 kHz

3.0 3.5 4.0 4.5 5.0 0.1 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 0.1 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0

24.44 23.71 22.99 22.33 22.53 30.64 30.42 29.52 28.29 27.02 25.77 24.74 23.81 22.98 22.12 21.40 25.16 25.19 25.23 25.28 25.30 25.36 25.43 25.48 25.53 25.62 25.78

6.09 6.60 7.16 7.63 7.68 3.65 3.73 4.08 4.59 5.24 5.97 6.70 7.48 8.19 9.02 9.74 6.18 6.23 6.30 6.32 6.34 6.37 6.39 6.42 6.46 6.62 6.64

3.18 3.54 3.92 4.28 4.30 0.98 1.06 1.18 2.15 2.60 3.09 3.63 4.07 4.44 5.07 5.41 1.82 1.96 2.07 2.13 2.20 2.21 2.26 2.33 2.40 2.42 2.48

Apart from quantitative analysis, qualitative measure of denoised ECG signal is also very important for true quality assessment. Semi-blind mean opinion score (MOS) test [32] of ten evaluators (eight researchers who are working in digital signal, bio-image and bio-signal processing domains and two cardiologists) has been carried out. In doing so, a web-based survey-form has been created containing the original and the corresponding denoised ECG signals, and the evaluators were invited to take part in the survey. Six ECG features (QRS complex, T and P waves, ST segment, PR and QT intervals) are considered for the MOS test. Eleven ECG signals of MIT-BIH arrhythmia database are taken (record #100, #106, #112, #117, #123, #202, #208, #214, #220, #228, 234), and added with synthetic noises (PLI noise is generated using Equation 6, BWN noise is generated using Equation 7, 30dB WGN and RN with μ=0 and σ=0.1 to 5.0). Then, these noisy signals, 121 noisy signals in total (each ECG record is added with WGN and RN at 11 different values of σ along with PLI and BWN), are denoised using the proposed technique setting γ at 99.999%. In the survey-form the evaluators were asked to compare the original ECGs with their denoised counterparts, and give a score in accordance

25

with their similarities. The quality ratings are 1 (completely different), 2, 3, 4, and 5 (identical). MOS of the denoised ECG signals at different values of σ are calculated using the following equation. ∑∑ where Total number of evaluators, the feature given by the evaluator.

Total number of features,

Quality rating of

The MOS value of a particular ECG feature is expressed as ∑ Gold standard MOS error of each feature and the overall denoised ECG signal are calculated using (14) and (15), respectively. ( (

) )

Table 5 shows the MOS errors of overall denoised ECG signals at different values of σ on MITBIH arrhythmia database.

Table 5: MOS errors of denoised ECG signals (overall). σ (%) σ (%)

0.1 3.33 3.0 6.00

0.5 3.33 3.5 11.33

1.0 4.67 4.0 10.67

1.5 3.33 4.5 11.33

2.0 2.5 5.00 4.67 5.0 13.33

From Table 5 it can be noted that the MOS error values increase with increasing values of σ. However, according to the gold standard MOS error criteria [32], which says that an ECG signal can be considered as of „very good‟ quality if the value of

is found to be ≤15, the denoised

ECG signals (overall) fall under the category „very good‟. Values of

at different values of

σ are also calculated, and are found to be „very good‟ too as per the criteria.

26

Run time of an algorithm or technique, which is the function of the input data file size, is an important figure of merit, and is used to gauge the time complexity. The proposed ECG denoising technique is implemented on MATLAB platform with a computer having 64-bits Windows 7 Professional operating system, 16 GB RAM and Intel(R) Xeon(R) CPU E3-1225 v3 3.20 GHz. Run times of the proposed ECG denoising technique as a function of the length of the signal at different sampling rates are evaluated and shown in Figure 14. From this figure it is to be noted that: (i) the time complexity of the proposed technique vary almost linearly with the length of the signal, and (ii) rate of change of time-complexity with the length of the signal is proportional to the sampling rate. Since the technique is primarily implemented on software platform, the space complexity and the power requirement are not considered in this research. However, the run time of proposed technique can further be reduced by efficient codeoptimization and implementation.

Figure 12: Variation of rates.

,

and

with σ of WGN and RN at different sampling

27

Figure 13: (a) Performance of the proposed technique on ECG signal #s0002 taken from PTBDB, (b) performance of the technique on ECG signal #105 taken from MIT-BIH arrhythmia database.

Figure 14: Run times of the proposed ECG denoising technique as a function of length of the signal at different sampling rates.

4. Performance comparisons Performance of an ECG denoising technique depends on a number of important parameters such as the considered bandwidth, sampling frequency at which the signal is acquired, and the added noise level. Different ECG denoising techniques reported in the literature are tested on different 28

types of ECG signals and at different noise levels. Therefore, a straightforward comparison of denoising efficiency among different techniques is not possible. It has already been mentioned before that the performance of SSA in denoising ECG signals is studied in [28]. However, the numerical results of the denoising performance are not provided in detail in [28], and hence, it is not possible to directly compare the performance of the proposed method with that of [28]. Based on the spectral analysis of the principal components, it has been suggested

in

[28]

to

factor

choose

an

L

value,

which

is

at

least

twice

the

. The clinical bandwidth of an ECG signal lies between 0.5 Hz

to 100 Hz [29]. Therefore,

, where

. In

[28], the ECG signals are recorded at 250 Hz and 1250 Hz sampling rates, and then the acquired ECG signals are added with synthetically generated power-line interference noise. Then, the noisy ECG signals are processed with an SSA-based technique for different values of , and the denoising performance at different values of

are studied. It has been shown in [28] that the (1)

power-line interference noise removal efficiency improves with increasing the value of , and (2) the power-line interference noise removal efficiency gets particularly better at are close to the integer multiples of the factor denoising technique, the value of

values, which

. In this proposed SSA-based ECG

is calculated using Equation 2. For an example, if the

sampling rate of the ECG signal is 1250 Hz, then our proposed technique uses a 313

value of

for denoising the ECG signal, which well-satisfies the condition for choosing

the value of

given in [28]. However, the denoising efficiency of the technique proposed in [28]

is not tested on ECG signals corrupted with any other types of noises except the power-line interference noise. The following steps are performed in order to do a fair comparison of performance between the proposed ECG denoising technique with that of the one, which is proposed in [28]. 1) The algorithm proposed in [28] is reproduce. 2) All the 48 ECG signals (lead ML II) of MIT-BIH arrhythmia database are mixed with PLI noise (generated using Equation 6), BWN (generated using Equation 7), and both 30dB WGN and RN at σ =1.5 and μ=0 (the same amount of noises, which are used to testify the performance of the proposed algorithm; please see Table 2).

29

3) The noisy ECG signals are denoised using the reproduced-algorithm at L values of 14, 19, 22, 26, and 90, and finally, the denoising-performance of the reproduced-algorithm is compared to that of the proposed one and shown in Table 6. Table 6: Performance comparison of the proposed ECG denoising technique with that of the one proposed in [28].

MITDB (all 48-ECG records)

Database

Method The reproduced algorithm of [28] Proposed

Value of L 14 19 22 26 90 90

Average

Average

Average

0.00 0.00 0.00 0.00 0.00 42.70

682.67 682.67 682.67 682.67 682.67 5.28

93.53 93.53 93.52 93.52 93.52 2.51

In [28], the SSA-based denoising technique is tested on ECG signals recorded at 250 Hz and 1250 Hz sampling rates with four different values of . The values of , which are used at 250 Hz sampling rate are 10, 13, 15 and 18, and the values of , which are used at 1250 Hz sampling rate are 50, 65, 75 and 90. Following this rule, the value of L for an ECG signal recorded at 360 Hz would be 14, 19, 22 and 26, and the performance of the reproduced algorithm is tested on these values of L, and are shown in Table 6. On the other hand, following Equation 2, our proposed algorithm calculates the value of L=90. The performance of the reproduced algorithm is also tested using an L value of 90. Table 6, clearly demonstrate that the algorithm, which is proposed in [28], might able to expel out the powerline-line-interference noise, but its performance is not satisfactory in presence of other types of noises. The EMD-based ECG denoising technique proposed in [7] is tested only on 800 samples of the ECG data file #115 taken from mitdb. In [7], first, the amplitude of the ECG signal is normalized to ±1, and then the signal is added with four different types of synthetic noises (PLI, WGN, RN and BWN) taking only one type of noise at-a-time. Next, the noisy ECG signal is denoised, and one additional non-diagnostic distortion measure metric, namely, mean square error ( ∑

), is calculated along with PRD to show the efficiency of the

technique. The lower the value of

, the better the performance. For the sake of a fair

comparison, our proposed ECG denoisng technique is also applied on that particular ECG signal at the same level and type of noise contamination. The result is shown and compared to [7] in Table 7.

30

Table 7: Performance comparison of the proposed ECG denoising technique with that of the EMD-based method proposed in [7]. MITDB record #115 (800 samples) Noise PLI types WGN BWN

[7]

12.0526 7.53 7.6487

0.00096 0.00270 0.00260

Proposed

31.90 13.33 29.19

0.000007 0.000761 0.000307

From Table 7, it can be noted that the performance of the proposed ECG denoising technique outperforms the technique reported in [7] in removing PLI, WGN and BWN. The RN generation technique, which is used in [7] is not mentioned properly, and hence, it is not possible to compare the RN removal efficiency of the proposed technique with that of [7]. The EMD and adaptive switching mean filter (ASMF)-based ECG denoising technique proposed in [8] is tested on eight ECG data files (#100, #101, #103, #105, #115, #200, #215, #230) taken from mitdb. In [8], first, the amplitude of the ECG signal is normalized to ±1, and then they are added with synthetic RN in the frequency band 0.3 Hz - 2 Hz at different levels of input SNR. Finally, the noisy signals are denoised. For the sake of fair comparison, our proposed denoising technique is also applied on these ECG signals after adding synthetic levels.

at same input SNR

removal efficiency of the proposed technique is compared to that of [8] in Table 8.

From this Table it is interesting to note that the average values of method proposed in [8] in removing

and

of the

is about 3 and 2 times, respectively, better than that of

the proposed technique. However, the average value of

obtained using the proposed

technique is about 3 times better than that of [8]. Table 8: Performance comparison of the proposed ECG denoising technique with that of the EMD and ASMF-based method proposed in [8]. [8]

Proposed

Noise type: RN

in dB 0 5 10 15 20

9.2980 9.1351 8.7879 8.0516 5.6733

34.32 18.95 11.53 7.12 5.30

0.02022 0.00655 0.00232 0.00088 0.00050

3.0388 3.0125 2.9413 2.7238 1.9713

70.22 39.55 22.06 12.80 7.83

0.00903 0.00285 0.00090 0.00030 0.00011

The fractional zero-phase filtering (FZP)-based ECG denoising technique, which is proposed in [33] is tested only on the ECG record #115 taken from mitdb. In [33], first, the ECG signal is

31

added with PLI and 15dB WGN taking one type of noise at-a-time and then the noisy signal is denoised. Frequency and peak-to-peak amplitude of the added synthetic PLI noise, which is considered in [33], are 50 Hz and 0.15 mV, respectively. In [33], the power of the added WGN (i.e., 15 dB) is mentioned, but the standard deviation value, which is used to generate the WGN is not mentioned. Therefore, it is difficult to directly compare the WGN removal efficiency of the proposed technique with that of [33]. For the sake of fair comparison, the proposed technique is modified so as to obtain a

value, which is close to that of [33], and then the

performance is compared and shown in Table 9. Table 9: Performance comparison of the proposed ECG denoising technique with that of the FZP-based method proposed in [33]. MITDB record #115 Noise PLI types WGN

[33]

Proposed

14.2565 13.4716

0.0128 0.0153

34.0500 13.4576

0.000006 0.000733

It can be noted from Table 9, that the PLI removal efficiency of the proposed technique outperforms [33]. In removing WGN, our proposed technique achieves a 13.4576 by setting σ at 0.039, and the obtained

value of

value is very close to that in [33]. Now,

from Table 9, it can be seen that at almost same level of

the

value obtained using

the proposed technique is approximately 21 times better than that of [33]. In [17], an adaptive Fourier decomposition-based (AFD)-based ECG denoising technique is proposed. In the performance evaluation phase of [17], 16 ECG signals are taken from mitdb, and are mixed with WGN at an input SNR level of 10dB. Then, those noisy signals are filtered using the AFD-based technique. Performance of our proposed ECG denoising technique is also tested on these 16 ECG signals at the same level of noise contamination, and the obtained values are shown and compared to [17] in Table 10. Table 10 clearly demonstrate that the proposed ECG denoising technique outperforms [17]. Table 10: Performance comparison of the proposed ECG denoising technique with that of the AFD-based method proposed in [17].

Noise type: WGN

MIT-BIH record # 101 102 103 104

[17] 36.0 32.6 76.0 55.7

MSE Proposed 0.002654 0.004392 0.003356 0.003058 32

105 106 107 108 109 201 202 203 205 207 208 209

73.6 98.9 572.6 24.0 112.1 38.3 28.4 321.3 29.5 93.9 199.2 60.8

0.004169 0.003828 0.012570 0.002472 0.009362 0.005129 0.003801 0.007476 0.003213 0.008841 0.007246 0.005104

A dictionary-learning-based ECG denoising technique is proposed by Sajita et al. in [15], and the performance of the technique is evaluated on a number of ECG signals taken from mitdb with different combinations of various synthetically generated noises. An additional non-diagnostic distortion measure metric, namely, maximum absolute amplitude ( ), is calculated along with value of

to show the efficiency of the technique. The lower the

, the better the performance. In Table 11, the efficiency of the proposed ECG

denoising technique is compared to that of [15] in removing a combination of BWN and PLI. It can be seen from this table that our proposed technique performs better than that of [15] in most of the cases. Table 11: Performance comparison of the proposed ECG denoising technique with that of the dictionary-learning-based method proposed in [15].

Noise type: BWN+PLI

MIT-BIH record # 102 109 111 119 124 203 208

[15] 20.90 23.75 22.09 23.83 22.43 23.61 23.33

Proposed 0.21 0.35 0.21 0.30 0.33 0.27 0.34

21.23 33.39 24.61 23.24 37.88 23.95 25.65

0.2278 0.2302 0.2367 0.2242 0.2326 0.2300 0.2301

5. Conclusion and discussion A reliable, robust and high performance singular spectrum analysis-based ECG denoising technique has been proposed in this paper. The performance of the technique is tested on a large number of ECG data files taken from five different publicly available ECG databases. The reconstructed ECG signal is the noise-reduces and smoothed version of the original ECG signal.

33

Since the reconstructed components are filtered using zero-phase bi-directional bandpass and notch filters, the phase difference, i.e., time-lag between the original and reconstructed signal is zero. The smoothing and denoising performance of the proposed ECG denoising technique on real-time ECG signals can be visually examined in Figure 13. ECG tracings are generally scrutinized for calculating the heart-rate and heart-rate-variability, identifying the presence of different types of arrhythmias as well as other cardiovascular disorders such as myocardial infraction, bundle branch block. Therefore, the reconstructed, i.e., the denoised ECG signal must be of the same amplitude and time scales of the original ECG signal. As the singular spectrum analysis technique takes the original time-domain ECG signal as input, and provides the reconstructed components at the output in the same domain, there happens to be neither time nor amplitude-scale variation with the denoised ECG signal compared to the original one. The novelty of the proposed ECG-denoising technique lies in its performance. From Table 2 and also from Figures 4 to 7 it can be noted that proposed technique is robust enough to bring out the ECG signal from a combination of different types of noises including white Gaussian and random types of noises. Further, Tables 6 to 11 show the competency of the proposed technique in removing each of the four different types of noises (baseline wander, power-line interference, white Gaussian and random noises), and also its superiority in comparison with other state-of-the art methods. Not only the proposed technique performs well on synthetically-added noises, it is also ace in expelling out the real noises. Figures 8 and 13 portray the efficiency of the technique in removing real noises. Both the quantitative and qualitative distortion measure metrics show that the proposed ECG denoising technique is robust enough to filter various noises present in the signal without jeopardizing the clinical content. It is also possible that the proposed technique could be adapted for denoising other biomedical signals exhibiting periodic or quasi-periodic nature such as photoplethysmogram, esophageal pressure signal. From Table 1 it can be noted that the proposed technique is tested on all the available leads of ECG signal of all the five databases, except MITBIH arrhythmia database. There are two leads of ECG signals present in MIT-BIH arrhythmia database: ML II and V5. However, the performance of other ECG denoising techniques ([7], [8], [15], [17] and [33]), which are mentioned in Section IV, are reported only on lead ML II of MIT-

34

BIH arrhythmia database. Hence, in this research the lead ML II is taken to expedite the comparison of performance to other techniques. It could be said that the (i) optimum value of L, (ii) enhanced covariance among the eigenvector of the matrix

and also (iii) carrying out the bandpass and notch filtering operations at the

reconstruction-component level are the key reasons for getting such an excellent ECG denoising performance. The main advantages of using this proposed technique are: (i) the singular spectrum analysis parameters are made adaptive to the sampling rate of the ECG signal, and therefore the technique does not require any manual intervention, and (ii) the proposed technique is a model-free/function-independent method of ECG denoising. It captures the periodicity of the oscillatory modes in ECG signals better than that of the fixed filtering-based approaches. Also in the case of wavelet transform-based methods the mother wavelet function and the number of decomposition levels are to be selected manually, and moreover, the number of decomposition levels are also to be adjusted manually for variations in sampling rate. The threshold value (γ) could be set to either 99.999% or 100%, and the quality of the denoised signal is identical for both the threshold values. However, the optimum number of filtered reconstruction-components,

, required to reconstruct the signal is much higher at γ=100%

than that of at 99.999%. This is the reason for which γ is set to 99.999% in this research. Figures 4 to 7 show that the proposed technique is able to efficiently process different types of ECG signals irrespective of their QRS-pattern and clinical morphologies. From Figure 2, it can be seen that the coefficients of certain filtered reconstruction-components clearly indicate the QRS complexes regions. Therefore, the proposed technique, alongside denoising the ECG, can also be adapted for detecting its QRS complexes. Although a few reliable quantitative ECG distortion measure metrics, such as percent root mean square difference, are there, nonetheless, the result of subjective assessment test provides better insight into the quality of the denoised ECG signals as the test is conducted by the cardiologists. The more the ECG signals could be included in the MOS test, the better the performance analysis would be. In this research the denoising technique of a single-lead ECG signal using SSA is proposed. The 1D ECG signal is converted to 2D through proper embedding operation. The technique can also be revamped and extended in order to be applicable on tensor data so as to denoise multi-lead ECG signals. This will be the topic of our future research work. 35

Acknowledgements The authors convey their sincere thankfulness to Dr. Vijay S Chauhan, MD, FRCPC (Cardiologist, Toronto General Hospital, Toronto, Canada), Dr. Uma Mandal, M.B.B.S. (College of Medicine and JNM Hospital, Kolkata, India), and the graduate students of the Department of Electrical, Computer and Biomedical Engineering, Ryerson University, for their precious time and effort in carrying out the subjective quality estimation. We would also thank the funding provided by the Natural Sciences and Engineering Research Council (NSERC) Discovery Grant program of Canada.

Conflict of interest statement The authors report no conflicts of interest. The authors alone are responsible for the content and writing of the paper.

References [1] [2]

[3]

[4]

[5]

[6]

[7]

A.L. Goldberger, Clinical Electrocardiography, a Simplified Approach, 7 th ed., New York, NY, Elsevier, USA, 2006. M. L. Ahlstrom and W. J. Tompkins, “Automated High-Speed Analysis of Holter Tapes with Microcomputers”, IEEE Transactions on Biomedical Engineering, vol. 30, iss. 10, pp. 651-657, Oct. 1983. Z. A. A. Alyasseri, A. T. Khader, M. A. Al-Betar and M. A. Awadallah, “Hybridizing βhill climbing with wavelet transform for denoising ECG signals”, Information Sciences, vol. 429, iss. 8, pp. 229-246, March 2018. H. Liu, Y. Li, Y. Zhou, X. Jing and T. K. Truong, “Joint power line interference suppression and ECG signal recovery in transform domains”, Biomedical Signal Processing and Control, vol. 44, iss. 6, pp. 58-66, July 2018. A. Kumar, R. Komaragiri and M. Kumar, “Design of wavelet transform based electrocardiogram monitoring system”, ISA Transactions, vol. 80, iss. 9, pp. 381-398, Sep. 2018. P. Singh, G. Pradhan and S. Shahnawazuddin, “Denoising of ECG signal by non-local estimation of approximation coefficients in DWT”, Biocybernetics and Biomedical Engineering, vol. 37, iss, 3, pp. 599-610, 2017. S. Jain, V. Bajaj and A. Kumar, “Riemann Liouvelle Fractional Integral Based Empirical Mode Decomposition for ECG Denoising”, IEEE Journal of Biomedical and Health Informatics, vol. 22, no. 4, pp. 1133-1139, July 2018. 36

[8]

[9] [10]

[11]

[12]

[13]

[14]

[15]

[16] [17]

[18] [19]

[20]

[21]

M. Rakshit and S. Das, “An efficient ECG denoising methodology using empirical mode decomposition and adaptive switching mean filter”, Biomedical Signal Processing and Control, vol. 40, iss. 2, pp. 140-148, Feb. 2018. S. Jain, V. Bajaj and A. Kumar, “Effective de-noising of ECG by optimised adaptive thresholding on noisy modes”, IET Sci. Meas. Technol., vol. 12, iss. 5, pp. 640-644, 2018. P. Nguyen and J. M. Kim, “Adaptive ECG denoising using genetic algorithm-based thresholding and ensemble empirical mode decomposition”, Information Sciences, vol. 373, iss. 35, pp. 499-511, Dec. 2016. S. Kumar, D. Panigrahy and P.K. Sahu, “Denoising of Electrocardiogram (ECG) signal by using empirical mode decomposition (EMD) with non-local mean (NLM) technique”, Biocybernetics and Biomedical Engineering, vol. 38, iss, 2, pp. 297-312, 2018. D. P. Tobon and T. H. Falk, “Adaptive Spectro-Temporal Filtering for Electrocardiogram Signal Enhancement”, IEEE Journal of Biomedical and Health Informatics, vol. 22, no. 2, pp. 421-428, March 2018. H. D. Hesar and M. Mohebbi, “ECG Denoising Using Marginalized Particle Extended Kalman Filter With an Automatic Particle Weighting Strategy”, IEEE Journal of Biomedical and Health Informatics, vol. 21, no. 3, pp. 635-644, July 2018. H. D. Hesar and M. Mohebbi, “An Adaptive Particle Weighting Strategy for ECG Denoising Using Marginalized Particle Extended Kalman Filter: An Evaluation in Arrhythmia Contexts”, IEEE Journal of Biomedical and Health Informatics, vol. 21, no. 6, pp. 1581-1592, Nov 2017. U. Satija, B. Ramkumar and S. Manikandan, “Noise-aware dictionary-learning-based sparse representation framework for detection and removal of single and combined noises from ECG signal”, Healthcare Technology Letters, vol. 4, iss. 1, pp. 2-12, Feb. 2017. J. Zhu and X. Li, “Electrocardiograph signal denoising based on sparse decomposition”, Healthcare Technology Letters, vol. 4, iss. 4, pp. 134-137, Aug. 2017. Z. Wang, F. Wan, C. M. Wong and L. Zhang, “Adaptive Fourier decomposition based ECG denoising”, Computers in Biology and Medicine, vol. 77, iss. 10, pp. 195-205, Oct. 2016. J. Wang, Y. Ye, X. Pan, X. Gao and C. Zhuang, “Fractional zero-phase filtering based on the Riemann-Liouville integral”, Signal Processing, vol. 98, iss. 6, pp. 150-157, May 2014. Y. Lee and D. Hwang, “Periodicity-based nonlocal-means denoising method for electrocardiography in low SNR non-white noisy conditions”, Biomedical Signal Processing and Control, vol. 39, iss. 1, pp. 284-293, Jan. 2018. S. Cuomo, R. Farina and F. Piccialli, “An inverse Bayesian scheme for the denoising of ECG signals”, Journal of Network and Computer Applications, vol. 115, iss. 1, pp. 48-58, Aug. 2018. M. Seo, M. Choi, J. S. Lee and S. W. Kim, “Adaptive Noise Reduction Algorithm to Improve R Peak Detection in ECG Measured by Capacitive ECG Sensors”, Sensors, vol. 18, iss. 7, pp. 1-16, July 2018.

37

[22] N. Golyandna and A. Zhigljavsky, Singular Spectrum Analysis for Time Series, 1 st ed., New York, NY, Springer, USA, 2013. [23] L. Zotov and E. Scheplova, “MSSA of globally gridded OAM from ECCO, AAM from ECMWF, and gravity from GRACE”, in Proc. of 3rd International Conference on Digital Information Processing, Data Mining, and Wireless Communications (DIPDMWC), 6-8 July 2016, Moscow, Russia, pp. 127-132. [24] M. L. de Menezes, R. C. Souza and J. F. M. Pessanha, “Electricity consumption forecasting using singular spectrum analysis”, Dyna, vol. 82, no. 190, pp. 138-146, 2015. [25] S. Wang et al., “Localizing Microaneurysms in Fundus Images Through Singular Spectrum Analysis”, IEEE Transactions on Biomedical Engineering, vol. 64, no. 5, pp. 990-1002, May 2017. [26] S. Ferdowsi and V. Abolghasemi, “Semiblind Spectral Factorization Approach for Magnetic Resonance Spectroscopy Quantification”, IEEE Transactions on Biomedical Engineering, vol. 65, no. 8, pp. 1717-1724, Aug 2018. [27] D. Jarchi, C. Wong, R. M. Kwasnicki, B. Heller, G. A. Tew and G. Z. Yang, “Gait Parameter Estimation From a Miniaturized Ear-Worn Sensor Using Singular Spectrum Analysis and Longest Common Subsequence”, IEEE Transactions on Biomedical Engineering, vol. 61, no. 4, pp. 1261-1273, April 2014. [28] B. Yang, Y. Dong, C. Yu and Z. Hou, “Singular Spectrum Analysis Window Length Selection in Processing Capacitive Captured Biopotential Signals”, IEEE Sensors Journal, vol. 16, no. 19, pp. 7183-7193, Oct. 2016. [29] S. K. Mukhopadhyay, S. Mitra and M. Mitra, “ECG Feature Extraction Using Differentiation, Hilbert Transform, Variable Threshold and Slope Reversal Approach”, Journal of Medical Engineering and Technology, vol. 36, no. 7, pp. 372-386, Oct. 2012. [30] PhysioNet. [Online]. Available: https://physionet.org/cgi-bin/atm/ATM, accessed on: June. 2019. [31] M.S. Manikandan and S. Dandapat, “Wavelet energy based diagnostic distortion measure for ECG”, Biomedical Signal Processing and Control, vol. 2, iss. 2, pp. 80-96, April 2007. [32] Y. Zigel, A. Cohen and A. Katz, “The weighted diagnostic distortion (WDD) measure for ECG signal compression,” IEEE Transactions on Biomedical Engineering, vol. 47, no. 11, pp. 1422-1430, Nov. 2000. [33] J. Wang, Y. Ye, X. Pan and X. Gao, “Parallel-type fractional zero-phase filtering for ECG signal denoising”, Biomedical Signal Processing and Control, vol. 18, iss. 4, pp. 36-41, April, 2015. [34] N. Golyandina and A. Zhigljavsky, Singular Spectrum Analysis for Time Series. Berlin, Germany: Springer-Verlag, 2013, pp. 25-32. [35] T. J. Harris and H. Yuan, “Filtering and frequency interpretations of singular spectrum analysis,” Physica D, Nonlinear Phenomena, vol. 239, iss. 20–22, pp. 1958-1967, Oct. 2010.

38

[36] D. Skillicorn, “Understanding Complex Datasets: Data Mining with Matrix Decompositions”, 1st ed., New York, NY, Chapman and Hall/CRC, USA, 2007. [37] S. W. Chen, H. C. Chen, and H. L. Chan, “A real-time QRS detection method based on moving-averaging incorporating with wavelet denoising”, Computer Methods and Programs in Biomedicine, vol. 82, iss. 3, June 2006, pp. 187-195. [38] O. Haavisto, “Detection and analysis of oscillations in a mineral flotation circuit”, Control Engineering Practice, vol. 18, iss. 1, Jan 2010, pp. 23-30. [39] D. G. Manolakis, V. K. Ingle, S. M. Kogon, Statistical and Adaptive Signal Processing: Spectral Estimation, Signal Modeling, Adaptive Filtering and Array Processing. Norwood, MA, USA: Artech House, 2005. [40] J. Harmouche, D. Fourer, F. Auger, P. Borgnat, P. Flandrin, “The Sliding Singular Spectrum Analysis: A Data-Driven Nonstationary Signal Decomposition Tool” IEEE Transactions on Signal Processing, vol. 66, No. 1, Jan 2018.

39