Submit manuscript...
MOJ
eISSN: 2374-6920

Proteomics & Bioinformatics

Research Article Volume 7 Issue 2

Analysis of heart rate variability based on quantitative approach

Hla Myo Tun,1 Zaw Min Naing,1 Win Khaing Moe,2 Maung Maung Latt3

1Department of Electronic Enginering, Yangon Technological University, Myanmar
2Department of Research and Innovation, Ministry of Education, Myanmar
3Technological University (Toungoo), Ministry of Education, Myanmar

Correspondence: Hla Myo Tun, Department of Electronic Enginering, Yangon Technological University, Myanmar

Received: February 26, 2018 | Published: April 30, 2018

Citation: Tun HM, Naing ZM, Moe WK, et al. Analysis of heart rate variability based on quantitative approach. MOJ Proteomics Bioinform. 2018;7(2):131 ?1 41. DOI: 10.15406/mojpb.2018.07.00223

Download PDF

Abstract

Heart rate variability (HRV) is a measure of variations of heart rate between two successive heart beats and it is a relatively new method for assessing the effects of stress on our body. It is measured as the time gap between our heart beats that varies as we breathe in and out. Simple measures of the small changes in each beat of our heart can provide a wealth of information on the health of our heart and nervous system; such measures are called heart rate variability or HRV. It is usually calculated by analyzing the time series of beat-to-beat intervals from ECG signal. HRV is measured based on variation of time in milliseconds between two heartbeats also known as RR interval where R is a point corresponding to the peak of the QRS complex of the ECG wave, and RR is the interval between successive Rs. The evaluation of HRV can provide an indication of cardiovascular health. To accomplish this evaluation, the raw ECG signal is firstly de-noised by the application of DWT by automatically determining the optimal order of decomposition. After the purification, the wavelet filtering is used for R-peak detection since this method is efficient and accurate in the compute of the R peaks positions without changing of the shape or position of the original signal. Finally, the RR interval is analyzed because it consists in studying the HRV. The values of RR intervals are then plotted versus time, giving a curve called RR tachogram. After all, the Lomb-Scargleperiodogram of frequency-domain method (spectral analysis) is used to investigate the sympathovagal balance of HRV from RR tachogram. By using the Lomb-Scargleperiodogram for power spectral density estimation, we have no need to make de-trending and re-sampling. As a result, all signal shows as arrhythamia by comparing with the normal value of standard measurement. And then, the traditional time-domain method (statistical analysis) and spectral analysis of frequency-domain method are applied to analyze the variation of heart rate in arrhythmia database. The MATLAB programming is used to implement the algorithm for HRV analysis and MIT/BIH arrhythmia databases is used as data inputs.

Keywords: heart rate variability, arrhythmia database, autonomic nervous system, sinoatrial, parasympathetic

Abbreviations

HRV, Heart rate variability; ANS, autonomic nervous system; SA, sinoatrial; LF, low-frequency; HF, high-frequency;

Introduction

Heart rate variability (HRV) is the analysis of variations in the instantaneous heart rate in time series using beat-to-beat RR interval. HRV is considered as an indicator of the activity of autonomic regulation of circulatory function and as the one of the most significant methods of analyzing the activity of the autonomic nervous system (ANS). There are many methods that are used to detect heartbeats, such as ECG and blood pressure. Among them, ECG is considered as the best method to detect HRV signals since it provides a clear waveform, which makes it easier to exclude heartbeats not originating in the sinoatrial (SA) node.

The heart rate may be increased by slow acting sympathetic activity or decreased by fast acting parasympathetic (vagal) activity. The balance between the effects of the sympathetic and parasympathetic systems, the two opposite acting branches of the autonomic nervous system, is referred to as the sympathovagal balance and the reflected in the beat-to-beat changes of the cardiac cycle. Spectral analysis of the RR tachogram is typically used to estimate the effect of the sympathetic and parasympathetic modulation of the RR intervals. The two main frequency bands of interest are referred to as the low-frequency (LF) band (0.04 to 0.15 Hz) and the high-frequency (HF) band (0.15 to 0.4Hz).The ratio of the power contained in the LF and HF components has been used as a measure of the sympathovagal balance. The analysis of HRV signals is an important tool for studying the autonomic nervous system, as it allows the evaluation of the balance between the sympathetic and parasympathetic influences on heart rhythm.

Time-domain method is the simplest method to perform. Simple time-domain variables that can be calculated include the mean NN interval, the mean heart rate, the difference between the longest and shortest NN interval. In frequency domain methods, various spectral methods are used for the analysis of RR tachogram. Since RR intervals in time series are non-stationary and they are spaced unevenly, it needs to make re-sampling. The Lomb-Scargleperiodogram can overcome this weakness and it can estimate the PSD directly from irregularly sampled RR interval series without resampling. The most common classifications for HRV spectrum are: very low frequency (VLF: 0.003-0.04Hz), low frequency (LF: 0.04-0.15Hz) and high frequency (HF: 0.15-0.4Hz). The ratio of LF to HF appears to be a sensitive measure of the autonomic nervous system’s (ANS) response to a sudden change in cardiovascular control and this ratio represents an evaluation of the ANS balance.1‒10

Overview of heart rate variability

In this research, there are four major steps to implement the analysis of heart rate variability namely preprocessing, R-peak detection, and HRV analysis in time-domain and frequency-domain method.Raw ECG is taken from online MIT-BIH database, each data contains different ECG wave shape; each portion of the ECG waveform carries information that about the cardiac conditions. Preprocessing means the de-noising of raw ECG signal. During the recording of ECG signals, different types of noise from various sources can be superimposed to the original signal [12Nil]. In this reason, wavelet method is applied for the removal of baseline drift. Once the signal is free from noise the discrete wavelet transform is applied to extract the R-peak location and the amplitude in mV (millivolts) [16Sav, 12Say]. The accuracy of the determined location of R-peak is essential for HRV analysis.

After detecting the R points, the RR interval analysis is performed. HRV metrics are calculated from an analysis of the RR tachogram, the time series of RR intervals. There are severalapproaches for HRV analysis, which can be subdivided into linear analysis (time-domain and frequency-domain) and non-linear analysis (wavelet). In this paper, the time-domain and frequency-domain method of linear analysis are used. The standard measurements of time-domain method are SDNN (ms), SDANN (ms), RMSSD (ms), NN50 (count), and PNN50 (%). These variables directly derived from the beat-to-beat intervals, such as the mean heart rate (HR) and the standard deviation (SD) for the entire record. The idea for the frequency domain methods is density (PSD) by using Lomb-Scargleperiodogram. After all, these results are verified with the normal values of standard measures of HRV.11‒24

Implementation of HRV analysis

Autonomic nervous system (ANS) plays an important role in the regulation of the physiological processes of the human organism during normal and pathological conditions. Among the techniques used in its evaluation, the heart rate variability (HRV) has arising as a simple and non-invasive measure of the autonomic impulses, representing one of the most promising quantitative markers of the autonomic balance. The HRV describes the oscillations in the interval between consecutive heart beats (RR interval), as well as the oscillations between consecutive instantaneous heart rates. It is a measure that can be used to assess the ANS modulation under physiological conditions, such as wakefulness and sleep conditions, different body positions, physical training and also pathological conditions. Changes in the HRV patterns provide a sensible and advanced indicator of health involvements. Higher HRV is a signal of good adaptation and characterizes a health person with efficient autonomic mechanisms, while lower HRV is frequently an indicator of abnormal and insufficient adaptation of the ANS, provoking poor patient’s physiological function. Because of its importance as a marker that reflects the autonomic nervous system activity on the sinus node and as a clinical instrument to assess and identify health involvements, this study reviews conceptual aspects of the HRV, measurement devices, filtering methods, indexes used in the HRV analyses [13 Hem], limitations in the use and clinical applications of the HRV.

Procedure of HRV analysis

Analysis of Heart rate variability (HRV) provides a non-invasive method to assess the neuronal influences on the cardiaoregulatory function. Since as defined before HRV is the fluctuation of RR intervals, these physiological fluctuations reflect the nonlinear feedback control systems created by the interaction between sympathetic and parasympathetic activities. The HRV analysis can be processed as shown in Figure 1. It includes four major steps namely signal de-noising, R-peak detection and quantification of HRV.

Figure 1 Patient 1: Pre operative imaging showing caudate lobe abscess.

Procedure for ECG signal de-noising

A signal is often corrupted by noise during its acquisition or transmission. The de-noising process is to remove the noise while retaining and not distorting the quality of the processed signal. The traditional way of signal de-noising is filtering. From few years a lot of research about non-linear method of signal de-noising has been developed. These methods are mainly based on thresholding the Discrete Wavelet transform (DWT) coefficients. Simple de-noising algorithm that used DWT consists of three steps as shown in Figure 2.

Figure 2 Block Diagram for ECG Signal De-noising.

Firstly, DWT is adopted to decompose the noisy signal and get the wavelet coefficients. And these wavelet coefficients are de-noised with wavelet threshold. Finally, inverse transform is applied to the modified coefficients and get de-noised signal. The equations of ECG signal de-noising with DWT method is as follows.

The DWT of a signal x(t) is given by:

X DW T K = x( t ) 2 m/2 ψ( 2 m tk )dt MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGybWdamaaBaaajuaibaWdbiaadseacaWGxbGaamivaKqb a+aadaWgaaqcfasaa8qacaWGlbaapaqabaaajuaGbeaapeGaeyypa0 ZaaybCaeqajuaipaqaa8qacqGHsislcqGHEisPa8aabaWdbiabg6Hi Lcqcfa4daeaapeGaey4kIipaaiaadIhadaqadaWdaeaapeGaamiDaa GaayjkaiaawMcaaiaaikdapaWaaWbaaeqajuaibaWdbiaad2gacaGG VaGaaGOmaaaajuaGcqaHipqEdaqadaWdaeaapeGaaGOma8aadaahaa qabKqbGeaapeGaamyBaaaajuaGcaWG0bGaeyOeI0Iaam4AaaGaayjk aiaawMcaaiaadsgacaWG0baaaa@57BE@  (1)

where ψ m,k MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacqaHipqEpaWaaSbaaKqbGeaapeGaamyBaiaacYcacaWGRbaa juaGpaqabaaaaa@3C10@ is the wavelet function.

The approximation and details coefficients are respectively defined by the following Equations 4.2 and 4.3.

c A j+1 ( k )= n= h( n2k )c A j ( k ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGJbGaamyqa8aadaWgaaqcfasaa8qacaWGQbGaey4kaSIa aGymaaqcfa4daeqaa8qadaqadaWdaeaapeGaam4AaaGaayjkaiaawM caaiabg2da9maawahabeqcfaYdaeaapeGaamOBaiabg2da9iabgkHi Tiabg6HiLcWdaeaapeGaeyOhIukajuaGpaqaa8qacqGHris5aaGaam iAamaabmaapaqaa8qacaWGUbGaeyOeI0IaaGOmaiaadUgaaiaawIca caGLPaaacaWGJbGaamyqa8aadaWgaaqcfasaa8qacaWGQbaajuaGpa qabaWdbmaabmaapaqaa8qacaWGRbaacaGLOaGaayzkaaaaaa@551C@                                                                           (2)

c D j+1 ( k )= n= g( n2k )c A j ( k ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGJbGaamira8aadaWgaaqcfasaa8qacaWGQbGaey4kaSIa aGymaaqcfa4daeqaa8qadaqadaWdaeaapeGaam4AaaGaayjkaiaawM caaiabg2da9maawahabeqcfaYdaeaapeGaamOBaiabg2da9iabgkHi Tiabg6HiLcWdaeaapeGaeyOhIukajuaGpaqaa8qacqGHris5aaGaam 4zamaabmaapaqaa8qacaWGUbGaeyOeI0IaaGOmaiaadUgaaiaawIca caGLPaaacaWGJbGaamyqa8aadaWgaaqcfasaa8qacaWGQbaajuaGpa qabaWdbmaabmaapaqaa8qacaWGRbaacaGLOaGaayzkaaaaaa@551E@                                                                           (3)

where j is the decomposition level.

The reconstruction phase begins with the oversampling data followed by two reconstruction filters  and  to obtain the original signal as expressed by the Equation 4.

X IDWT ( t )= m= k= X k m 2 m/2 Ψ( tk )dt MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiwamaaBa aaleaacaWGjbGaamiraiaadEfacaWGubaabeaakmaabmaabaGaamiD aaGaayjkaiaawMcaaiabg2da9maapehabaWaaabCaeaadaaeWbqaai aadIfadaqhaaWcbaGaam4Aaaqaaiaad2gaaaaabaGaam4Aaiabg2da 9iabgkHiTiabg6HiLcqaaiabg6HiLcqdcqGHris5aaWcbaGaamyBai abg2da9iabgkHiTiabg6HiLcqaaiabg6HiLcqdcqGHris5aaWcbaGa eyOeI0IaeyOhIukabaGaeyOhIukaniabgUIiYdGccaaIYaWaaWbaaS qabeaadaWcgaqaaiaad2gaaeaacaaIYaaaaaaakiabfI6aznaabmaa baGaamiDaiabgkHiTiaadUgaaiaawIcacaGLPaaacaWGKbGaamiDaa aa@6197@                                                                         (4)

Estimation of threshold value and noise level

In the Wavelet transform method, the estimation of threshold value and noise level includes in the essential rule. There are various rules for thresholding such as soft thresholding and hard thresholding. In the research, global threshold calculating rule is usedand soft thresholding is applied for de-noising. The mean  from the absolute values of the first detail coefficients c D 1 ( k ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGJbGaamiramaaBaaajuaibaGaaGymaaqabaqcfa4aaeWa aeaacaWGRbaacaGLOaGaayzkaaaaaa@3C66@ then the standard deviation can be calculated by the Equations 5 and 6.
μ= 1 N Abs((c D 1 )) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacqaH8oqBcqGH9aqpdaWcaaWdaeaapeGaaGymaaWdaeaapeGa amOtaaaadaGfWbqabKqbG8aabaWdbiabgkHiTiabg6HiLcWdaeaape GaeyOhIukajuaGpaqaa8qacqGHris5aaGaamyqaiaadkgacaWGZbGa aiikaiaacIcacaWGJbGaamira8aadaWgaaqcfasaa8qacaaIXaaapa qabaqcfa4dbiaacMcacaGGPaaaaa@4B35@                                                                                                              (5)
Where N is the number of the detail coefficients c D 1 MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGJbGaamiramaaBaaajuaibaGaaGymaaqabaaaaa@395F@ .

σ= μ 0.6745 MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacqaHdpWCcqGH9aqpdaWcaaWdaeaapeGaeqiVd0gapaqaa8qa caaIWaGaaiOlaiaaiAdacaaI3aGaaGinaiaaiwdaaaaaaa@3FDC@                                                                                                             (6)

Where 0.6745 is an empirical value used to calibrate the mean with standard deviation for a Gaussian process.

The threshold value is obtained by Equation 7.

S=σ 2Ln( N ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGtbGaeyypa0Jaeq4Wdm3aaOaaa8aabaWdbiaaikdacaWG mbGaamOBamaabmaapaqaa8qacaWGobaacaGLOaGaayzkaaaabeaaaa a@3F70@  (7)
Finally, the IDWT can be computed using the new coefficients to reconstruct the de-noised signal.

Procedure for R-peak detection

For R-peak detection, the signal is decomposed into shifted and scaled versions of the original mother wavelet Db (6) and decomposed into level-8 because the wavelet Db6 closely matches with the shape of ECG QRS complex. An analysis of signal using ECG includes decomposition of the signal, thresholding of the wavelet coefficients and reconstruction of the signal using modified wavelet coefficients. And the following steps are managed for R-peak detection.

Step 1: ECG signal is read and the length is calculated.

Step 2: The signal is decomposed using db6 wavelet.

Step 3: 3rd, 4th and 5th detail coefficients are selected, as most energy of the QRS complex is concentrated in these coefficients.

Step 4: The wave is reconstructed using detail coefficient 3, 4, 5 ( D1=d3+d4+d5 ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfa4aaeWaae aaqaaaaaaaaaWdbiaadseajuaicaaIXaqcfaOaeyypa0JaamizaKqb GiaaiodajuaGcqGHRaWkcaWGKbqcfaIaaGinaKqbakabgUcaRiaads gajuaicaaI1aaajuaGpaGaayjkaiaawMcaaaaa@446F@ .

Step 5: A function d4*(d3+d5) /2^n is defined to reduce the oscillatory nature of the signal where d3, d4, d5 are the ,  , detail coefficients and n is the level of decomposition.

Step 6: Derivate up to level 5 is made using the transfer function

H(z)=(T/8)( z 2 2 z 1 +2 z 1 + z 2 ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGibGaaiikaiaadQhacaGGPaGaeyypa0Jaaiikaiaadsfa caGGVaGaaGioaiaacMcacaGGOaGaeyOeI0IaamOEa8aadaahaaqabK qbGeaapeGaeyOeI0IaaGOmaaaajuaGcqGHsislcaaIYaGaamOEa8aa daahaaqabKqbGeaapeGaeyOeI0IaaGymaaaajuaGcqGHRaWkcaaIYa GaamOEa8aadaahaaqabKqbGeaapeGaaGymaaaajuaGcqGHRaWkcaWG 6bWdamaaCaaabeqcfasaa8qacaaIYaaaaKqbakaacMcaaaa@519D@  (8)

Step 7: The differential equation:

y( nT )=( T 8 )( x( nT2T )2x( nTT )+2x( nT+T )+x( nT+2T ) ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWG5bWaaeWaa8aabaWdbiaad6gacaWGubaacaGLOaGaayzk aaGaeyypa0ZaaeWaa8aabaWdbmaalaaapaqaa8qacaWGubaapaqaa8 qacaaI4aaaaaGaayjkaiaawMcaamaabmaapaqaa8qacqGHsislcaWG 4bWaaeWaa8aabaWdbiaad6gacaWGubGaeyOeI0IaaGOmaiaadsfaai aawIcacaGLPaaacqGHsislcaaIYaGaamiEamaabmaapaqaa8qacaWG UbGaamivaiabgkHiTiaadsfaaiaawIcacaGLPaaacqGHRaWkcaaIYa GaamiEamaabmaapaqaa8qacaWGUbGaamivaiabgUcaRiaadsfaaiaa wIcacaGLPaaacqGHRaWkcaWG4bWaaeWaa8aabaWdbiaad6gacaWGub Gaey4kaSIaaGOmaiaadsfaaiaawIcacaGLPaaaaiaawIcacaGLPaaa aaa@60AA@ (9)

is applied to the signal, by using the transfer function and taking the amplitude response as

| H( ωT ) |=( T 4 )[ sin( 2ωT )+2sin( ωT ) ] MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qadaabdaWdaeaapeGaamisamaabmaapaqaa8qacqaHjpWDcaWG ubaacaGLOaGaayzkaaaacaGLhWUaayjcSdGaeyypa0ZaaeWaa8aaba Wdbmaalaaapaqaa8qacaWGubaapaqaa8qacaaI0aaaaaGaayjkaiaa wMcaamaadmaapaqaa8qacaWGZbGaamyAaiaad6gadaqadaWdaeaape GaaGOmaiabeM8a3jaadsfaaiaawIcacaGLPaaacqGHRaWkcaaIYaGa am4CaiaadMgacaWGUbWaaeWaa8aabaWdbiabeM8a3jaadsfaaiaawI cacaGLPaaaaiaawUfacaGLDbaaaaa@564D@ (10)

Step 8: The signal is squared point by point using the equation
 (11)

to emphasize R wave from the ECG signal.

Step 9: A moving window is integrated using the equation<

Y= ( 1 N ) * [ x( nT( N1 )T )+x( nT( N2 )T )++x( nT ) ] MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGzbGaeyypa0ZaaeWaa8aabaWdbmaalaaapaqaa8qacaaI Xaaapaqaa8qacaWGobaaaaGaayjkaiaawMcaa8aadaahaaqabeaape GaaiOkaaaadaWadaWdaeaapeGaamiEamaabmaapaqaa8qacaWGUbGa amivaiabgkHiTmaabmaapaqaa8qacaWGobGaeyOeI0IaaGymaaGaay jkaiaawMcaaiaadsfaaiaawIcacaGLPaaacqGHRaWkcaWG4bWaaeWa a8aabaWdbiaad6gacaWGubGaeyOeI0YaaeWaa8aabaWdbiaad6eacq GHsislcaaIYaaacaGLOaGaayzkaaGaamivaaGaayjkaiaawMcaaiab gUcaRiabgAci8kabgUcaRiaadIhadaqadaWdaeaapeGaamOBaiaads faaiaawIcacaGLPaaaaiaawUfacaGLDbaaaaa@5C67@ (12)

to obtain the waveform feature information.

Step 10: The threshold value is calculated corresponding to the product of max and mean of the signal to locate the end points of the moving window.

Step 11: The PQRST peaks are located based on the amplitude of the signal within each moving window.

Step12: The time intervals are calculated considering the positions of two consecutive same labeled peaks and stored.

Step 13: Diagnosis of various cardiac diseases is done by comparing ground truth conditions with the data.

HRV Metrics from the RR tachogram

The RR interval analysis consists in studying the Heart Rate Variability (HRV). The HRV is a measure of the variation of heart rate. It is usually calculated by analyzing the time series of beat to beat intervals from ECG (RR interval) or traces of blood pressure. It is obtained by measuring the time between RR intervals on the electrocardiogram. The values of RR intervalsare then plotted versus time, giving a curve called tachogram of HRV or RR tachogram as shown in Figure 3. This tachogram is a combination of sinusoidal waves of different frequencies.

Figure 3 RR Tachogram for HRV Analysis.

HRV metrics are calculated from an analysis of the RR tachogram; the time series of RR intervals. It is important to note that this is an unusual time series in that both axes are time intervals, one being related to the other. Furthermore, since the variability in HR occurs on a beat-to-beat basis, the time series is inherently unevenly spaced along the horizontal axis. Figure 4 illustrates this concept; each star indicates the location of a beat in time (along the horizontal axis). The horizontal distance between each point (time stamp) is different for each adjacent pair, with the difference recorded on the vertical axis. The fact that the RR tachogram is inherently unevenly sampled leads to complications and errors in metrics that utilize interpolation.

Figure 4 Area Calculation by Using Trapezoidal Rule.

HRV analysis in statistical indices of timedomain

Simple time–domain variables that can be calculated include the mean NN interval, the mean heart rate, the difference between the longest and shortest NN interval. For the analysis of RR intervals, we interest primarily in the extraction of beat frequency (F) in beats per minute (beats / min or rpm).

The first parameter is derived from the developed R peaks detection algorithm by the following formula:

RR Distance: R R mean ( S )= L. T samp N QRS MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaamOuai aadkfadaWgaaqcfasaaiaad2gacaWGLbGaamyyaiaad6gaaKqbagqa amaabmaabaGaam4uaaGaayjkaiaawMcaaiabg2da9maalaaabaGaam itaiaac6cacaWGubWaaSbaaKqbGeaacaWGZbGaamyyaiaad2gacaWG WbaajuaGbeaaaeaacaWGobWaaSbaaKqbGeaacaWGrbGaamOuaiaado faaKqbagqaaaaaaaa@4B6E@  (13)

Beat frequency: F( beats/mn )= 60 R R mean ( S ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaamOram aabmaabaGaamOyaiaadwgacaWGHbGaamiDaiaadohacaGGVaGaamyB aiaad6gaaiaawIcacaGLPaaacqGH9aqpdaWcaaqaaiaaiAdacaaIWa aabaGaamOuaiaadkfadaWgaaqcfasaaiaad2gacaWGLbGaamyyaiaa d6gaaKqbagqaamaabmaabaGaam4uaaGaayjkaiaawMcaaaaaaaa@4B49@ or F( beats/mn )= 60. N QRS L. T samp MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaamOram aabmaabaGaamOyaiaadwgacaWGHbGaamiDaiaadohacaGGVaGaamyB aiaad6gaaiaawIcacaGLPaaacqGH9aqpdaWcaaqaaiaaiAdacaaIWa GaaiOlaiaad6eadaWgaaqcfasaaiaadgfacaWGsbGaam4uaaqcfaya baaabaGaamitaiaac6cacaWGubWaaSbaaKqbGeaacaWGZbGaamyyai aad2gacaWGWbaajuaGbeaaaaaaaa@4E8D@  (14)

Where L=total length of the recording (sample number), T samp MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaamivam aaBaaajuaibaGaam4CaiaadggacaWGTbGaamiCaaqcfayabaaaaa@3BFF@ =the sampling period(s), N QRS MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaamOtam aaBaaajuaibaGaamyuaiaadkfacaWGtbaajuaGbeaaaaa@3AB9@ = number of detect QRS.

The root mean square of successive differences (RMSSD) is calculated for the purposes of HRV analysis by

RMSS D x = 1 N t=1 N1 ( x( t )x( t+1 ) ) 2 MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGsbGaamytaiaadofacaWGtbGaamira8aadaWgaaqcfasa a8qacaWG4baajuaGpaqabaWdbiabg2da9maakaaapaqaa8qadaWcaa WdaeaapeGaaGymaaWdaeaapeGaamOtaaaadaGfWbqabKqbG8aabaWd biaadshacqGH9aqpcaaIXaaapaqaa8qacaWGobGaeyOeI0IaaGymaa qcfa4daeaapeGaeyyeIuoaamaabmaapaqaa8qacaWG4bWaaeWaa8aa baWdbiaadshaaiaawIcacaGLPaaacqGHsislcaWG4bWaaeWaa8aaba WdbiaadshacqGHRaWkcaaIXaaacaGLOaGaayzkaaaacaGLOaGaayzk aaWdamaaCaaabeqcfasaa8qacaaIYaaaaaqcfayabaaaaa@5597@ (15)

HRV analysis by using lomb-scargle periodogram in frequency-domain

Since RR intervals in time series are non-stationary and they are spaced unevenly, it needs to make re-sampling. There are different interpolation methods for re-sampling. However, if the time series contains missing samples in heart rate time series, PSD estimates can be severely artificial and in such cases resampling is complicated by the need to infer probable values as replacements. The Lomb-Scargleperiodogram can overcome this weakness and it can estimate the PSD directly from irregularly sampled RR interval series without resampling.
This method is based on the definition of (DWT) as shown in Equation 16.

DWT= n=1 N x( t n ) e jω t n MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGebGaam4vaiaadsfacqGH9aqpdaGfWbqabKqbG8aabaWd biaad6gacqGH9aqpcaaIXaaapaqaa8qacaWGobaajuaGpaqaa8qacq GHris5aaGaamiEamaabmaapaqaa8qacaWG0bWdamaaBaaabaWdbiaa d6gaa8aabeaaa8qacaGLOaGaayzkaaGaamyza8aadaahaaqabKqbGe aapeGaeyOeI0IaamOAaiabeM8a3jaadshajuaGpaWaaSbaaKqbGeaa peGaamOBaaWdaeqaaaaaaaa@4DE4@  (16)
Where, x( t n ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqk0Jf9crFfpeea0xh9v8qiW7rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiEaiaacI cacaWG0bWaaSbaaSqaaiaad6gaaeqaaOGaaiykaaaa@3C81@ is unevenly sampled signal, for(n =1, 2,…N) and ω=2πf MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqk0Jf9crFfpeea0xh9v8qiW7rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyYdCNaey ypa0JaaGOmaiabec8aWjaadAgaaaa@3E40@ is the angular frequency. This equation has been used to define the transform for unevenly sampled series.

The equation for the transformation of unevenly sampled signals to the power spectral density is

P x ( f )= 1 2 σ 2 { [ n=1 N ( x( t n ) x ¯ )cos( 2πf( t n τ ) ) ] n=1 N cos 2 ( 2πf( t n τ ) ) 2 + [ n=1 N ( x( t n ) x ¯ )sin( 2πf( t nτ ) ) ] n=1 N sin 2 ( 2πf( t n τ ) ) 2 } MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGqbWdamaaBaaajuaibaWdbiaadIhaaKqba+aabeaapeWa aeWaa8aabaWdbiaadAgaaiaawIcacaGLPaaacqGH9aqpdaWcaaWdae aapeGaaGymaaWdaeaapeGaaGOmaiabeo8aZ9aadaahaaqabKqbGeaa peGaaGOmaaaaaaqcfa4aaiWaa8aabaWdbmaalaaapaqaa8qadaWada WdaeaapeWaaubmaeqajuaipaqaa8qacaWGUbGaeyypa0JaaGymaaWd aeaapeGaamOtaaqcfa4daeaapeGaeyyeIuoaamaabmaapaqaa8qaca WG4bWaaeWaa8aabaWdbiaadshapaWaaSbaaeaapeGaamOBaaWdaeqa aaWdbiaawIcacaGLPaaacqGHsislceWG4bGbaebaaiaawIcacaGLPa aaciGGJbGaai4BaiaacohadaqadaWdaeaapeGaaGOmaiabec8aWjaa dAgadaqadaWdaeaapeGaamiDa8aadaWgaaqcfasaa8qacaWGUbaaju aGpaqabaWdbiabgkHiTiabes8a0bGaayjkaiaawMcaaaGaayjkaiaa wMcaaaGaay5waiaaw2faaaWdaeaapeWaaubmaeqajuaipaqaa8qaca WGUbGaeyypa0JaaGymaaWdaeaapeGaamOtaaqcfa4daeaapeGaeyye IuoaaiGacogacaGGVbGaai4Ca8aadaahaaqabKqbGeaapeGaaGOmaa aajuaGdaqadaWdaeaapeGaaGOmaiabec8aWjaadAgadaqadaWdaeaa peGaamiDa8aadaWgaaqcfasaa8qacaWGUbaajuaGpaqabaWdbiabgk HiTiabes8a0bGaayjkaiaawMcaaaGaayjkaiaawMcaaaaapaWaaWba aeqajuaibaWdbiaaikdaaaqcfaOaey4kaSYaaSaaa8aabaWdbmaadm aapaqaa8qadaqfWaqabKqbG8aabaWdbiaad6gacqGH9aqpcaaIXaaa paqaa8qacaWGobaajuaGpaqaa8qacqGHris5aaWaaeWaa8aabaWdbi aadIhadaqadaWdaeaapeGaamiDa8aadaWgaaqcfasaa8qacaWGUbaa juaGpaqabaaapeGaayjkaiaawMcaaiabgkHiTiqadIhagaqeaaGaay jkaiaawMcaaiGacohacaGGPbGaaiOBamaabmaapaqaa8qacaaIYaGa eqiWdaNaamOzamaabmaapaqaa8qacaWG0bWdamaaBaaajuaibaWdbi aad6gacqGHsislcqaHepaDaKqba+aabeaaa8qacaGLOaGaayzkaaaa caGLOaGaayzkaaaacaGLBbGaayzxaaaapaqaa8qadaqfWaqabKqbG8 aabaWdbiaad6gacqGH9aqpcaaIXaaapaqaa8qacaWGobaajuaGpaqa a8qacqGHris5aaGaci4CaiaacMgacaGGUbWdamaaCaaabeqcfasaa8 qacaaIYaaaaKqbaoaabmaapaqaa8qacaaIYaGaeqiWdaNaamOzamaa bmaapaqaa8qacaWG0bWdamaaBaaajuaibaWdbiaad6gaaKqba+aabe aapeGaeyOeI0IaeqiXdqhacaGLOaGaayzkaaaacaGLOaGaayzkaaaa a8aadaahaaqabKqbGeaapeGaaGOmaaaaaKqbakaawUhacaGL9baaaa a@BA7D@  (17)

Where x ¯ MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqk0Jf9crFfpeea0xh9v8qiW7rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaWaa0aaaeaaca WG4baaaaaa@3917@  and σ 2 MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqk0Jf9crFfpeea0xh9v8qiW7rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4Wdm3aaW baaSqabeaacaaIYaaaaaaa@3AB5@  are the mean and variance of the series { x( t n ) } MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfa4aaiWaae aacaWG4bWaaeWaaeaacaWG0bWaaSbaaKqbGeaacaWGUbaajuaGbeaa aiaawIcacaGLPaaaaiaawUhacaGL9baaaaa@3E04@ , and τ( f ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaeqiXdq 3aaeWaaeaacaWGMbaacaGLOaGaayzkaaaaaa@3ABD@ is a frequency dependent time delay, defined to make the transform insensitive to time shift, and is computed as

tan( 4πft )= n=1 N sin( 4πf t n ) n=1 N cos( 4 πf t n ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qaciGG0bGaaiyyaiaac6gadaqadaWdaeaapeGaaGinaiabec8a WjaadAgacaWG0baacaGLOaGaayzkaaGaeyypa0ZaaSaaa8aabaWdbm aavadabeqcfaYdaeaapeGaamOBaiabg2da9iaaigdaa8aabaWdbiaa d6eaaKqba+aabaWdbiabggHiLdaaciGGZbGaaiyAaiaac6gadaqada WdaeaapeGaaGinaiabec8aWjaadAgacaWG0bWdamaaBaaajuaibaWd biaad6gaaKqba+aabeaaa8qacaGLOaGaayzkaaaapaqaa8qadaqfWa qabKqbG8aabaWdbiaad6gacqGH9aqpcaaIXaaapaqaa8qacaWGobaa juaGpaqaa8qacqGHris5aaGaci4yaiaac+gacaGGZbWaaeWaa8aaba WdbiaaisdacaGGGcGaeqiWdaNaamOzaiaadshapaWaaSbaaKqbGeaa peGaamOBaaqcfa4daeqaaaWdbiaawIcacaGLPaaaaaaaaa@64EC@  (18)

From this method, the power spectral density can be calculated several frequency bands of interest have been defined in humans.

Area calculation within frequency bands

The trapezoidal rule is a numerical method that approximates the value of a definite integral. We consider the definite integral as

a b f( x )dx MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qadaGfWbqabKqbG8aabaWdbiaadggaa8aabaWdbiaadkgaaKqb a+aabaWdbiabgUIiYdaacaWGMbWaaeWaa8aabaWdbiaadIhaaiaawI cacaGLPaaacaWGKbGaamiEaaaa@415A@ (19)

Assume that f(x) is continuous on [a, b] and divide [a, b] into n subintervals of equal length

Δx= ba n MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacqGHuoarcaWG4bGaeyypa0ZaaSaaa8aabaWdbiaadkgacqGH sislcaWGHbaapaqaa8qacaWGUbaaaaaa@3E0A@                (20)

using the n + 1 point, x 0  =a,  x 1 =a+Δx,   x 2 =a+2Δx, ,  x n =a+nΔx=b MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWG4bWdamaaBaaajuaibaWdbiaaicdacaGGGcaajuaGpaqa baWdbiabg2da9iaadggacaGGSaGaaiiOaiaadIhapaWaaSbaaKqbGe aapeGaaGymaaqcfa4daeqaa8qacqGH9aqpcaWGHbGaey4kaSIaeyiL dqKaamiEaiaacYcacaGGGcGaaiiOaiaadIhapaWaaSbaaKqbGeaape GaaGOmaaqcfa4daeqaa8qacqGH9aqpcaWGHbGaey4kaSIaaGOmaiab gs5aejaadIhacaGGSaGaaiiOaiabgAci8kaacYcacaGGGcGaamiEa8 aadaWgaaqcfasaa8qacaWGUbaajuaGpaqabaWdbiabg2da9iaadgga cqGHRaWkcaWGUbGaeyiLdqKaamiEaiabg2da9iaadkgaaaa@626D@

The value of f(x) can be calculated at these points.

y 0 =f( x 0 ),  y 1 =f( x 1 ),  y 2 =f( x 2 ), ,  y n =f( x n ) MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWG5bWdamaaBaaajuaibaWdbiaaicdaaKqba+aabeaapeGa eyypa0JaamOzamaabmaapaqaa8qacaWG4bWdamaaBaaajuaibaWdbi aaicdaaKqba+aabeaaa8qacaGLOaGaayzkaaGaaiilaiaacckacaWG 5bWdamaaBaaajuaibaWdbiaaigdaaKqba+aabeaapeGaeyypa0Jaam Ozamaabmaapaqaa8qacaWG4bWdamaaBaaajuaibaWdbiaaigdaaKqb a+aabeaaa8qacaGLOaGaayzkaaGaaiilaiaacckacaWG5bWdamaaBa aajuaibaWdbiaaikdaaKqba+aabeaapeGaeyypa0JaamOzamaabmaa paqaa8qacaWG4bWdamaaBaaajuaibaWdbiaaikdaaKqba+aabeaaa8 qacaGLOaGaayzkaaGaaiilaiaacckacqGHMacVcaGGSaGaaiiOaiaa dMhapaWaaSbaaKqbGeaapeGaamOBaaqcfa4daeqaa8qacqGH9aqpca WGMbWaaeWaa8aabaWdbiaadIhapaWaaSbaaKqbGeaapeGaamOBaaqc fa4daeqaaaWdbiaawIcacaGLPaaaaaa@64F3@

We approximate the integral by using n trapezoids formed by using straight line segments between the points

( x i1 , y i1 ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfa4aaeWaae aacaWG4bWaaSbaaKqbGeaacaWGPbGaeyOeI0IaaGymaaqcfayabaGa aiilaiaadMhadaWgaaqcfasaaiaadMgacqGHsislcaaIXaaajuaGbe aaaiaawIcacaGLPaaaaaa@419E@ and ( x i , y i ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfa4aaeWaae aacaWG4bWaaSbaaKqbGeaacaWGPbaajuaGbeaacaGGSaGaamyEamaa BaaajuaibaGaamyAaaqcfayabaaacaGLOaGaayzkaaaaaa@3E4E@ for 1in MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaOaaGymai abgsMiJkaadMgacqGHKjYOcaWGUbaaaa@3C8A@ as shown in the figure below.

The area of a trapezoid is obtained by adding the area of a rectangle and a triangle.

A= y 0 Δx+ 1 2 ( y 1 y 0 )Δx= ( y 0 + y 1 )Δx 2 MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGbbGaeyypa0JaamyEa8aadaWgaaqcfasaa8qacaaIWaaa juaGpaqabaWdbiabgs5aejaadIhacqGHRaWkdaWcaaWdaeaapeGaaG ymaaWdaeaapeGaaGOmaaaadaqadaWdaeaapeGaamyEa8aadaWgaaqc fasaa8qacaaIXaaajuaGpaqabaWdbiabgkHiTiaadMhapaWaaSbaaK qbGeaapeGaaGimaaqcfa4daeqaaaWdbiaawIcacaGLPaaacqGHuoar caWG4bGaeyypa0ZaaSaaa8aabaWdbmaabmaapaqaa8qacaWG5bWdam aaBaaajuaibaWdbiaaicdaaKqba+aabeaapeGaey4kaSIaamyEa8aa daWgaaqcfasaa8qacaaIXaaajuaGpaqabaaapeGaayjkaiaawMcaai abgs5aejaadIhaa8aabaWdbiaaikdaaaaaaa@5794@ (21)

By adding the area of the n trapezoids, we obtain the approximation

a b f( x )dx ( y 0 + y 1 )Δx 2 + ( y 1 + y 2 )Δx 2 + ( y 2 + y 3 )Δx 2 ++ ( y n1 + y n )Δx 2 MathType@MTEF@5@5@+= feaagKart1ev2aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qadaGfWbqabKqbG8aabaWdbiaadggaa8aabaWdbiaadkgaaKqb a+aabaWdbiabgUIiYdaacaWGMbWaaeWaa8aabaWdbiaadIhaaiaawI cacaGLPaaacaWGKbGaamiEaiabgIKi7oaalaaapaqaa8qadaqadaWd aeaapeGaamyEa8aadaWgaaqcfasaa8qacaaIWaaajuaGpaqabaWdbi abgUcaRiaadMhapaWaaSbaaKqbGeaapeGaaGymaaqcfa4daeqaaaWd biaawIcacaGLPaaacqGHuoarcaWG4baapaqaa8qacaaIYaaaaiabgU caRmaalaaapaqaa8qadaqadaWdaeaapeGaamyEa8aadaWgaaqcfasa a8qacaaIXaaajuaGpaqabaWdbiabgUcaRiaadMhapaWaaSbaaKqbGe aapeGaaGOmaaqcfa4daeqaaaWdbiaawIcacaGLPaaacqGHuoarcaWG 4baapaqaa8qacaaIYaaaaiabgUcaRmaalaaapaqaa8qadaqadaWdae aapeGaamyEa8aadaWgaaqcfasaa8qacaaIYaaapaqabaqcfa4dbiab gUcaRiaadMhapaWaaSbaaKqbGeaapeGaaG4maaqcfa4daeqaaaWdbi aawIcacaGLPaaacqGHuoarcaWG4baapaqaa8qacaaIYaaaaiabgUca RiabgAci8kabgUcaRmaalaaapaqaa8qadaqadaWdaeaapeGaamyEa8 aadaWgaaqcfasaa8qacaWGUbGaeyOeI0IaaGymaaqcfa4daeqaa8qa cqGHRaWkcaWG5bWdamaaBaaajuaibaWdbiaad6gaaKqba+aabeaaa8 qacaGLOaGaayzkaaGaeyiLdqKaamiEaaWdaeaapeGaaGOmaaaaaaa@78BC@ (22)

which simplifies to the trapezoidal rule formula.

a b f( x )dx Δx 2 ( y 0 +2 y 1 +2 y 2 ++2 y n1 + y n ) MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qadaGfWbqabKqbG8aabaWdbiaadggaa8aabaWdbiaadkgaaKqb a+aabaWdbiabgUIiYdaacaWGMbWaaeWaa8aabaWdbiaadIhaaiaawI cacaGLPaaacaWGKbGaamiEaiabgIKi7oaalaaapaqaa8qacqGHuoar caWG4baapaqaa8qacaaIYaaaamaabmaapaqaa8qacaWG5bWdamaaBa aajuaibaWdbiaaicdaaKqba+aabeaapeGaey4kaSIaaGOmaiaadMha paWaaSbaaKqbGeaapeGaaGymaaqcfa4daeqaa8qacqGHRaWkcaaIYa GaamyEa8aadaWgaaqcfasaa8qacaaIYaaajuaGpaqabaWdbiabgUca RiabgAci8kabgUcaRiaaikdacaWG5bWdamaaBaaajuaibaWdbiaad6 gacqGHsislcaaIXaaajuaGpaqabaWdbiabgUcaRiaadMhapaWaaSba aKqbGeaapeGaamOBaaqcfa4daeqaaaWdbiaawIcacaGLPaaaaaa@6088@ (23)

Tests and results

ECG de-noising is the initial procedure for heart rate variability analysis. Since the raw ECG loaded from the MIT-BIH database is normally corrupted with different types of noises. Therefore, it is needed to purify before R-peak detection stage. Figure 5 shows original ECG signal with base-line wander drift.

Figure 5 Measuring the Area of a Trapezoid.

Wander baseline drift removal

Wander baseline drift and motion artifact belong to low frequency in which the wander baseline drift frequency is lower than 1 Hz. Since the ECG is a non-stationary signal, normal filters cannot be effective to remove the noise; so, several techniques are used to do so for such types of signals. Baseline wandering is a major source of noise which can be efficiently removed using wavelet based approach. In the original ECG signal, it needs to make de-trending and de-noising. Figure 6 shows the result of base line de-trending.

Figure 6 Original ECG Signal with Base-line Drift of MIT-BIH Record-100.

procedure with DWT

Generally, there are three steps for ECG signal de-noising by using Wavelet transform method. They are decomposition, thresholding detail coefficients and reconstruction. The decomposition of detail levels and the reconstruction of approximation and detail levels are shown in Figures 7-9 shows the base-line drift de-noised signal. The comparison of original and de-noised ECG signal is shown in Figure 10.

Figure 7 Baseline De-trending of Original ECG Signal of MIT-BIH Record-100.

Figure 8 Decomposition of Detail-level for De-noising.

Figure 9 Reconstruction of Approximation and Detail Levels for De-noising.

Figure 10 Base-line Drift De-noising ECG Signal ofMIT-BIH Record-100.

R-peak detection

In wavelet decomposition, the signal is decomposed into shifted and scaled versions of the original mother wavelet. An analysis of signal using ECG includes decomposition of the signal, thresholding of the wavelet coefficients and reconstruction of the signal using modified wavelet coefficients. The wavelet Db6 closely matches with the shape of ECG QRS complex. The ECG signals are decomposed into eight levels using Db6. Addition of third, fourth and fifth detail levels is taken as reconstructed signal as most energy of the QRS complex is concentrated in these coefficients. The peaks in reconstructed signal are observed in the vicinity of ECG R-peaks with maximum tolerance of six samples. All samples in reconstructed signal are scanned for peak. Figure 11 shows the R-waves that find in the restrictions of Wavelet Transform method and Figure 12 shows the detection of R-peak.

Figure 11 Comparison for Original and De-noising Signal of Record-100.

Figure 12 R-wave Detection in ECG Record-100.

Overall performance of R peak detection

The performance of the R peak detection is essential for HRV analysis in detection of arrhythmia. Accuracies of R peak detection can be estimated according to the output results of 48 records. The detected R peaks are compared with MIT/BIH annotation records which represents the actual time intervals from patient. In matching with each annotation files, it is not needed that matching annotations have exactly equal times. While comparing with annotation, the performance of discrete wavelet transform (DWT) algorithm can be checked with the accuracy. From this result, we can see that DWT method for R peak detection is enough for HRV analysis.

Accuracy= TP+TN TP+TN+FP+FN 100% MathType@MTEF@5@5@+= feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaqcfaieaaaaaa aaa8qacaWGbbGaam4yaiaadogacaWG1bGaamOCaiaadggacaWGJbGa amyEaiabg2da9maalaaabaGaamivaiaadcfacqGHRaWkcaWGubGaam OtaaqaaiaadsfacaWGqbGaey4kaSIaamivaiaad6eacqGHRaWkcaWG gbGaamiuaiabgUcaRiaadAeacaWGobaaaiabgEHiQiaaigdacaaIWa GaaGimaiaacwcaaaa@504E@  (24)

Where, True positive (TP) = the number of beats correctly classified as normal.

True negative (TN)= the number of beats correctly classified as abnormal.

False positive (FP)= the number of beats correctly classified as normal when actually abnormal.

False negative (FN)= the number of beats incorrectly classified as abnormal when actually normal.

From the Table 1 the some results for R-peak detection signals are not good such as record 102, 104, 108 and 207 because the characteristics of these signals are abnormal. The features of some of these ECG signals are presented by the Figures 13-15.

Record

TP

TN

FP

FN

Accuracy

100

2272

0

0

1

99%

101

1864

0

0

9

99%

102

442

0

0

162

73%

103

2084

0

0

6

99%

104

1687

0

180

443

73%

105

2602

0

0

89

97%

106

1997

0

10

77

96%

107

2079

0

0

60

97%

108

1185

0

0

638

65%

109

2523

0

0

11

99%

111

2012

0

0

120

94%

112

2536

0

0

13

99%

113

1795

0

0

1

99%

114

1880

0

2

7

99%

115

1953

0

0

7

99%

116

2393

0

3

28

99%

117

1535

0

0

3

99%

118

2281

0

0

19

99%

119

1989

0

12

92

95%

121

1740

0

0

353

83%

122

2476

0

0

2

99%

123

1517

0

0

3

99%

124

1604

0

16

13

98%

200

2533

0

0

258

91%

201

1869

0

0

169

90%

202

2125

0

0

24

99%

203

2829

0

95

183

91%

205

2625

0

26

20

98%

207

1747

0

0

638

73%

208

2876

0

0

163

94%

209

3015

0

0

36

99%

210

2665

0

0

120

94%

212

2751

0

0

11

99%

213

3246

0

0

47

98%

214

2259

0

0

39

98%

215

3362

0

0

37

98%

217

2207

0

0

73

97%

219

2154

0

0

157

93%

220

2048

0

0

20

99%

221

2406

0

0

55

98%

222

2252

0

0

240

90%

223

2603

0

0

39

99%

228

2026

0

0

65

97%

230

2262

0

0

203

92%

231

1141

0

0

9

99%

232

1770

0

0

45

98%

233

3068

0

0

83

97%

234

2752

0

0

11

99%

Table 1 Accuracies of R-peaks detection

Figure 13 R-peak Detection in ECG Record-100.

Figure 14 R-peaks Detection of ECG Record-102.

Figure 15 R-peaks Detection of ECG Record-207.

Results for HRV analysis using time-domain method

Simple time–domain variables that can be calculated include the mean NN interval, the mean heart rate, the difference between the longest and shortest NN interval. The root mean square of successive differences (RMSSD) is calculated. Table 2 shows the results of statistical analysis in time-domain method.

Signal

SDNN
(ms)

SDANN
(ms)

RMSSD
(ms)

NN50
(count)

pNN50
(%)

100

38.7

28.4

55.7

24

6.5

105

74.6

73.9

119.5

43

10.3

106

243.7

200

385

182

55.8

112

16

14.6

18.4

2

0.5

113

242.7

240.5

287.7

222

69.4

114

144.4

130.7

234.4

107

39.2

115

67

66.6

71.6

154

49

116

61.4

55.2

106.1

30

7.6

117

30

27.8

33.7

32

12.9

118

76.4

76.6

106.7

48

13.3

122

38.1

31.2

20

7

1.7

123

117.9

114.8

115

162

65.5

124

71.7

45.3

40.5

53

21.2

201

153.9

153.5

207.7

335

76.3

202

87.5

79.7

143.1

79

30

205

44.2

26.7

48.4

10

2.2

209

49.1

36.9

48.1

52

10.7

212

35.5

33

27

32

6.9

213

27.9

23.2

41.2

36

6.6

217

64.4

61.4

91

75

20.8

219

141.2

112

211.6

309

81.5

220

29

22.1

34.7

10

2.8

221

216.3

214.8

322.8

337

85.1

222

294

267.1

378.5

211

63.9

223

64

57.7

104.1

69

17.1

230

70

52.7

25.6

17

4.3

231

224.9

171

135.6

144

49.5

232

579.2

581.7

843.8

169

59.1

233

141.7

140.5

244.9

349

68.2

234

14.4

10.7

17.3

1

0.2

103

38.5

37.2

29.4

2.4

8.2

107

3104.6

1405.1

5377.4

22

6.3

Table 2 Results of Statistical Analysis in Time-domain Method

RR intervals analysis

The RR interval analysis consists in studying the heart rate variability. It is usually calculated by analyzing the time series of beat to beat intervals from ECG (RR interval) or traces of blood pressure. It is obtained by measuring the time between RR intervals on the electrocardiogram. The values of RR intervals are then plotted versus time, giving a curve called tachogram of HRV or RR tachogram as shown in Figure 16. This tachogram is a combination of sinusoidal waves of different frequencies. From this RR-tachogram, this is an unusual time series in that both axes are time intervals, one being related to the other. Since the variability in HR occurs on a beat to beat basic, the time series is inherently unevenly spaced along the horizontal axis.

Figure 16 R-peaks Detection of ECG Record-108.

Conversion of RR iterval in time-series into PSD

Since RR intervals in time series are non-stationary and they are spaced unevenly, it needs to make re-sampling. There are different interpolation methods for resampling. However, if the time series contains missing samples in heart rate time series, PSD estimates can be severely artificial and in such cases resampling is complicated by the need to infer probable values as replacements. The Lomb-Scargleperiodogram can overcome this weakness and it can estimate the PSD directly from irregularly sampled RR interval series without resampling Figures 17-19 show the result of power spectrum for HRV analysis.

Figure 17 RR interval for a Five minutes segment of ECG Signal of Record-100.

Figure 18 Power Spectrum of ECG Signal of Record-100.

Figure 19 Power Spectrum of ECG Signal of Record-105.

Figure 20 Power Spectrum of ECG Signal of Record-234.

Results for HRV analysis using frequency-domain method

Since RR intervals in time series are non-stationary and they are spaced unevenly, it needs to make re-sampling. There are different interpolation methods for re-sampling. However, if the time series contains missing samples in heart rate time series, PSD estimates can be severely artificial and in such cases resampling is complicated by the need to infer probable values as replacements. The Lomb-Scargleperiodogram can overcome this weakness and it can estimate the PSD directly from irregularly sampled RR interval series without resampling. The results of spectral analysis in frequency-domain method are shown in Table 3.

Signal

Peak
VLF
(Hz)

Peak
LF
(Hz)

Peak
HF
(Hz)

aVLF
(m)

aLF
(m)

aHF
(m)

LF/HF

100

0.03

0.13

0.17

0.001

0.001

0.009

0.099

105

0.03

0.12

0.17

0.001

0.002

0.014

0.144

106

0.04

0.11

0.39

0.006

0.006

0.015

0.392

112

0.02

0.08

0.27

0.014

0.003

0.002

1.334

113

0.02

0.09

0.17

0.004

0.01

0.036

0.278

114

0.01

0.13

0.17

0.001

0.001

0.026

0.052

115

0.02

0.06

0.3

0.009

0.011

0.015

0.705

116

0.03

0.12

0.36

0.001

0.001

0.009

0.074

117

0.02

0.08

0.34

0.007

0.019

0.02

0.926

118

0.01

0.07

0.17

0.009

0.017

0.027

0.619

122

0.02

0.05

0.36

0.013

0.004

0.001

2.975

123

0.04

0.08

0.25

0.013

0.019

0.019

1

124

0.01

0.06

0.36

0.008

0.001

0.002

0.511

201

0.03

0.13

0.17

0.004

0.024

0.044

0.543

202

0.01

0.05

0.17

0.001

0.002

0.014

0.113

205

0.03

0.08

0.28

0.009

0.015

0.028

0.516

209

0.03

0.07

0.17

0.004

0.002

0.006

0.255

212

0.02

0.09

0.17

0.004

0.004

0.008

0.425

213

0.02

0.08

0.24

0.003

0.009

0.023

0.416

217

0.03

0.14

0.17

0.005

0.007

0.02

0.342

219

0.03

0.13

0.39

0.004

0.018

0.044

0.398

220

0.02

0.06

0.3

0.007

0.003

0.008

0.443

221

0.01

0.12

0.18

0.003

0.012

0.032

0.376

222

0.01

0.05

0.39

0.007

0.003

0.014

0.251

223

0.01

0.09

0.29

0.002

0.005

0.013

0.351

230

0.02

0.09

0.19

0.005

0.002

0.001

3.401

231

0.03

0.07

0.37

0.006

0.001

0.001

0.689

232

0.02

0.13

0.33

0.004

0.009

0.041

0.209

233

0.01

0.13

0.39

0

0.001

0.003

0.288

234

0.01

0.05

0.3

0.002

0.001

0.005

0.143

103

0.01

0.05

0.17

0.005

0.002

0.007

0.28

107

0.04

0.15

0.39

0

0.006

0.097

0.067

Table 3 Results of HRV Analysis using Frequency-Domain Method

Discussions on results

Heart rate variability (HRV) is a measure of variations of heart rate between two successive heart beats. It is usually calculated by analyzing the time series of beat-to-beat intervals from ECG signal. The detection of cardiac arrhythmias is a crucial point in the cardiac diseases diagnosis. An arrhythmia is characterized by the irregularity of the heart rate. A heart rate is regular if it is of the order of 60 beats per minute; otherwise; it's called bradycardia or tachycardia. The most commonly used modality for the arrhythmia diagnosis is the ECG. The detection of this cardiac irregularity is based on the R peaks detection and analysis of their regularity (RR intervals).

During the recording of ECG signals, different types of noise from various sources can be superimposed to the original signal. To purify the signal, many methods have been proposed: adaptive filtering, digital filtering and wavelet filtering. In this research, DWT algorithm is used for base-line wander drift de-noising. Using the knowledge that wavelet filtering is efficient and accurate in the compute of the R peaks positions without changing of the shape or position of the original signal, this technique is used to filter the ECG signal. For the choice of the mother function, a comparative study has achieved using various types of mother functions: coiflet, symlet and Daubechies. The best results are obtained by the use of the Db6 mother wavelet and decomposition level-8 is appropriate for de-noising. After de-noising, the purified ECG signal is used as input for R-peak detection. For the R peak detection, DWT method is also used. Firstly, de-noised ECG signal is read and length is calculated. ECG signals are decomposed into eight levels using Db6. After all,, and  detail coefficients are selected and the wave is reconstructed using these coefficients. The peaks in reconstructed signal are observed in the vicinity of ECG R-peak with maximum tolerance of six samples. All samples in reconstructed signals are scanned for peaks.

The variation of heart rate is analyzed by using the linear method such as time-domain (statistical analysis) and frequency-domain (PSD) method. Simple time–domain variables that can be calculated include the mean NN interval, the mean heart rate, the difference between the longest and shortest NN interval. Since RR intervals in time series are non-stationary and they are spaced unevenly, it needs to make re-sampling. There are different interpolation methods for re-sampling. However, if the time series contains missing samples in heart rate time series, PSD estimates can be severely artificial and in such cases resampling is complicated by the need to infer probable values as replacements. The Lomb-Scargleperiodogram can overcome this weakness and it can estimate the PSD directly from irregularly sampled RR interval series without resampling. After all, the results of both methods are verified by normal value of standard measures of HRV. HRV plays an important role in an indication of cardiovascular health since it is considered as an indicator of the activity of autonomic regulation of circulatory function and as the one of the most significant method of analyzing the activity of the autonomic nervous system.

Conclusion

The quantification of HRV has been shown to afford a suggestion of cardiovascular health and can predict the hopeful different arrhythmia conditions such as Myocardial Infarction (MI), heart failure, angina, and sudden death. By applying the Lomb-Scargleperiodogram for HRV analysis, the variation of heart rate can be described more exactly. Moreover, the benefit of frequency-domain method (spectral analysis) is that it can make identification of different arrhythmia. The traditional time-domain method also yields the same results as abnormal. The performance of HRV analysis strongly depends upon the raw ECG signal de-noising and R-peak detection. By using DWT for these processes, the more accurate HRV analysis can be performed. And the HRV analysis also depends upon intra-patient basis (depending upon activity) and on an inter-patient basis (depending on their cardiovascular fitness).

Acknowledgements

None

Conflict of interest

Authors declare there is no conflict of interest.

References

  1. Savita Chandel. A Review on Wavelet Techniques for Different Noises Removal from ECG Signal. International Journal of Advanced Research in Computer Science and Software Engineering. 2016;6(5):558–564.
  2. Laxmi S. Sargar, Manasi M. et al. Automated Detection of R–peaks in Electrocardiogram. International Journal of Scientific & Engineering Research. 2015;6(8):1265–1269.
  3. Joshi AK, Tomar A, Tomar M. Review Paper on Analysis of Electrocardiograph (ECG) Signal for the Detection of Arrhythmia Abnormalities. International Journal of Advanced Research in Electrical, Electronics and Instrumentation Engineering. 2014;3(10):2466–12476.
  4. Subramanian B, Ramasamy A, Rangasamy J. Performance Comparison of Wavelet and Multiwavelet Denoising Methods for an Electrocardiogram Signal. Hindawi Publishing Corporation Journal of Applied Mathematics. 2014:241540.
  5. Sawant C, Patil HT. ECG Signal De–noising using Discrete Wavelet Transform. International Journal of Electronics Communication and Computer Engineering. 2014;5(4):23–28.
  6. Dhir Er JS, Panag Er NK. ECG Analysis and R Peak Detection Using Filters and Wavelet Transform. International Journal of Innovative Research in Computer and Communication Engineering. 2014;2(2):2883–2890.
  7. Ara L, Hossain N, Yahea Mahbub SM. Baseline Drift Removal and De–Noising of the ECG Signal using Wavelet. International Journal of Computer Applications. 2014;95(16):0975–8887.
  8. Transform. International Journal of Computer Applications. 2014;95(16):15–17.
  9. Mohammed Al Mahamdy, Riley HB. Performance Study of Different Denoising Methods for ECG Signals. The 4th International Conference on Current and Future Trends of Information and Communication Technologies in Healthcare (ICTH–2014). Procedia Computer Science. 2014;37:325–332.
  10. Reema S Kalda, Pramod J Deore. ECG Denoising using Wavelet Transform. International Journal of Innovative Research In Electrical, Electronics, Instrumentation And Control Engineering. 2014;2(7);1689–1692.
  11. Janani R. Analysis of Wavelet Families for Baseline Wander Removal in ECG signals. International Journal of Advance Research in Computer Science and Management Studies. 2014;2(2):169–178.
  12. Pandey S, Ayub S. Wavelet Based R Peak Detection in ECG Signals Using Matlab. Journal of Basic and Applied Engineering Research. 2014;1(2):101–110.
  13. Zoltan German–Sallo. Wavelet transform based HRV analysis. The 7th International Conference Inter disciplinarily in Engineering (INTER–ENG 2013). Procedia Technology. 2014;2:105–111.
  14. Tsaneva GG, Tcheshmedjiev J. Denoising of Electrocardiogram Data with Methods of Wavelet Transform. International Conference on Computer Systems and Technologies – CompSysTech’13. Ruse, Bulgaria. 2013;9–16.
  15. Kumar VH, Sudarshan BG, Prasannakumar SC. A Review of Different Analyzing Methods to Extract HRV Data. International Journal of Computer Trends and Technology. 2013;4(4):951–954.
  16. Nouira I, Abdallah AB, Bedoui MH, et al. A Robust R Peak Detection Algorithm Using Wavelet Transform for Heart Rate Variability Studies. International Journal on Electrical Engineering and Informatics. 2013;5(3):270–284.
  17. Hari Mohan Rai, Anurag Trivedi. De–noising of ECG Waveforms based on Multi–resolution Wavelet Transform. International Journal of Computer Applications. 2012;5(18):25–30.
  18. Parihar N, Chouhan VS. Removal of Baseline wander and detection of QRS complex using wavelets. International Journal of Scientific & Engineering Research. 2012;3(4):1–5.
  19. Mukhopadhyay S, Biswas S, Roy AB, et al. Wavelet Based QRS Complex Detection of ECG Signal. International Journal of Engineering Research and Applications (IJERA). 2012;2(3):2361–2365.
  20. Khaing AS, Naing ZM. Quantitative Investigation of Digital Filters in Electrocardiogram with Simulated Noises. International Journal of Information and Electronics Engineering. 2011;1(3):210–214.
  21. https://waset.org/publications/6930/improved–power–spectrum–estimation–for–rr–interval–time–series
  22. Alfaouri M, Khaled Daqrouq. ECG Signal Denoising By Wavelet Transform Thresholding. American Journal of Applied Sciences. 2008;5(3):276–281.
  23. Reed MJ, Robertson CE, Addison PS. Heart Rate Variability Measurements and the Prediction of Ventricular Arrhythmias. Q J Med. 2005;98(2):87–95.  
  24. Moody GB, Mark RG. The Impact of the MIT–BIH Arrhythmia Database. IEEE Eng Med Biol Mag. 2001;20(3):45–50.
Creative Commons Attribution License

©2018 Tun, et al. This is an open access article distributed under the terms of the, which permits unrestricted use, distribution, and build upon your work non-commercially.