1. Introduction

Over the last few decades, revolutionary technological progress provided significant breakthroughs to diminish restrictions on health state monitoring. New techniques and various approaches have been adapted with still-ongoing modifications and researches. Such advancements allowed developing new devices and computer-aided tools that added remarkable evolution in clinical cardiology and also enabled exploiting the opportunity to monitor and analyze ECG signals on a beat-to-beat basis.

This present effort offers an extensive body of information concerning the techniques most utilized in the field, thus granting the ability of extracting ECG features, which in turn enables the determination of the Heart Rate Variability (HRV). The latter is the milestone that permits dealing with various physiological disorders and also offers physiologically good solutions and control mechanisms. Examples include drug dosage control mechanisms such as drugs that affect directly the stability of the autonomic nervous system which in turn alters whole body homeostasis. Ultimately, this study aims at providing a state-of-the-art update on QRS detection algorithms and heart rate variability analysis techniques based on different medical and engineering research perspectives.

In terms of organization, the paper is structured as shown in Figure 1, which presents the general block diagram that demonstrates the main components used for HRV evaluation. As shown in the Figure, the block diagram starts with the way the data acquisition process is performed, which could be easily carried out by considering MIT-BIH database as a reference for ECG signals with different specifications. The reason behind using this particular database is that real-time acquired signals by hardware are significantly affected by noise, which will cause difficulties with validating accuracy of the results. Therefore, this database simplifies preprocessing of entered signals. The second section in the block diagram is the preprocessing stage, which consists mainly of signal filtration used to remove the noise and any existing artifact. The third block is used to suppress the waves comprising the ECG signal except for the QRS complex while the following subsequent stage (the fourth stage) deals with labeling the R-peaks with their time of occurrence. It should be noted here that the first four sections are combined altogether to emphasize on QRS detection techniques. In contrast, the 5^th stage forms the time series describing the R-R interval variations with respect to time. The significance of extracting the R-R interval is to determine the length of heart beat intervals, since HRV is the physiological phenomenon at which the time interval between each heart beat varies inconsistently and thus is evaluated by measuring the interbeat time interval fluctuations. Consequently, a reliable and accurate QRS detection algorithm must be exploited to obtain trustworthy HRV data. Finally, the 6^th component is related to the employment of this time series in different mathematical and statistical approaches to assess the HRV. Note that the 5^th and the 6^th steps combine to highlight the HRV evaluation processes.

2. Physiological Origins

2.1. ECG Morphology

ECG signals can be resolved into heartbeats, each of which represents, from a physiological point of view, a cardiac cycle. A heartbeat is segmented into smaller functional standard waves, designated as P, Q, R, S and T. Generally, these waves declare depolarization and repolarization phases of the heart muscle. A heartbeat also contains interwave time segments, specifically, PR interval, QRS complex, QT interval and ST segment in addition to the refractory period. The latter extends from the end of a T wave of the preceding heartbeat to the beginning of next heartbeat P wave. This period guarantees enough time for heart muscle to relax. When dealing with multiple heartbeats, the R-R interval, which is the time taken between each ventricular blood ejection, or the time between each two consecutive heartbeats, is a significantly important feature. These waves and time segments are illustrated in Figure 2 while the values of such morphological features in terms of their amplitude and time for both normal and abnormal conditions are listed in Table 1 which also gives the physiological definition of each feature ^{[1, 2, 3]}.

2.2. HRV Definition

HRV, which has been for long used as a screening tool for diagnostic purposes, is the physiological phenomenon that reveals the state cases of having non-consistent R-R durations/intervals over a number of cardiac cycles per unit time. The origin of heart beating is the electric initiation caused by the SA node that generates about 100-120 pulses per minute during rest. However, when subject is in equilibrium, Heart Rate (HR) will never reach this rate because of the continuous non-voluntary control provided by the Autonomic Nervous System (ANS). The ANS has two main branches, namely, the Sympathetic Nervous System (SNS) that triggers and escalates body responses and the Parasympathetic Nervous System (PNS), which contradicts or reverses actions caused by the sympathetic ANS in order to reach an optimized stable state.

Sympathetic and parasympathetic systems control all body organs. This review inspects the significant impact of such systems on HR, which basically reveals the HRV contribution to body homeostasis.

A wide variety of HRV analysis methodologies are available including time domain, statistical methods, geometric methods, frequency domain, and power spectral density analysis ^{[4, 5]}. When compared to other methods such as Poincare plot, frequency or spectral analysis, geometric analysis and non-linear analysis, time-domain analysis is characterized by the simplicity of calculating statistical parameters. The bottom line is that HRV time-domain analysis methods have higher procedural simplicity to perform than other assessment methods. Nevertheless, it is important to state that this method is not informative enough when implemented by itself as will be shown in the following sections.

It is established that HRV describes the ability of the heart to adapt and respond to changing circumstances of different types and cases of stimulations and is predominantly dependent on the extrinsic regulation of HR. The aim of HRV analysis is to evaluate the health conditions and the state of ANS that is mainly in charge of cardiac activity regulation. Moreover, HRV is in general a numeric description of HR fluctuation around a baseline that represents an average heart rate of a subject. Consequently, when an average HR is 60 bpm, R peaks are not equally spaced from each other with exactly one second rather R-R intervals may fluctuate from 0.5 second and extend up to 2 seconds. This fluctuation is due to the fact that many physiologically related aspects affect the change of heart rate variability. Figure 3 presents a visual illustration of the fluctuations in HRV per unit sample.

2.3. Physiological Relevance of HRV

HRV is considered a powerful tool to judge the predominant impact of ANS on HR, which in turn affects body homeostasis. Some disorders reflect major fluctuations on HRV and thus, body systems gather to revert such actions back into balance.

Literature surveys indicate the presence of a set of influential factors that could significantly affect HRV. These factors are age [6-11]^[6], gender ^{[10, 11, 12, 13]}, heart disease ^[14], neurological disease ^{[15, 16, 17]} and exercise ^{[11, 14, 18, 19, 20, 21]}.

2.3.1. HRV and Autonomic Nervous System

As previously mentioned, the ANS consists of both SNS and PNS, both of which work in correlation with each other to fulfill the task of controlling the cardiac activity ^[22]. Sympathetic stimulation occurs as a response to stress, exercise and heart disease. It causes an increase in HR by increasing the firing rate of pacemaker cells in the SA node. Parasympathetic activity, resulting primarily from the function of internal organs, trauma, allergic reactions and the inhalation of irritants and anesthetics, decreases the firing rate of pacemaker cells and the HR creating a regulatory balance in physiological autonomic function. Separate contributions of both the SNS and PNS integrate the ability of modulating heart rate and regulating QRS intervals at distinct frequencies. Sympathetic activity, with respect to HR modulation frequencies, falls within the range of low frequencies (0.04–0.15 Hz) while parasympathetic activity lies within the higher frequency range (0.15–0.4 Hz). This difference in frequency ranges allows HRV analysis to distinguish between sympathetic and parasympathetic participation events.

Studies dealing with this topic proved that when PNS activity dominates the activity of SNS, heart rate will decrease and HRV will increase. The opposite also holds true ^{[23, 24]}.

2.3.2. HRV and Blood Pressure

Blood Pressure (BP) alterations may increase the risks of having abnormalities in the cardiovascular system, especially in hypertensive individuals. Electrocardiographic evidences show that HRV is significantly decreased in cases with cardiovascular abnormalities such as Left Ventricular Hypertrophy (LVH), aortic valve disease and hypertension. A strong correlation was found between HRV and arterial baroreflex, which is reduced in the case of diabetes and hypertension ^[25].

2.3.3. HRV and Myocardial Infarction

Dominance of sympathetic activity with respect to reduction in parasympathetic cardiac control has been found in patients with acute Myocardial Infarction (MI). Sympathetic activity dominance over the parasympathetic activity reduces thresholds of fibrillations and shows a significant tendency for Ventricular Fibrillation (VF). Vagal activity, known also as the parasympathetic activity, increases thresholds of fibrillations and seems to provide protection mechanism against malignant cardiac tachyarrhythmia. It was also shown that HRV decreases with the presence of MI and that extra training and exercise does not alter or improve HRV after having MI. Left ventricular dysfunction is diagnosed on the basis of having reduced numeric values of the standard deviation between each R-R interval (SDRR). This approach is obtained by continuous observations of the imbalance and dominance of parasympathetic activity over sympathetic activity after the occurrence of MI ^{[26, 27, 28, 29]}.

2.3.4. HRV and Diabetes

Diabetes can cause severe autonomic dysfunction. Decreased beat-to-beat variability causes lowered HRV as a complication of poorly controlled diabetes. It was also concluded that cardiac PNS activity is diminished in diabetic patients ^{[30, 31, 32]}.

2.3.5. HRV and Renal Failure

Many tests were run on patients with renal failure. Five-minute HRV analysis experiments on patients with chronic renal failure showed a negative correlation between the mean of R-R intervals and calcium ions concentration after hemodialysis. In renal failure cases, HR power spectral density analysis exhibited a strong reduction of power spectral density in all frequency ranges, which dictates reduction of HRV ^{[33, 34, 35, 36]}.

2.3.6. HRV and Gender, and HRV and Age

It is a proven fact that gender and age are also parameters that affect HRV. Emese et al. ^[37] conducted a study the results of which proved that in newborns, HRV was lower for boys than girls. The HR variation for healthy subjects from 20 to 70 years was studied by Bonnemeir et al. ^[38]. Authors found that HRV decreases with age and that HRV in females is higher than in males. Previous studies revealed that day and night hours have different HRVs when evaluated using time and frequency domain analysis methods. Amounts of HRV are affected by maturity of sympathetic and parasympathetic (vagal) nerves. Maturation of these nerves in the early stages of a person life increases HRV while as age increases HRV decreases ^[39].

2.3.7. HRV and Drugs/Medications

HRV is significantly influenced by administering drugs and medications and as such, HRV helps to quantify the impact of certain drugs on ANS. Since some drugs and anesthetics are aimed to suppress SNS or PNS activity, there will be significant changes in HRV response, which assists in quantifying the amounts of drugs that should be administered. In addition, it provides positive predictions of health disorders and prevents life threatening events ^{[40, 41, 42]}.

2.3.8. HRV and Smoking

Studies using HRV analysis have shown that smoking increases sympathetic and reduces vagal activity. Studies also suggested that smoking reduces HRV ^{[43, 44]} and further indicated that HRV decrements are highly related to exposure to tobacco outlets or smoking since this affects the control activity of ANS ^{[45, 46]}.

3. ECG Feature Extraction and Cassification – State of the Art

ECG parameter extraction has been investigated for several decades. Lots of methodologies as well as transformations have been proposed for accurate and fast ECG feature extraction. Various techniques related to this area of research are reviewed and discussed in the following sections.

Feature extraction approaches may be classified into the following categories:

•  Time-domain thresholding “statistical approach”.

•  Spectral analysis.

•  Geometry analysis.

•  Principal component analysis.

•  Fuzzy logic systems.

•  Artificial neural networks.

•  Neuro-Fuzzy networks “hybrid systems”.

In this section, some of the QRS detection methods as categorized above, are briefly discussed and compared. Simplicity, accuracy, predictivity, consistency, mathematical complexity, and hardware demands are the main characteristics of concern in the comparison performed herein.

3.1. Time-Domain Thresholding
3.1.1. Amplitude Threshold

In an attempt to horizontally cut the ECG signal, early algorithms used amplitude thresholds prior to the step of differentiation. This practice was preferred mainly to prevent the potent undesired impacts of P and T waves over R waves. The first derivative was then applied to increase the steepness of the QRS complex slope. The amplitude threshold was computed according to certain criteria as a fraction of the measured ECG signal.

The robustness of this technique to noise was judged by relative deficiency due to a set of defined technical downsides. This is attributed to the insufficient attenuation of the noise signal. In addition, the parameter selection requires the attainment of ECG data segments of fixed lengths that are exactly the same as in the first-derivative technique (discussed below). Moreover, in this algorithm, these processed segments are commonly determined empirically. In contrast, in order to obtain the value of the amplitude threshold, it could be adjusted only once prior to analysis stages of ECG signals. Hence, the value of the threshold is held constant through the entire ECG analysis. Friesen et al. ^[47] employed ECG data with fixed segment lengths of 33 seconds. Although this algorithm provided relatively high accuracy, it is strongly believed that the choice of small segments helped obtain such accuracies. Accordingly, the performance of this algorithm on longer ECG signals is expected to be poor unless these longer signals are broken down into smaller segments. Furthermore, there is still a possibility of data loss in the form of missing beats at the beginning and at the end of each processed ECG segment ^[48].

3.1.2. First Derivative

The algebraic derivative, which was first developed by ^[49], has proved its robustness to noise in image processing ^[50], signal processing as well as in biomedical engineering. Most commonly, the differentiator is represented by a high-pass filter. Benefits provided by applying a differentiator to ECG signals are modifying signal phase, and generating zero crossings at R peak locations. Pan and Tompkins ^[51] are considered pioneers in developing the technique that utilized derivation as one of its prime stages.

Literature exhibits diverse and extensive usage of the first derivative to detect QRS complexes in numerous formulae. For example, it has been used in ^{[52, 53, 54]} to derive ECG signals followed by a threshold while in ^{[55, 56]}, an amplitude threshold was applied to ECG signals followed by a first derivative and then by a threshold. Also, a combination of first and second derivatives (presented in the next section) followed by thresholding was implemented in ^{[57, 58]}. First derivative that precedes digital filtration was found in ^[59] as well as in digital filters applied to ECG signals followed by a first derivative and subsequently by a threshold ^[51]. Mathematical morphology filtering applied to ECG signals followed by a first derivative tailed by a threshold was also employed ^{[60, 61]}. The literature also reports using the first derivative prior to applying either Hilbert transform as in ^{[62, 63, 64]} or Wavelet transform as in ^[65], both tailed by thresholding.

Technically, comparative studies show that the first derivative provides remarkable assistance in reducing motion artifacts and baseline drifts, but failed to reduce high-frequency noise ^[65]. In terms of parameter selection, researchers have previously utilized first derivative in various combinations without declaring the reasons behind their choices. In addition, the processed ECG segments all had equally fixed lengths and thresholds ^{[52, 53, 54]}. In an attempt to assess the numerical efficiency of this technique, it is noteworthy that amplitude and first derivative class of algorithms is simple and encloses a single equation for feature extraction. In this regard, the majority of cases encountered the Okada’s equation. However, in this class of algorithms, complexity becomes apparent if segmentation takes place. Therefore, the order of complexity depends primarily on the count of processed segments per record.

3.1.3. Combined First and Second Derivatives

Literature reports numerous research efforts that have employed a combined first-and-second-derivative approach in an attempt to enhance the QRS complex detection. Although this technique is not fairly justified, it is believed that the second derivative contributes to emphasizing higher frequencies of the QRS complex but unfortunately including those due to noise, which requires substantially more significant smoothing. This approach involves simply linear combination of the magnitudes of the first and second derivative and is utilized to highlight the QRS area with respect to the rest of the ECG components. Due to its relatively low computational complexity, the technique was extensively employed in the early development of hardware and software approaches ^{[66, 67, 68]}. Therefore, the approach is still in active use in real or quasi-real time analysis as well as in inspecting long-term recordings. The second derivative technique could be used before applying Hilbert transform, followed by a threshold ^[69]. Also, it can be combined with the first derivative of ECG signal followed by a threshold ^[70].

From a technical perspective, this method is characterized by simplicity since it encloses up to four equations dedicated to feature extraction that have, however, to be experimentally conducted ^{[71, 72]}. In the process, all selected parameters must be held fixed through the entire analysis procedure. In addition, like all derivation-based techniques, this approach operates well when ECG data segmentation takes place. Consequently, as the number of segments increases, complexity increases.

3.1.4. Matched Filters

A robust ECG feature extraction algorithm that introduces the possible benefits of using a matched filter was proposed in ^{[72, 73]}. This method is roughly built on signal-to-noise ratio (SNR) maximization basis by using this type of filters. The experimental data in matched filters are obtained from the MIT-BIH database and the signals were acquired from four persons of different ages. One of them had normal sinus rhythm while each of the remaining three had cardiac abnormality of various degrees. All acquired signals were re-sampled at a sampling frequency of 500 Hz. A peak detection threshold was set at two-thirds of the averaged output of the matched filter for a specific feature of interest. On the one hand, ST segment extraction was hard to perform because of amplitude variability and noise and, on the other hand, the detection of other features was harder than expected not because of the implementation of the matched filter, but due to difficulty in discriminating the R-peak from the ST-peak. Consequently, the need arises to set a modified threshold detection algorithm and to use a different form of adaptive matched filters. The best results using this method were obtained only when normal cases were considered.

In linear systems, the matched filters demonstrate themselves as optimized linear filters for maximizing the SNR in the presence of additive stochastic noise. Technically, a matched filter aligns predefined signal properties like shape or magnitude altogether. Matched filters are commonly used in radar applications, where a signal is sent out and the reflected signals are measured and examined for something in close similarity to the sent out signal. In medical sciences, two-dimensional matched filters are commonly utilized in image processing to improve the SNR in x-ray images. Nonetheless, an ECG signal is a nonlinear signal produced by a nonlinear system, i.e., the human body. Thus, it is problematic and relatively difficult when trying to adapt a nonlinear signal through a linear model. Simplistically, this explains the debates over claiming that the performance of the matched filters is not optimized in this field of research.

Literature indicates the existence of a set of QRS complex detectors that operate well in the presence of moderate noise ^[51]. These detectors usually employ a band-pass filter with a center frequency in the range of 10-17 Hz. The signal passes through multiple processing stages to estimate the local energy in the defined pass band. However, these techniques typically face two main technical obstacles: (1) the signal pass band of the QRS complex is variant not only for different subjects, but also for different beats of the same subject, and (2) pass bands of both the QRS complex and the noise keep overlapping. In cases of having a familiar signal in the noise to be detected, a matched filter could maximize the SNR. Otherwise, the operating principle of this method declares dealing with both desired and undesired quantities of signals. Therefore, the design of this filter has to be optimized to minimize the effects of noise at the filter output in some statistical sense to enhance the detection of the desired signal. Subsequently, to optimize the filter design, the identification of both the signal and correlation statistics of the noise becomes inevitable. It may be added here that the non-stationary properties of ECG signals and noises represent a major hindrance in the application of matched filters to detect QRS complexes. A remedial action for detecting QRS complexes under the impact of extremely noisy environments that entailed the development of a linear adaptive matched filter was developed in ^[74]. Although matched filters in ^[51] are linear adaptive filters, they intend to adapt to compensate for changes in signal shapes and noise corruption. In addition, this filter demonstrated better performance than local energy estimation method cited in ^[51], particularly, in the presence of motion artifacts.

Real-time computations of matched filters are reported in ^{[75, 76]}. The literature reports on the application of matched filters in QRS detection as they are independently applied to the ECG signal ^[77], follow digital filters ^{[78, 79]}, or follow neural networks ^[80]. Regarding the extent of immunity to noise, the matched filters proved to offer remarkable improvement by increasing the SNR ^[81]. Nevertheless, the method seems to be characterized by mathematical complexity, particularly as related to parameter selection such as the fixed template length and the filter design. In terms of numerical efficiency, the method is computationally time-consuming because of the sample-by-sample moving comparison with the template along the ECG signals.

3.1.5. Digital Filters

A simple digital filter with cut-off frequencies at extrema of the to-be-detected ECG feature would attenuate other wave components and some artifacts ^[82]. In addition, digital-differentiation-based techniques explicate the circumstance that QRS complexes enclose the maximal slope change in samples. Most of the real-time QRS detectors use a differentiator at the initial stages of preprocessing. More specifically, in an attempt to suppress the sensitivity to noise, a band-pass filter is usually optimized and is fundamentally based on high- and low- pass filters cascaded and designed accordingly.

Algorithms based on more refined designs of digital filters were proposed in literature. Mehta and Lingaya ^[83] employed digital filters to reduce the power-line interference and drifts, with satisfactory success. Two years later, the same authors developed a digital differentiator preceding a Support Vector Machine (SVM) classifier to depict QRS complexes and scored an accuracy of 99.75% ^[84]. In ^[52], the proposed algorithm involved ECG filtering in parallel using two different low-pass filters designed with two distinct cut-off frequencies. The difference between the filter outputs represents the band pass filtered ECG signal. This tactic is considered a nonlinear operation that regularly issues an effective suppression of small values and grants a slight smoothing of the peaks. The principles of digital filtering in conjunction with wavelet transform and dyadic wavelet transform basics were also utilized ^[85]. The design situated the digital filters before R peaks detection stage. The filters were constructed using a low-pass filter based on discrete wavelet transform, and a high-pass filter based on the dyadic wavelet transform. This algorithm recorded a percentage of detected R peaks of 99.26%. Remarkable evolution in this path has been pioneered by Pan and Tompkins (PT) ^[51]. They introduced a processor based on real-time QRS classification mechanism that utilized digital analysis of slope, amplitude, and width of QRS complexes after preprocessing. Thereafter, it was modified by Hamilton and Tompkins ^[86] that added an adaptive threshold to the same algorithm. Later, Savitsky-Golay Smoothing Filters (SGSF) came across to replenish the original high-pass filter and the differentiator with a SGSF and a subtractor, which significantly improved the efficiency of ECG processing ^{[87, 88]}. Das and Chakraborty ^[90] applied PT by using SGSF modification and a multiplication of backward difference (MOBD) algorithm, which is basically an AND-combination of adjacent values of the derivative, was projected as digital filters to detect QRS complexes ^[91]. Digital filters have been used in various combinations that used a band-pass filter followed by matching filter followed by threshold ^[92]. Moreover, they could be combined with derivative-based techniques ^[93], applied before performing Hilbert transform ^[51], or before applying Wavelet transform ^[65].

Digital filters seem to show satisfactory robustness to noise due to their ability to uplift the SNR according to the nature of the implemented filter and its order. With respect to parameter selection, the processed segments have equal and fixed lengths as well as all selected parameters. Besides, selecting a digital differentiator grants operation of a notch filter. The mathematical operations in digital filters class of algorithms like squaring, difference, and multiplication are not generalized as must-perform procedures. This means that they are set and justified by researchers according to the foreseen output of the algorithm.

3.1.6. Hidden Markov Model

Hidden Markov Models (HMMs) have proved their potential and flexibility as an effective class of statistical models for describing sequential data sets ^[94]. In ^[95], Coast and Cano investigated the application of HMMs to QRS and ECG waveform detection and were the first to demonstrate the usefulness of HMM as a promising technique. The main objective of this method is to deduce the observed data sequence by a probability function that differs according to the state of the underlying Markov chain ^[96]. In the case of ECG signals, potent states are P-wave, QRS complex, and T-wave. The added value of this technique is that it detects not only QRS complexes, but also the P and T waves.

Studies in ^{[95, 96]} investigate implementing HMM following the application of band-pass filtering to ECG signals. In ^[98], the authors applied wavelet to ECG signals prior to applying HMM. Problems encountered in this method comprise necessary manual segmentation for training prior to the analysis of the record. This manual segmentation includes determining the number of states, transition probabilities and output functions. In addition, the underlying (hidden) state sequence of the produced data is unknown, and hence, the parameters of HMM cannot be simply obtained by applying likelihood estimation formulae ^[95]. In addition, HMM method requires large numbers of parameters to be evaluated. Specifically, 15 to 50 parameters all of which are to be held fixed through the entire procedure. Nevertheless, HMM is susceptible to noise, sensitive to baseline wander, DC shift, and HRV.

3.1.7. Pan and Tompkins

Pan and Tompkins (PT), also known as the low-pass differentiation algorithm (LPD), introduced a major evolution in ECG signal processing in 1985 ^{[51, 99]}. Slope, amplitude and width information were used to detect QRS complexes in PT algorithm. The QRS detection was implemented through three detection steps: linear digital filtering, non-linear transformations, and decision rule algorithms.

The PT method does not consume significant power rates. In the process, a raw ECG signal is passed into an A/D converter at a sampling frequency of about 200 Hz after filtration by an analog band-pass filter to limit the band of the ECG signal at about 50 Hz before digitization. Pattern recognition tasks involve a preprocessing step, which consists of a band-pass filter cascaded into a low-pass and a high-pass filter configuration. The low-pass filter is used to limit the operating range of an ECG signal and also to reduce higher frequency noise effects while the high-pass filter is used to highlight the onset of each QRS complex.

A digital filter design having integer coefficients allows for real time processing. The overall band-pass filter significantly reduces all kinds of undesired interferences and frequency noise impacts. The resultant signal is then passed through a local peak detection algorithm, which identifies and marks all R peaks of the ECG signal. This algorithm uses a set of thresholds to select the required QRS complexes. The threshold adaptation is based on the amplitude of the identified peak.

The output of the band-pass filter constitutes the input for a differentiation element. The latter is used to find the slope of the QRS complex. Initially, a five-point derivative operator is performed. If detection is carried out properly, the output is moved on to the next stage. If there is a missing beat or R peak, a two-point derivative operator is performed, and the highest peak of them is considered as an R peak. Generally speaking, the detection of such features requires a precise knowledge of locations of the boundaries “onset and offset” of each interval. The most important task of the decision rule is to determine the threshold. Problems that might be faced is a missing beat, noise interference, non-linearity, non-stationary variations and human source of errors while setting the threshold database.

The following squaring process intensifies the slope of the frequency response to help to detect false peaks such as those related to T-waves. A point-by-point squaring function is applied to the output of the differentiator to guarantee that all parts of the signal are within the positive polarity to emphasize the higher frequencies present due to the QRS complex according to Equation (1)

(1)

This represents the non-linear processing part of the signal.

Subsequently, a moving window integrator obtains information about the width of the QRS complex. Time averaging is performed by adding together the 32 most recent values from the squaring function and dividing the total by 32. This result is then passed through the same local peak detection and threshold setting algorithms as the original band-pass filtered signal to identify QRS slope information. All potential QRS peaks found in both filtered and transformed waveforms are then compared. Only those appearing in both processed waveforms are classified to be true QRS complexes. The output is a stream of pulses indicating the locations of the QRS complexes ^[100]. Such an algorithm relies not only on slope information, but also on the QRS complex amplitude and width information. In the Pan Tompkins algorithm, human factor plays a role in the choice of thresholds. Operator experience in thresholds settings might be the main source of error in this method. Table 2 summarizes the performance of the algorithm compared to matched filters and indicates that PT is superior. Consequently, the focus in the next section will be on this algorithm. The discussion will briefly touch on some of the proposed PT modifications done through the past few years.

3.1.8. Hamilton and Tompkins

Patrick S. Hamilton and Willis J. Tompkins developed a modification for the QRS detection by adding new decision rules ^[86]. In this method, signals are filtered to yield a set of QRS complex and noise events. Decision rules are then applied to identify the QRS complex events from noisy ones. Here, determining the threshold is the most important task.

At the pre-processing stage, ECG signals are filtered using low-pass and high-pass filters and the peak is detected by the time averaged window applied to the QRS complex. The authors in ^[86] developed an algorithm that eliminates detection of both ripples on large waves as well as very small noise peaks. The decision rule section uses an eight-point median peak level estimator, a detection threshold coefficient of 0.1825, a search back threshold that is 50% of the normal threshold, an eight-point median R-R interval estimator, and a 200 ms refractory blanking period. The accuracy of the system turned out to be 99.77%.

3.1.9. Savitsky-Golay Smoothing Filters

Savitsky-Golay low-pass smoothing filters (SGFs), which have been used for ECG noise reduction and compression, represent the first modification to PT ^{[87, 89, 101]}. Such filters are also named least-squares or Digital Smoothing Polynomial (DISPO) filters. The simplest of SGFs are the ones with a zero order and a predefined frame length. The special characterisitc feature illustrating the significance of SGFs utilization resides in the ability to deal with any present data in the time-domain directly, which allows for faster real-time data processing.

The SGF is mainly based on least-square, multiple line regression and more complex numerical methods. The basic idea is that for some fixed nL = nR, each g_i is computed as the average of the data points from f_i - n_L to f_i + n_R with the constant C_n = 1/(n_L +n_R +1). The core of Savitsky-Golay filtering is to find filter coefficients C_n such that

(2)

where g_i is a linear combination, f_i represents data values, i is the index of data points, n_L is the number of data to the left, and n_R is the number of data to the right.

The benefits of adding SGFs to PT lie in the replacement of the high-pass and the differentiator components by the Savitsky-Golay smoothing LPF and a subtractor. Still, however, high mathematical and numerical complexity could be faced while finding the filter correlated coefficients, an issue that could be resolved by a software tool.

A comparison between PT and SGFs as related to their application in ECG signal processing is illustrated in Figure 4 through Figure 6. The Figures show that the SGF is superior to PT in terms of the quality of setting the thresholds resulting in a better R peak amplitude detection. The general block diagrams of PT and SGFs are given in Figure 7 and Figure 8.

3.1.10. Mathematical Morphology

ECG signal analysis processes is usually not a trivial task due to the presence of bulky ECG-ridden noises such as impulsive noises caused by muscle contraction, power-line interference, and baseline drift. Mathematical morphology has been adopted mainly because it provides the eligibility of removing these kinds of noise sources using morphology-based filtering. In principle, mathematical morphology is a prevailing tool for quantitative interpretation and analysis of geometrical structures using nonlinear signal operators. This apparatus aims at extracting data related to the shape and size of geometrical objects. The mathematical morphology technique was first emerged from image processing, and then it was considered as a dedicated ECG feature extraction technique ^[102].

The mathematical morphology filter utilizes two basic morphological operators which are erosion and dilation. To enhance the algorithm immunity to noise, a low-pass filter composed of the accumulation modulus is applied. Finally, peak detection is performed using adaptive thresholding and a set of decision rules ^[103]. Literature reports on the application of mathematical morphology to ECG signals, followed by the first derivative, then by a threshold ^[104]. However, opening and closing were found to outperform the multiscale dilation and erosion filtering as extended morphological operators, seemingly because they subscribe to space independence. These two operators could work as morphology filters enhanced with clipping properties to depict an ECG signal valleys and peaks from its within.

The mathematical morphology approach has attracted remarkable attention due to its rich theoretical framework, low computational complexity, and simple hardware implementation. Moreover, prior awareness of the frequency spectrum is no longer necessary ^[105]. That is, this approach enables avoiding frequency band overlapping of QRS complexes and other signal components such as P and T waves. Nevertheless, the calculations of erosion and dilation operators for ECG signal are time-consuming ^[106].

3.2. Spectral Analysis
3.2.1. Fast Fourier Transform

The idea behind using fast Fourier, discrete wavelet and Hilbert-Huang transforms in signal processing is to enhance accuracy. Although many previous investigations proved that Fast Fourier Transform (FFT) is not the appropriate technique when processing non-stationary and nonlinear signals, FFT has attracted many researchers over the last decades due to its ability to break down the signal into many simpler components, i.e., a mother frequency and its harmonics as stated in Equation (3) below. The main obstacle to this technique is that it resolves the signal along its time domain into an infinite frequency span. This span is represented with a complex exponential-enveloped sinusoidal signal. As a result, the transformed signal frequency range extends from -∞ up to +∞, which prevents finding a certain frequency content in a finite window ^{[107, 108, 109]}.

(3)

3.2.2. Wavelet Transform

The infinite range of the Fourier integral leads to time-averaged analysis results making it incapable of determining specific features within the signal. This technical obstacle could be alleviated by adding a sliding time window of fixed length to confine the analysis in time. Therefore, a number of time-frequency signal analysis techniques presently exist that provide the ability of high resolution in the time-frequency space. Such techniques include short time Fourrier transform (STFT), Winger-Ville transform (WVT), Choi-Williams Distribution (CWD), and wavelet transform.

Among these, the wavelet transform has recently emerged as one of the most dominant tools for analyzing challenging signals across a variety of areas in engineering and medicine ^[110]. The wavelet transform is mainly classified into two separate varieties: Continuous Wavelet Transform (CWT), and Discrete Wavelet Transform (DWT).

The CWT is a time-frequency analysis method that grants a time-scale representation similar to the time-frequency representation of the STFT ^[65]. However, the CWT performs signal analysis using a set of functions permitting a variable time and frequency resolution for various frequency bands. This is achieved as CWT implements a variable width window, which allows for proper isolation of the high frequency features contained within signals. Another distinction from the STFT is that the CWT is not constricted to sinusoidal series and analyzing functions. Instead, a large selection of localized wavelets that satisfy mathematical criteria is used. The wavelet transform of a continuous time signal, is defined as:

(4)

where ^* is the complex conjugate of the analyzing wavelet function, a is the dilation parameter of the wavelet, and b is the location parameter of the wavelet.

The DWT, however, is about choosing an integer power of two scaling in a and b. The selection of these parameters leads to the dyadic DWT (dyDWT). The transform integral remains continuous for the DWT, but is identified only on discretized values of a and b. The sampling of these parameters is performed using logarithmic discretization of the scale a, and link this to the step size of steps taken between locations.

(5)

where m and are integers that control the wavelet dilation and translation respectively; a₀ is a specified fixed dilation step parameter set at a value greater than 1, and is the location parameter, which must be greater than zero. The power of two logarithmic scaling of both parameters is commonly known as the dyadic grid arrangement, which guarantees obtaining orthonormal wavelet basis. Both of the discrete wavelet parameters; the dilation and the translation parameters, are frequently given values of 2 and 1, respectively. The discrete wavelet transform could be written as follows:

(6)

where stands for the wavelet or detail coefficient at scale and location indices m and n. Nevertheless, it is significant to clarify the difference between DWT and the discretized approximations of CWT used in practice. The discretization of the CWT comprises a discrete approximation of the transform integral computed on a discrete grid of a scales and locations, and the same goes for the inverse CWT. Degrees of closeness of these approximations to the original signal mainly depends on the resolution of discretization used. On the other side, the transform integral of the DWT remains continuous, but computed only on a discretized grid of a scales and locations. Then, summation of DWT coefficients to infinity over and could take place in order to restore the original signal reliably. However, a few years later, a set of DWT-derived methods have emerged for reasons such as improving analytic techniques used for compression systems using linear combinations of discrete wavelets like Wavelet Packet Transform (WPT), and to solve problems of translation invariance that arose in both DWT and WPT. The latter is solved by introducing the Stationary Wavelet Transform (SWT) ^[111].

Several research efforts based on the WT method concerning the QRS detection were reported in literature. A multi-resolution wavelet-based system for detecting P, Q, R, S, and T peaks is reported in ^[112]. The authors also detected RRI, PP, QQ, SS, and TT time lapses. The peaks are perceived by the composition of Daubeheissub-bands wavelet of original ECG signal and the accuracy was up to 100%. Martinez et al. ^[113] developed a robust ECG delineation system which was tested on several manually annotated databases, such as MIT-BIH arrhythmia, QT, European ST-T, and CSE databases. The QRS detector obtained a sensitivity of 99.66% and a positive predictivity of 99.8% for signals taken from the MIT-BIH database. Baas et al. ^[114] examined prolongation of QT interval as a major cause of ventricular tachyarrhythmia. Three automatic T-end delineation methods based on wavelet filter banks (WAM), correlation (CORM) and Principal Component Analysis PCA (PCAM) have been developed and applied to Physionet QT database. All algorithms showed good results.

A method based on discrete wavelet transform (DWT) was suggested in ^[115], where DWT was used to extract the ECG features to complete proper classifications. All wavelet transforms share the same principle which is recovering the main signal from a mother wavelet and the wavelet coefficients. The work included data acquisition, signal preprocessing, beat detection, feature extraction, and classification. The method was mainly designed to trace non-stationary ECG signals. The DWT may lead to an optimal frequency resolution since it uses a varying window size that is narrow at higher frequencies and broad at lower frequencies. Despite the fact that DWT provides stable features to morphology variations, it does not follow each physiological temporal variation. Also, DWT offers implementation simplicity, moderate accuracy and consistency, but suffers from low predictivity and high mathematical complexity ^{[116, 117, 118]}.

The disadvantage of using DWT in bio-signal processing is that once the basic wavelet is determined, this basic wavelet must be used to perform all signal analyses. In addition, the window size chosen plays a vital role in affecting signal resolution and accuracy, which restricts the applicability of DWT to analyzing physiological signals ^{[119, 120, 121]}. Nevertheless, to classify any given system as linear, two main requirements must be satisfied: superposition and homogeneity, which is not the case with most of the biopotential signals including PPG, ECG, EMG, and EEG. The reason behind this is that such signals contain more complex physical natures and variations that could not be treated as linear systems. DWT and FFT showed limited follow ups to every single instantaneous variation in the signal ^[122]. The governing equation of DWT is given by

(7)

3.2.3. Empirical Mode Decomposition and Hilbert-Huang Transform

The Empirical Mode Decomposition (EMD) method has been first proposed by Huang et al. ^[123]. This technique has been mainly introduced to improve the analysis of nonlinear and non-stationary signals. The key role of this technique is that it breaks down any complicated data set into small and finite counts of Intrinsic Mode Functions (IMFs), all according to certain criteria. A remarkable benefit of this approach is that when a raw ECG signal gets decomposed into a number of IMFs, the grouping of these IMFs produces a signal with more distinguished QRS complexes. It also reduces ECG-ridden noise by omitting out all undesired decomposed fragments that do not meet the method main conditions. And hence, this procedure can be thought of as an adaptive filtering.

Literature shows that EMD algorithm has been used in several formulae. In ^[124], the authors applied EMD to ECG signals followed by threshold. The algorithm was also applied to ECG signals followed by singularity, followed by threshold in ^[125]. Noting that in ^[125], the researchers applied a high-pass filter to ECG signals prior to EMD filtering, and then followed by threshold.

The EMD algorithm is characterized by being robust to noise, since the first IMFs can sift out the noise and preserve the QRS content with respect to other signal components ^[124]. Thus, the first IMFs usually improve the SNR. On the other side, though the number of produced IMFs is directly proportional to the signal length, the length of these segments is not determined experimentally. The processed segments are all of the same length, and the selection criteria of IMFs are the result of applying trial-and-error methodology. Moreover, the algorithm class is numerically simple and comprises at least nine steps with more than a few specific equations for extraction. Complexity is obviously higher than derivative-based and digital filters algorithms, but, EMD algorithm reduces noise more efficiently. However, the degree of complexity is certainly higher for increasing numbers of processed ECG segments.

All physiological signals change with time due to the physiological status and hence are classified as non-linear and non-stationary. One of the best techniques for the analysis of such signals that do not usually follow regular pattern or stay still is the Hilbert-Huang Transform (HHT). HHT overcomes all of FFT and DWT limitations mentioned earlier. HHT is a combination of EMD and Hilbert transform (HT).

Once the EMD decomposes the signal into IMFs, HT then comes in to transform all selected IMFs into frequency domain. The analytic signal z(t) of the real signal x(t) is obtained by equations (8) through (10) ^{[126, 127]}.

(8)

(9)

(10)

where a(t) and ɸ(t) are the instantaneous amplitude and phase of z(t), respectively.

The instantaneous frequency ω(t) is the first derivative of the phase detected, as given by equation (11)

(11)

The representative governing equation of HT is

(12)

There are two main conditions for classifying the data set as an IMF:

1. The numbers of all detected maxima must be equal to or differs at most by one from the number of zero crossings.

2. The mean value of all points that are within both of maxima envelope and minima envelope must be equal to zero.

Researchers in ^[128] implemented HHT technique as a tool for ECG QRS detection. The outcome was a sensitivity of 99.84% and specificity of 99.92%.

3.3. Geometry Analysis
3.3.1 Poincare` Plot

A new approach of ECG signal analysis was proposed in ^[129] that entails the use of Poincare` plot to identify cardiac diseases. P-P and R-R intervals were precisely detected and combined by Poincare’ plot analysis. The basic principle of disease identification involves the determination of the geometric differences between normal and abnormal signals. The authors focused on two ECG signal categories; normal sinus rhythm and arrhythmia.

Poincare` plot analysis technique is based on a set of numeric descriptors that are useful for arrhythmia identification. The descriptors include mean square, standard deviation and geometric median positions of RR or PP Poincare` plot.

The P-wave is selected since it is closely related to arrhythmia syndrome. Compared to R peaks, P-waves are generally more challenging to detect due to low amplitude, adjacent QRS complexes, inconsistently varying shape and low SNR. Quantitative measures like amplitude, local derivative, and spatial velocity have been widely used for P-wave detection. However, these measures don’t take into account the energy content of the P-wave as a whole portion of a signal. Therefore, they fail in presenting accurate data when noise or small amplitudes are present. Consequently, Gabor Wavelet Transform, (GWT), which is governed by equation (13) below, has been used as a possible solution for P-wave detection and localization.

(13)

where c is a fixed constant that stands for the number of oscillations, σ and t are scale and time shift of function g, respectively.

Strong relationship has been reported between cardiac diseases and Poincare` plot results. The output shape of Poincare` curve was subjected to multiple linear regression models and curve fitting to describe this relationship. Moreover, this analysis technique could classify two ECG signal categories, successively ^{[130, 131]}.

3.4. Principal Component Analysis

Principal component analysis (PCA) linearly transforms the data into a new coordinate system. The technique is mainly used for pattern recognition and dimensionality reduction, while preserving informative data. PCA breaks down the data into orthogonal components that are independent of each other. All components are separated by means of the greatest variance such that each extracted component lies on its own coordinate ^{[132, 133]}.

PCA possesses a set of equivalent series expansions that help it to figure out its own calculation procedures. These series expansions are truncated orthogonal series expansion, Karhunen-Loeve series expansion, and orthogonal and orthonormal series expansions ^[134].

PCA calculations can be carried out with one of a variety of methods, such as the covariance matrix method, the singular value decomposition (SVD), the correlation matrix method, and Holelling’s T squared statistics ^[135]. For ECG analysis purposes, PCA is contemplated as an influential tool for pattern recognition, ST-T morphology, and QRS complexes detection ^{[136, 137, 138, 139]}. Fatemian et. al. in ^[140] has put forward an approach for ECG feature extraction utilizing PCA. Robust preprocessing system was implemented and noise and outliers were removed reliably. The method is capable of dealing with numerous cardiac rhythms regardless of HR and stress levels. The PCA characteristic of dimensionality reduction contributed to accelerate remarkably the classification process and thus superseded other conventional ECG classification techniques.

Researchers in ^[141] have endorsed another approach using PCA with assistance of Gaussian mixture model (GMM) as a classifier. Preprocessing and R peaks labeling succeeded using extended PT algorithm, which consists of resampling by applying FFT. Feature segmentation was achieved by putting data through linear predictive model estimation. Orthogonal and independent components were determined and considered features of interest. Features are then passed to GMM for classification. Prognostic performance is evaluated in terms of accuracy and probability of classification error. This technique resulted in accuracy of 94.29%.

3.5. Fuzzy logic systems

Fuzzy logic systems (FLS) exploit the concept of partial truth. Compared to digital binary logic systems, FLS investigates degrees of truth that lie between zero and one. Therefore, FLS has improved decision rules of judgment, since it enlarges the grey scale area between partly true and partly false. FLS have introduced remarkable evolution in computerized clinical applications. Fuzziness concepts have enrolled the depiction of possibilities among “yes” and “no” decisions through membership functions and decision rules ^[142].

Authors in ^[143] utilized FLS for the diagnosis of cardiac arrhythmia. Fuzzy logic (FL) includes fuzzy sets definitions, fuzzy rules, fuzzy interference engine design and defuzzification. Four input variables were defined as QRS duration, QTP interval “the time duration between Q and T in a QRS complex”, ratio of R-R intervals, and area of R, S and T waves. Consequently, five output variables were defined as “heartbeat case”, namely, normal, left bundle branch block (LBBB), right bundle branch block (RBBB), ventricular premature contraction (VPC), and atrial premature contractions (APC). This technique ensued with sensitivities of 95.06%, 91.03%, 90.5%, 92.63% and 93.77% for each classified case, respectively. The total classification accuracy was 93.78%.

Qidawi and Shakir in ^[144] built another fuzzy-based classification system. Raw ECG signals are normalized, preprocessed and then disintegrated into smaller frequency components each of which is related to ECG signal features. Finally, each decomposed feature is classified into a set of pre-defined categories. Fuzzy interference system (FIS) has been operated using MATLAB fuzzy logic tool box for categorization purposes. A selection of recursive statistical measures such as mean, standard deviation, median, energy “sum-squared values”, skewness, kurtosis, harmonic mean and mean deviation were optimized as classifiers. Four combinations of skewness and kurtosis were defined as input variables. Ventricular fibrillation, second degree AV block, and premature ventricular contraction are three classification cases defined as membership output variables. Five linguistic logical membership decision rules have adjudged the fate of system input variables. Results have revealed almost 100% correct detection for three diseases that are defined as system output variables.

3.6. Artificial Neural Networks

In theory, artificial neural networks (ANN) are trained to classify target outputs based on particular inputs. That is, the chosen network structure compares the emerging output with the desired target. Until a match occurs, the network keeps readjusting weights that store data gained from training sets.

An ANN stores empirical expertise based on training sets and makes knowledge available for use whenever needed. As such, ANNs simulate human brain abilities in two aspects, namely, taught familiarity by the network through supervised training, and inter-neuron connection known as synaptic weights that store such acquired knowledge.

ANNs include several neuron models and architectures such as single-layer feed-forward network, matrix-vector input, and multi-layer feed-forward network. Such network structures can be trained by several training algorithms including backpropagation ^[145], conjugate gradient algorithm ^{[146, 147]}, and Levenberg-Marquardt algorithm ^{[148, 149]}.

Jadhav et. al. ^[150] introduced a new approach for cardiac arrhythmia disease categorization. A modular neural network (MNN) was employed to separate normal from abnormal cases. The constructed neural network is characterized by altering the number of hidden layers from one to three. Classification performance is assessed by means of six measures: classification accuracy, sensitivity, specificity, mean squared error (MSE), receiver operating characteristics (ROC) and area under the curve (AUC). The obtained classification accuracy was 82.22%.

Qi Gao in ^[151] discussed computerized detection and classification of five cardiac conditions using ANNs. Executed techniques in this study were signal preprocessing, QRS detection, feature extraction and neural networks for signal classification to conduct pattern recognition. The study focused on the optimum neural network structure for ECG classification. A set of statistical features that included maximum amplitude of signal histogram, autocorrelation coefficient, QRS energy and mean or expectation vector were selected as classifiers. The core enhancement produced a three-layer network with 25 inputs, 5 neurons in the output layer and 5 neurons in its hidden layers that resulted in 91.8% recognition rate of the five cardiac conditions and an average accuracy of 84.93%. Lower complexity of neural network structures and more effective ECG features are suggested for future works to increase accuracy rate and testing speed.

ECG signal analysis and arrhythmia classification using data mining and ANNs has been devised by Gupta and Chatur in ^[152]. In their work, attempts were made to classify four types of heartbeats which are atrial premature contraction, ventricular premature contraction, left branch bundle block beat, and normal beat. This study showed that the application of ANNs made significant assistance.in coping with nonlinearities associated with the beat types considered in this study.

3.7. Neuro-Fuzzy Networks

Golpayegani and Jafari introduced a novel methodology in ECG beat recognition using adaptive neural fuzzy filters ^[153]. The interference of adaptive neural fuzzy filters (ANFF) introduced remarkable improvements on ECG classification accuracy. The authors presented a comparative study between multilayer perceptron (MLP) with back-propagation training algorithm and ANFF for early diagnosis of ECG arrhythmia. In this study, ANFF has been deployed to classify ECG signals based on the use of 4 different features. Conclusions drawn from experimental results showed that ANFF classifier is superior to MLP architecture and that techniques implemented using ANFF learn better and faster than MLP. Results also showed that a combination of Daubeches, Symlets and Biorthogonal wavelet transforms with mid 4 order auto regression model coefficients, would give optimum results. Classification accuracy of 97.6%, 100% specificity, and high sensitivity were accomplished.

A multi-lead ECG classification tactic employing neuro-fuzzy technique has been suggested in ^[154]. Unlike other researches, the authors implemented random projection (RP) rather than PCA for dimensionality reduction. In their study, a neuro-fuzzy classifier was employed to refine health state monitoring and cardiac diagnosis. Signals with normal sinus rhythm, MI, and cardiacmyopathy were acquired and used. Each signal category containing 20 beats were collected from 60 subjects. Experiments gave a recognition rate of 100% with a 25-RP coefficient number. Lower-complexity features selected for memory usage enhanced real-time operation. It was concluded that implementing a wireless, wearable sensor platform are possible future research interests.

4. HRV Analysis Techniques – State of the Art

HRV assessment is mainly affected by several factors including the time window of the signal of interest, age of the subject, physical conditions, activity, sleep or awake cycle, disease infections, drug effects, and gender. Normalizing those factors for performing HRV evaluations is crucial and can significantly affect the results of such a study.

A variety of analysis techniques for HRV have been proposed in previous researches ^[155]. Time-domain, frequency-domain, geometry-based analysis, and nonlinear techniques are the most common approaches.

4.1. Time-Domain Analysis

So far, research in time-domain analysis has been aiming at obtaining roughly estimated ranges of each statistical measure, which might lead to better indicative results for various cardiovascular diseases. In general, time domain measures used for signal analysis include standard deviation of all NN intervals (SDNN) in seconds, standard deviation of the averages of NN intervals in all 5- min segments of the entire recording (SDANN) in milliseconds, the square root of the mean of the sum of the squares of differences between adjacent NN intervals (RMSSD) in milliseconds, mean of the standard deviations of all NN intervals for all 5-min segments of the entire recording (SDNN index) in milliseconds, standard deviation of differences between adjacent NN intervals (SDSD), and number of pairs of adjacent NN intervals differing by more than 50 ms in the entire recording. In this regard, three variants are possible, namely, counting all such NN interval pairs, only pairs in which the first or the second interval is longer (NN50 count), or NN50 count divided by the total number of all NN intervals (pNN50) in percentage ^[156]. The resulting numbers that fall out of the normal case may be due to either lower or higher number of heart beats. Previous studies showed apparent implications of SDNN and pNN50 in patients with chronic heart failure (CHF) and acute myocardial infarction (AMI) [157-160]^[157]. Although RMSSD has been preferred over pNN50 due to its robustness, pNN50 indicates cardiovascular risk levels more clearly.

Having a value of SDNN that is less than 50 ms or a pNN50 value lower than 3% is regarded as an implication of high risk. In contrast, if SDNN falls between 50 ms and 100 ms, implies moderate risk and having SDNN greater than 100 ms or pNN50 over 3% is considered normal. Most of the time-domain measures are related to the evaluation of cardiovascular high-risk levels. Still, however, their beneficial reference data could also be supportive in simplifying the judgment of cardiovascular health state via analyzing HRV at rest.

Time-domain measures are characterized by simplicity of calculation, but are not sufficiently informative when utilized in a stand-alone manner. However, frequency-domain and geometry-based analysis techniques do not have reference data such as those of time-domain. Therefore, to maintain a reliable HRV-related study, data interpretation should be performed by combining more than one analysis technique.

The validity of normal values of HRV at rest in a young (18 to 25 years old), healthy and active Mexican population was investigated ^[161]. The argument involved time-domain, frequency-domain analysis methods, and a Poincare plot. A thirty-minute time window of HR recordings is acquired for 200 individuals. Significant variations were found between athletes and active subjects. These major differences, however, were not based on gender. Percentile ranges and distributions were obtained for all predefined population categories that can be used as a reference for future researches.

4.2. Frequency-Domain Analysis

Time domain measures are easy to compute, but they lack the capability of distinguishing between sympathetic and parasympathetic contributions to HRV. In comparison, frequency domain methods evaluate HRVs by examining the frequency content of each acquired ECG signal. The two main frequencies that resemble ANS activity are the low frequency component (LF) ranging from 0.04 to 0.15 Hz, and the high frequency component (HF) falling between 0.15 and 0.4 Hz. It is established that LF indicates physiologically the sympathetic modulation of heart rate while HF exploits the vagal or parasympathetic activity of ANS and matches the respiratory activities.

The index ratio LF/HF describes the sympathovagal balance and shows which part of the ANS is dominating. Each parameter (LH or HF) is computed by integrating with respect to frequency and expressed in units of (ms)². The result is also expressed in terms of Power Spectral Density (PSD) ^{[162, 163, 164, 165]}. To obtain LF and HF components, a variety of approaches can be implemented one of which is FFT that transforms the signal to frequency domain in an attempt to simplify the computation of the PSD from the measured parameters. However, FFT is not an adequate approach for analyzing non-stationary signals as previously expressed ^[166].

As indicated earlier in the text, HHT has the capabilities to handle the FFT shortcomings related to the nonlinear and nonstationary nature of the signal. Initially, EMD is applied to the signal, which resolves the signal into smaller IMFs or components. Each obtained IMF is then transformed into its own frequency domains by applying Hilbert transform ^[167]. However, previous investigations demonstrated that there exists an approximate correlation between time domain and frequency domain measures such that both are related to each other and also contribute to HRV assessment as shown in Table 3.

4.3. Geometry Analysis

Geometrical analysis methods, also known as pattern-based analysis methods, depend on the graphical properties of the ECG curve in such a way that only tangents and curvatures are the informative parameters of the curve. Geometric methods translate the resulting pattern of the histogram representing the R-R interval variations into statistical descriptive measures and numbers to assess HRV. Such measures include sample density distribution of R-R intervals, and sample density distribution of differences between R-R intervals. The triangular index is a geometric measure with the length of R-R intervals as the x-axis and the number of R-R intervals as the y-axis. The base of the triangle represents the distribution of R-R intervals with respect to the peak R-R interval. The triangular interpolation of R-R intervals (TIRR), is the baseline width of R-R interval distribution measured as a triangle base. Triangular interpolation of the resulting pattern is then performed to assess HRV as an approximate linear function ^{[168, 169]}. TIRR is highly related to SDRR in time-domain analysis, but it is highly insensitive to instantaneous fluctuations like ectopic beats or artifacts since such events are omitted out of the triangular shape of the pattern, which means it could be used to generally assess HRV ^{[170, 171]}.

Poincare` plot, also named as Lorenz plot, is a geometric technique that mainly depends on the non-linear dynamics of ECG signals. It portrays each R-R interval as a function of the preceding interval. The resulting plot is then classified and segmented into smaller portions based on the resulting shape of the curve. This technique provides descriptive information on the overall behavior of the heart on a beat-to-beat basis and also indicative data of a certain heart degree failure in a subject ^{[172, 173, 174, 175, 176]}.

4.4. Non Linear Techniques of QRS Event Extraction

Although many researchers have made use of linear analytic techniques to interpret HRV information, the cardiovascular system is too complex to be considered as a linear system. The cardiovascular system in conjunction with the ANS is continuously modulated nonlinear fluctuations. Such non-stationary perturbations could result due to interactions with other physiological systems or even small impacts of cardiac health cases like PVC or atrioventricular block. Accordingly, linear techniques using time domain, spectral domain and geometrical measures are more susceptible to generate imprecise cardiac health state prognosis.

This concern was the subject of extensive debate that has articulated the ability of spectral analyses in quantifying the autonomic tone. For example, the predominance of HF and LF measures over each other introduced remarkable support to assessing MI. However, it is worth noting that reduced LF power could not get forecasted in high risk patients with MI. Moreover, it was found that these measures host less sensitivity, specificity and positive predictivity by 30% than anticipated. Such queries pointed many researchers’ attention to analyze RRI time series with consideration of HRV dynamics.

Nonlinearities are not entirely random; it follows a chaotic nature. Thus, examining HRV signals by means of nonlinear mathematics and chaos theory appears to generate superior efficacy to linear approaches when forecasting cardiac risk potentials. However, a single nonlinear approach is not sufficient to extract useful information of HRV signals because of attributes such as the numerous pattern configurations and length of signal recordings. This section introduces the most recent nonlinear practices utilized to inspect HRV dynamics.

A dynamic system can be branded as chaotic when it can cope with variations of its own initial conditions. That is, each point located in a certain phase space could be closely approximated by other points in future trajectories. Topologically speaking, the system is capable of evolving over time, so that any given region could overlap with any other region or any set of its own phase space. If this occurs, the system could be described by a random nature that behaves according to chaos theory. The current section overviews a brief definition of each recently used nonlinear technique to analyze HRV signals. Then it proceeds to mention the beneficial clinical significance provided by each, and offers a summary highlight of literature review over the most recent years.

4.4.1. Correlation Dimensions (CD)

Generally, the complexity of a system requires a transition between time domain and phase space. Mathematically, a phase space is a space in which all possible system states or configurations are represented analogous to a unique point in the phase space. Every parameter or degree of freedom of the system represents an axis of a multidimensional space keeping in mind that a phase space can cover numerous dimensions. For a combination of these parameters, a point is repeatedly plotted in the multidimensional space. This repetition eventually describes the system evolution over time.

Essential mathematical terms such as attractor, trajectory, and attractor reconstruction are mostly dealt with when determining CDs of a dynamic system. An attractor is a state toward which a dynamical system evolves over time. In other words, if points get close to the attractor, they persist close even if slightly disturbed. Geometrically, an attractor could be a point, a curve, a surface (called a manifold), or even a strange fractal structure. A trajectory of a dynamic system on the attractor does not have to meet specific constraints other than sustaining on the attractor, and it could be periodic, chaotic, or of any other type. Attractor reconstruction tactics have been utilized to reconstruct the phase space and generate new predictive prototypes. In this concern, RRI time series representing HR is recognized as a function of time, and then recycled to construct a general attitude of the observed states of the cardiovascular system ^[177]. Embedding dimension and delay time are the main deterministic parameters of the CD technique.

To judge whether the resulting correlations are linear or nonlinear, a Surrogate Data Analysis (SDA) ^{[178, 179]} for the time series must be computed. Significance levels of differences between the CD of the surrogate and the original time series implies the existence of nonlinear correlations present in the signal.

Efforts have been put forth in quantifying nonlinear parameters like CD for pathological signals. Such efforts have proved usefulness of CD as an indicator of cardiac pathology. HRV signal tracking and prediction of onset events like ventricular tachycardia for detecting congestive heart failure situations ^[180] have been iteratively practiced using many techniques based on chaos theory, such as Grassberger and Procaccia ^[181]. Judd ^[182] has developed a more robust method of appraising CD of a tachogram. Researchers in ^[183] analyzed HRV by estimating CD of healthy and morbid groups. They found that CD of morbid group is less than the healthy group, which is supported by the findings in ^[184]. In another context, HRV has been characterized using CD to study healthy and hypertensive control groups. The study resulted in having CD values for healthy subject that are more comparable to hypertension subjects ^[185].

4.4.2. Largest Lyapunov Exponents (LLE)

The sensitivity of any dynamic system to initial conditions is quantified by the Larget Lyapunov Exponent (λ). For a bounded dynamical system, and based on its initial conditions, two nearby phase spaces keep diverging and converging relative to each other at an average rate described by (λ). In this phase space, the solutions of such system will remain finite, which forms the common mechanism of creating deterministic randomness and unpredictability. Therefore, the presence of a positive (λ) for almost all initial conditions of a bounded dynamic system is regularly used to distinguish chaotic dynamics and periodic signals ^[186]. Accordingly, this study puts forth a definition of LLE as if it is the measure that describes the average rate of exponentially diverging two nearby trajectories relative to each other position. Subsequently, a negative exponent signifies that the orbits reach a common fixed point while a zero exponent means that the orbits maintain their relative positions. Finally, a positive exponent indicates chaotic orbits.

Authors in ^[184] have performed nonlinear analysis of RRI time series obtained from both healthy individuals and cardiac patients. All signals were converted into multidimensional phase space vectors by time-delay embedding. Nonlinearity was examined by checking both LLE and CD in the light of SDA ^[178]. Finite values of CD, which indicate complexity levels, and positive LE values that imply chaotic systems, were found for all the subjects. The findings of both measures signified low dimensional chaos and that HRV in cardiac patients is less complex than healthy subjects. It was found that these results are along some previous research releases ^{[187, 188]}. Nevertheless, SDA did not firmly support these outcomes for some of the control individuals, which calls for further investigations in this regard using more refined mathematical tools. The reason behind this is that some individuals displayed chaotic nature while the HRV of other subjects showed variations at specific time stamps. Accordingly, this has not added a significant diagnostic value of this approach.

4.4.3. Fractal Dimensions

The term ‘fractal’ was first presented by Mandelbrot ^[189]. Fractal Dimension (FD) is a measure that describes how fragmented a fractal object is. It can be understood as a characterization of its self-similarity. When seen at a smaller scale, a fractal is a set of points that resemble or represents the whole set. In other words, at a certain scale, these points look similar but not necessarily identical to the same set of points when considered at different scales. FD is a non-integer quantity that identifies the fractals in the metric space defining a fractal complexity compared to another. Recently, the FD of a signal has been considered an influential tool for transient detection. Therefore, this feature has been utilized in the analysis of ECG and EEG to recognize and discriminate specific states of physiological function ^{[190, 191]}. Plentiful fractal geometry approaches could be implemented to calculate FD of an object. These approaches can be categorized as belonging to the Hausdorff-Besicovitch Dimension (like the BoxCounting and Dividers methods) or to the Bouligand-Minkowski Dimension (Minkowski FD method). The latter is the one which produces the most accurate and consistent results for FD ^{[192, 193, 194]}. In addition to these algorithms, Higuchi and Katz algorithms could not be discarded out of consideration ^{[195, 196]}.

4.4.4. Hurst Exponent (HE)

The Hurst Exponent (HE) is a dimensionless estimator used to define self-similarity of a time-series. HE has been first generated by Harold Edwin Hurst ^[176]. Meaningful values of HE are contained within the range [0,1]. HE has been widely utilized to determine the correlation properties of fractional Brownian noise. However, due to the nonstationary nature of physiological signals, some other nonlinear techniques may be preferred over HE to identify self-similarities. The HE may be computed from

(14)

where T is the duration of sample of data, and R/S is the corresponding value of rescaled range. This expression represents the generalized formula of Hurst equation that is also valid for Brownian motion. HE is also related to CD by the expression (15)

(15)

where E refers to the Euclidean dimension. Noting that values of H = 0.5 and H > 0.5 indicates uncorrelated series, values of H < 0.5 imply anti-persistent time-series having higher ranges of correlations.

HRV series may contain an irregular singularity that could be identified as oscillations at that instant. This could be quantified using HE measure keeping in mind that HRV series might have many HE values, which indicates a multi-fractal behavior. It was found that HE is highly correlated to FD and CD values for the same signal. That is, for healthy subjects, FD values are higher than those having cardiomyopathy and AF while lower values of FD are mainly due to the low R-R variation as previously stated. The authors in ^[198] performed a nonlinear analysis of HRV signals in patients with Coronary Heart Disease (CHD) using different techniques including the HE. The study investigated correlations of dynamic examination between patients with CHD and healthy control group. As estimated by R/S, RRI dynamics were altered in patients with CHD. HE and other nonlinear techniques gave independent information that cannot be detected by linear analysis techniques. The study concluded that analyzing HRV with HE in conjunction with other measures can enhance detection of CHD.

4.4.5. Approximate Entropy and Sample Entropy (ApEn, SampEn)

In physics, the term ‘entropy’ is an indication of a system randomness, regularity and predictability. By means of entropy, information rates of loss or generation enable quantifying such systems. Most of the nonlinear analysis measures have a common purpose of reliable positive predictability. However, obstacles of not being able to forecast clear future recognition are present in some HRV signals. For example, LLE value keeps decreasing for slowly varying signals like dilated/Ischemic cardiomyopathy and will be higher for other cardiac cases due to higher R-R variation ^{[199, 200]}.

The Approximate Entropy (ApEn) offers a solution to such issues and has been efficaciously applied to short and noisy data in order to attain better future recognition of some cardiac health cases. ApEn is scale invariant, model independent, and can discriminate time series. It detects variations in underlying irregular behavior not reflected on peak occurrences or amplitudes ^[201]. ApEn distinguishes time series by labeling each one a nonnegative number with larger values that represent an implication of irregularity or complexity among the spread of data. This measure has been introduced first by Pincus ^[202], and is given by

(16)

Such that:

(17)

is the correlation integral with expressing the Heavyside step function. and are, respectively, the i^th and j^th R-R interval from the tachogram of length N (number of points). The input values m (the length of compared runs), and r (tolerance level), are suggested by Goldberger et al. ^[203] to have fixed values of 2 and 0.2, respectively. Having high ApEn values signifies high complexity and irregularity in time-series data.

It is noteworthy here that in ApEn, the comparison between template vectors and all remaining vectors also includes itself, which guarantees never having zero values of This self-matching quality lowers values of ApEn. As a result, the interpreted signals look more regular than they actually are. Nevertheless, Richman and Randall ^{[204, 205]} have developed Sample Entropy (SampEn), a similar statistic family to ApEn. SampEn improved examination of time-series regularity and complexity and added benefits over ApEn in assessing dynamic human cardiovascular physiology. Both statistic families share the same mathematical mechanism of action except for vector self-comparison elimination and the dependence on N and r variables. This leads to a new correlation integral of the form

(18)

Eventually, SampEn is defined as:

(19)

with a heavyside step function defined as

(20)

SampEn monotonically decreases as r increases and is, unlike ApEn, independent of N. However, SampEn showed a limitation with excessively large confidence intervals when analyzing time series that contains less than 200 data points. On the other hand, stationarity, which is scarcely found in physiological signals, is required for both entropy measures to attain proper performance. Furthermore, outliers such as missed beats and artifacts may affect both entropy values.

4.4.6. Detrended Fluctuation Analysis (DFA)

The DFA technique is mainly used to enumerate the fractal scaling properties of short RRI signals. This technique is a modification of root-mean-square analysis of random walks applied to nonstationary signals ^[206]. At various observation windows, the root-mean-square for integrated and detrended time series is measured, and then plotted on a log-log scale against each observation window. The following formula illustrates how R-R time series of length N is integrated:

(21)

Where represents the N^th value of the integrated series, RR(i) is the i^th inter-beat interval, and is the average inter-beat interval over the entire time series.

The integrated time series is divided into windows of equal length (n). In each window, a least-squares line is fitted to the R-R interval data, which represents the trend of data within this window. The integrated time series now is detrended, yn(k), which is also the symbol standing for the y-coordinate of the resulting straight line segments. Next, the fluctuating of the integrated and detrended time series could be computed by

(22)

This computation is repeated for the entire window sizes in order to find a relationship between F(n) and window size n. As F(n) will increase with window size, small fluctuations in a small window are related to the fluctuations that could be considered a scaling exponent or the self-similarity factor α. This scaling exponent defines the slope of this line relating log (fluctuation) to log (window size). This method has been developed and applied to physiological time series by Peng et al ^[207]. It validates existence or absence of fractal correlation properties that are contained within nonstationary time series records.

Physiological determinism provided by DFA technique shows that healthy subjects exhibit α values of approximately 1, which implies fractal-like behavior. On the contrary, patients diseased with acute MI tend to display a value of α that is less than 0.75 or 0.85 ^{[206, 208]}. However, these values signpost loss of fractal-like HR dynamics. Clinically, the slope value α is very low for highly varying signals that indicate cardiac abnormalities like PVC, Left Bundle Brach Block, AF, and VF. But for rhythmically varying signals like Sick Sinus Syndrome, Complete Heart Block (CHB), and Ischemic/dilated cardiomyopathy, this value is relatively higher and is close to unity. Therefore, a set of α ranges were the result of many studies on trial signals are as follows:

0 < α < 0.5: power-law anti-correlations are present such that large values are more likely to be followed by small values and vice versa.

α = 0.5: indicates white noise.

0.5 < α < 1: power-law correlations are present such that large values are more likely to be followed by large values and vice versa. The correlation is exponential.

α = 1: special case corresponding to 1/f noise.

α > 1: correlations exist, but cease to be of a power-law form.

α = 1.5: indicates Brownian noise.

The scaling exponent α is also considered as an indication of roughness of a time series. That is, having large α values means dealing with smoother time series. On the other hand, 1/f-noise reflects a trade-off between total unpredictability of white noise (rough landscape) and the much smoother landscape of Brownian noise.

4.4.7. Recurrence Plots (RP)

In order to effectively interpret time-series data, the set under study is considered valid and relevant only if data are stationary. RP is a graphical approach using mathematics that aims at judging nonstationarity of time series. These series were first introduced by Eckmann et al. ^[209] as a diagnostic tool to examine hidden periodicity and drifts in the time propagation. For healthy subjects, the RP plot graphically illustrates that normal cases show a diagonal line and less squares implying more variation, and hence, highly fluctuating HR. On the other hand, RP visualizes abnormalities like CHB and Ischemic/dilated cardiomyopathy cases by exhibiting more squares in the plot signifying inherent periodicity and lower HR fluctuations ^[210].

4.4.8. 1/f- Slope

Frequency dependence of RRI fluctuations power was reported by Kubayashi and Musha in ^[211]. This method utilizes logarithmic interpolation resulting with a slope of regression line that relates log (power) to log (frequency) in a 1/f relationship. Logarithmic scales are used mainly due to the uneven density of data points that mostly occurs when formulating plots describing high frequency ranges. The 1/f relation is usually calculated in the 10^-2 – 10^-4 Hz range corresponding to the negative exponent β, and offers an index for long-term scaling features ^[212]. This broadband spectrum, typifying primarily slow HR variations refers to a fractal-like behavior with a long-term dependency ^[213]. For healthy normal young men, Saul el al ^[212] stated that β is near -1, which signifies linearity of the regression line. This linearity means that the RRI plots for healthy young men versus time over 2 minutes (10^-2Hz), 20 minutes (10^-3Hz), and 3 hours (10^-4Hz) may appear alike. This phenomenon is commonly known as scale-invariance or self-similarity in fractal theory. In view of that, it has been suggested that the scale invariance may be a common feature of normal physiological function. The cessation of normal physiological functioning could cause either random or periodic behavior. This is indicated by steeper 1/f slopes. Supportively, Bigger et al ^[214] conveyed a regression line with a β of nearly -1.15 in post-MI patients.

Nevertheless, this technique’s foremost limitations are the need for large sets of data, and hence, stationarity is not guaranteed in such datasets. Furthermore, artifacts and patient movements substantially disturb spectral components.

5. Conclusions

Upon the completion of this extensive review, it may be concluded that appropriate QRS complex detection techniques with sustained accuracy must be used to obtain satisfactory heart rate variability interpretations. HRV has proved its significant improvements as a clinical application specially in assessing cardiac disorders. Pan-Tomkins algorithm achieves accurate R peak detection and could be simply implemented with respect to other methodologies, but the squaring function might enlarge the noise, if any, even after filtering and thus could be replaced with a rectification stage.

As far as the HRV measurements are concerned, time-domain measures are easy to compute but lack sufficient information and, therefore, time-frequency representation of the obtained data could be considered a better option. As for frequency-domain analysis, Hilbert-Huang Transform is a strong tool when compared with FFT and DWT techniques and is highly recommended in assessing low and high frequency contents of an ECG signal. In contrast, time- and frequency-domain analysis approaches are preferred over geometry-based analysis methods due to the fact that time-frequency representation is more informative. Thereby, to obtain reliable HRV results, more than one analysis technique must be combined and implemented.

References

[1]	Kovács, P. (2012, February). ECG signal generator based on geometrical features. In Annales Univ. Sci. Budapest., Sect. Comp (Vol. 37, pp. 247-260).
	In article
[2]	Mark, R. G. (2004). HST. 542J/2.792 J/BE. 371J/6.022 J Quantitative Physiology: Organ Transport Systems.
	In article
[3]	Reisner, A. T., Clifford, G. D., & Mark, R. G. (2007). The physiological basis of the electrocardiogram.
	In article
[4]	Aldersons, A., & Buikis, A. (2011, August). Mathematical algorithm for heart rate variability analysis. In Proceedings of the 11th WSEAS international conference on applied informatics and communications, and Proceedings of the 4th WSEAS International conference on Biomedical electronics and biomedical informatics, and Proceedings of the international conference on Computational engineering in systems applications (pp. 381-386).
	In article
[5]	Camm, A. J., Malik, M., Bigger, J. T., Breithardt, G., Cerutti, S., Cohen, R. J., & Singer, D. H. (1996). Heart rate variability: standards of measurement, physiological interpretation and clinical use. Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology. Circulation, 93(5), 1043-1065.
	In article		CrossRef
[6]	Kaplan, D. T., Furman, M. I., Pincus, S. M., Ryan, S. M., Lipsitz, L. A., & Goldberger, A. L. (1991). Aging and the complexity of cardiovascular dynamics.Biophysical Journal, 59(4), 945-949.
	In article		CrossRef
[7]	Pikkujämsä, S. M., Mäkikallio, T. H., Sourander, L. B., Räihä, I. J., Puukka, P., Skyttä, J., & Huikuri, H. V. (1999). Cardiac interbeat interval dynamics from childhood to senescence comparison of conventional and new measures based on fractals and chaos theory. Circulation, 100(4), 393-399.
	In article		CrossRef
[8]	Phyllis, K., Kleiger, M. D., Robert, E., Rottman, M. D., & Jeffrey, N. (1997). Differing effects of age on heart rate variability in men and women. The American journal of cardiology, 80(3), 302-305.
	In article		CrossRef
[9]	Sinnreich, R., Kark, J. D., Friedlander, Y., Sapoznikov, D., & Luria, M. H. (1998). Five minute recordings of heart rate variability for population studies: repeatability and age–sex characteristics. Heart, 80(2), 156-162.
	In article
[10]	Zhang, J. (2007). Effect of age and sex on heart rate variability in healthy subjects. Journal of manipulative and physiological therapeutics, 30(5), 374-379.
	In article		CrossRef
[11]	Carter, J. B., Banister, E. W., & Blaber, A. P. (2003). The effect of age and gender on heart rate variability after endurance training. Medicine and science in sports and exercise, 35(8), 1333-1340.
	In article		CrossRef
[12]	Huikuri, H. V., Pikkuja, S. M., Airaksinen, K. J., Ika, M. J., Rantala, A. O., Kauma, H., & Kesa, Y. A. (1996). Sex-related differences in autonomic modulation of heart rate in middle-aged subjects. Circulation, 94(2), 122-125.
	In article		CrossRef
[13]	Ramaekers, D., Ector, H., Aubert, A. E., Rubens, A., & Van de Werf, F. (1998). Heart rate variability and heart rate in healthy volunteers. Is the female autonomic nervous system cardioprotective?. European Heart Journal, 19(9), 1334-1341.
	In article		CrossRef
[14]	De La Cruz Torres, B., López, C. L., & Orellana, J. N. (2008). Analysis of heart rate variability at rest and during aerobic exercise: a study in healthy people and cardiac patients. British journal of sports medicine, 42(9), 715-720.
	In article		CrossRef
[15]	Bernardi, L., Ricordi, L., Lazzari, P., Solda, P., Calciati, A., Ferrari, M. R., & Fratino, P. (1992). Impaired circadian modulation of sympathovagal activity in diabetes. A possible explanation for altered temporal onset of cardiovascular disease. Circulation, 86(5), 1443-1452.
	In article		CrossRef
[16]	Guzzetti, S. T. E. F. A. N. O., Cogliati, C. H. I. A. R. A., Broggi, C. A. R. O. L. A., Carozzi, C. A. R. L. A., Caldiroli, D., Lombardi, F. E. D. E. R. I. C. O., & Malliani, A. L. B. E. R. T. (1994). Influences of neural mechanisms on heart period and arterial pressure variabilities in quadriplegic patients. American Journal of Physiology-Heart and Circulatory Physiology, 266(3), H1112-H1120.
	In article
[17]	Koh, J., Brown, T. E., Beightol, L. A., Ha, C. Y., & Eckberg, D. L. (1994). Human autonomic rhythms: vagal cardiac mechanisms in tetraplegic subjects.The Journal of physiology, 474(3), 483-495.
	In article
[18]	Garrido Esquivel, A., de la CruzTorres, B., Garrido Salazar, M. A., Medina Corrales, M., & Naranjo Orellana, J. (2009). Variabilidad de la frecuencia cardiaca en un deportista juvenil durante una competición de bádminton de máximo nivel. Revista Andaluza de Medicina del Deporte, 2(2), 70-74.
	In article
[19]	Saa, Y. D., Sarmiento, S., Martín González, J. M., Rodríguez Ruiz, D., Quiroga, M. E., & García Manso, J. M. (2009). Aplicación de la variabilidad de la frecuencia cardiaca en la caracterización de deportistas de élite de lucha canaria con diferente nivel de rendimiento. Archivos de medicina del deporte,2(4), 120-125.
	In article
[20]	Aubert, A. E., Seps, B., & Beckers, F. (2003). Heart rate variability in athletes.Sports Medicine, 33(12), 889-919.
	In article		CrossRef
[21]	Aubert, A. E., Beckers, F., & Ramaekers, D. (2000). Short-term heart rate variability in young athletes. Journal of cardiology, 37, 85-88.
	In article
[22]	Bonnemeier, H., Wiegand, U. K., Brandes, A., Kluge, N., Katus, H. A., Richardt, G., & Potratz, J. (2003). Circadian profile of cardiac Journal of cardiovascular electrophysiology, 14(8), 791-799.
	In article		CrossRef
[23]	Hall, M., Thayer, J. F., Germain, A., Moul, D., Vasko, R., Puhl, M., & Buysse, D. J. (2007). Psychological stress is associated with heightened physiological arousal during NREM sleep in primary insomnia. Behavioral sleep medicine, 5(3), 178-193.
	In article		CrossRef
[24]	Montano, N., Porta, A., Cogliati, C., Costantino, G., Tobaldini, E., Casali, K. R., & Iellamo, F. (2009). Heart rate variability explored in the frequency domain: a tool to investigate the link between heart and behavior. Neuroscience & Biobehavioral Reviews, 33(2), 71-80.
	In article		CrossRef
[25]	De Boer, R. W., Karemaker, J. M., & Strackee, J. (1985). Relationships between short-term blood-pressure fluctuations and heart-rate variability in resting subjects I: a spectral analysis approach. Medical and Biological Engineering and Computing, 23(4), 352-358.
	In article		CrossRef
[26]	Rothschild, M., Rothschild, A., & Pfeifer, M. (1988). Temporary decrease in cardiac parasympathetic tone after acute myocardial infarction. The American journal of cardiology, 62(9), 637-639.
	In article		CrossRef
[27]	Rossinen, M. D., Viitasalo, M. D., Partanen, M. D., Koskinen, M. D., Kupari, M. D., Nieminen, M. D., & Markku, S. (1997). Effects of acute alcohol ingestion on heart rate variability in patients with documented coronary artery disease and stable angina pectoris. The American journal of cardiology, 79(4), 487-491.
	In article		CrossRef
[28]	Carney, R. M., Blumenthal, J. A., Stein, P. K., Watkins, L., Catellier, D., Berkman, L. F., & Freedland, K. E. (2001). Depression, heart rate variability, and acute myocardial infarction. Circulation, 104(17), 2024-2028.
	In article		CrossRef
[29]	Carney, R. M., Blumenthal, J. A., Freedland, K. E., Stein, P. K., Howells, W. B., Berkman, L. F., & Jaffe, A. S. (2005). Low heart rate variability and the effect of depression on post-myocardial infarction mortality. Archives of internal medicine, 165(13), 1486.
	In article		CrossRef
[30]	Pfeifer, M. A., Cook, D., Brodsky, J., Tice, D., Reenan, A., Swedine, S., & Porte, D. (1982). Quantitative evaluation of cardiac parasympathetic activity in normal and diabetic man. Diabetes, 31(4), 339-345.
	In article		CrossRef
[31]	Singh, J. P., Larson, M. G., O’Donnell, C. J., Wilson, P. F., Tsuji, H., Lloyd-Jones, D. M., & Levy, D. (2000). Association of hyperglycemia with reduced heart rate variability (The Framingham Heart Study). The American journal of cardiology, 86(3), 309-312.
	In article		CrossRef
[32]	Wheeler, T., & Watkins, P. J. (1973). Cardiac denervation in diabetes. British Medical Journal, 4(5892), 584.
	In article		CrossRef
[33]	Forsström, J., Forsström, J., Heinonen, E., Välimäki, I., & Antila, K. (1986). Effects of haemodialysis on heart rate variability in chronic renal failure.Scandinavian journal of clinical & laboratory investigation, 46(7), 665-670.
	In article		CrossRef
[34]	Zoccali, C., Ciccarelli, M., & Maggiore, Q. (1982). Defective reflex control of heart rate in dialysis patients: evidence for an afferent autonomic lesion. Clinical Science, 63, 285-292.
	In article
[35]	Lerma, C., Minzoni, A., Infante, O., & José, M. V. (2004). A mathematical analysis for the cardiovascular control adaptations in chronic renal failure.Artificial organs, 28(4), 398-409.
	In article		CrossRef
[36]	Ewing, D. J., & Winney, R. (1975). Autonomic function in patients with chronic renal failure on intermittent haemodialysis. Nephron, 15(6), 424-429.
	In article		CrossRef
[37]	Nagy, E., Orvos, H., Bárdos, G., & Molnár, P. (2000). Gender-related heart rate differences in human neonates. Pediatric Research, 47(6), 778-780.
	In article		CrossRef
[38]	Spallone, V., Bernardi, L., Ricordi, L., Soldà, P., Maiello, M. R., Calciati, A., & Menzinger, G. (1993). Relationship between the circadian rhythms of blood pressure and sympathovagal balance in diabetic autonomic neuropathy.Diabetes, 42(12), 1745-1752.
	In article		CrossRef
[39]	van Ravenswaaij-Arts, C., Hopman, J. C., Kollée, L. A., van Amen, J. P., Stoelinga, G., & van Geijn, H. P. (1991). Influences on heart rate variability in spontaneously breathing preterm infants. Early human development, 27(3), 187-205.
	In article		CrossRef
[40]	Bekheit, S., Tangella, M., el-Sakr, A., Rasheed, Q., Craelius, W., & El-Sherif, N. (1990). Use of heart rate spectral analysis to study the effects of calcium channel blockers on sympathetic activity after myocardial infarction. American heart journal, 119(1), 79-85.
	In article		CrossRef
[41]	Coumel, P., Hermida, J. S., Wennerblöm, B., Leenhardt, A., Maison-Blanche, P., & Cauchemez, B. (1991). Heart rate variability in left ventricular hypertrophy and heart failure, and the effects of beta-blockade a non-spectral analysis of heart rate variability in the frequency domain and in the time domain. European heart journal, 12(3), 412-422.
	In article
[42]	Guzzetti, S., Piccaluga, E., Casati, R., Cerutti, S., Lombardi, F., Pagani, M., & Malliani, A. (1988). Sympathetic predominance an essential hypertension: a study employing spectral analysis of heart rate variability. Journal of hypertension, 6(9), 711-717.
	In article		CrossRef
[43]	Lucini, D., Bertocchi, F., Malliani, A., & Pagani, M. (1996). A controlled study of the autonomic changes produced by habitual cigarette smoking in healthy subjects. Cardiovascular research, 31(4), 633-639.
	In article		CrossRef
[44]	Hayano, J., Yamada, M., Sakakibara, Y., Fujinami, T., Yokoyama, K., Watanabe, Y., & Takata, K. (1990). Short-and long-term effects of cigarette smoking on heart rate variability. The American journal of cardiology, 65(1), 84-88.
	In article		CrossRef
[45]	Kamath, M. V., & Fallen, E. L. (1995). Correction of the heart rate variability signal for ectopics and missing beats. Heart rate variability. Armonk: Futura, 75-85.
	In article
[46]	Zeskind, P. S., & Gingras, J. L. (2006). Maternal cigarette-smoking during pregnancy disrupts rhythms in fetal heart rate. Journal of pediatric psychology,31(1), 5-14.
	In article		CrossRef
[47]	Webster, J. G. Medical instrumentation: application and design. 1998. John Wiley&Sons, NY.
	In article
[48]	Manolakis, D. G., Ingle, V. K., & Kogon, S. M. (2000). Statistical and adaptive signal processing: spectral estimation, signal modeling, adaptive filtering, and array processing (pp. 378-387). Boston: McGraw-Hill.
	In article
[49]	Li, Q., Mark, R. G., & Clifford, G. D. (2008). Robust heart rate estimation from multiple asynchronous noisy sources using signal quality indices and a Kalman filter. Physiological measurement, 29(1), 15.
	In article		CrossRef
[50]	Goutas, A., Ferdi, Y., Herbeuval, J. P., Boudraa, M., & Boucheham, B. (2005). Digital fractional order differentiation-based algorithm for P and T-waves detection and delineation. ITBM-RBM, 26(2), 127-132.
	In article		CrossRef
[51]	Pan, J., & Tompkins, W. J. (1985). A real-time QRS detection algorithm.Biomedical Engineering, IEEE Transactions on, (3), 230-236.
	In article		CrossRef
[52]	Okada, M. (1979). A digital filter for the ors complex detection. Biomedical Engineering, IEEE Transactions on, (12), 700-703.
	In article		CrossRef
[53]	Menrad, A. (1981). Dual microprocessor system for cardiovascular data acquisition, processing and recording. Proc. 1981 IEEE Inr. Con5 Industrial Elect. Contr. Instrument, 64-69.
	In article
[54]	Holsinger, W. P., Kempner, K. M., & Miller, M. H. (1971). A QRS preprocessor based on digital differentiation. Biomedical Engineering, IEEE Transactions on, (3), 212-217.
	In article		CrossRef
[55]	Morizet-Mahoudeaux, P., Moreau, C., Moreau, D., & Quarante, J. J. (1981). Simple microprocessor-based system for on-line ECG arrhythmia analysis.Medical and Biological Engineering and Computing, 19(4), 497-500.
	In article		CrossRef
[56]	Fraden, J., & Neuman, M. R. (1980). QRS wave detection. Medical and Biological Engineering and computing, 18(2), 125-132.
	In article		CrossRef
[57]	Balda, R. A., Diller, G., Deardorff, E., Doue, J., & Hsieh, P. (1977). The HP ECG analysis program. Trends in Computer-Processed Electrocardiograms, 197-205.
	In article
[58]	Ahlstrom, M. L., & Tompkins, W. J. (1983). Automated high-speed analysis of Holter tapes with microcomputers. Biomedical Engineering, IEEE Transactions on, (10), 651-657.
	In article		CrossRef
[59]	Engelse, W. A. H., & Zeelenberg, C. (1979). A single scan algorithm for QRS detection and feature extraction. Computers in cardiology, 6(1979), 37-42.
	In article
[60]	Trahanias, P., & Skordalakis, E. (1989). Bottom-up approach to the ECG pattern-recognition problem. Medical and Biological Engineering and Computing, 27(3), 221-229.
	In article		CrossRef
[61]	Zhang, F., & Lian, Y. (2007, August). Electrocardiogram QRS detection using multiscale filtering based on mathematical morphology. In Engineering in Medicine and Biology Society, 2007. EMBS 2007. 29th Annual International Conference of the IEEE (pp. 3196-3199). IEEE.
	In article
[62]	Arzeno, N. M., Deng, Z. D., & Poon, C. S. (2008). Analysis of first-derivative based QRS detection algorithms. Biomedical Engineering, IEEE Transactions on, 55(2), 478-484.
	In article		CrossRef
[63]	Benitez, D. S., Gaydecki, P. A., Zaidi, A., & Fitzpatrick, A. P. (2000). A new QRS detection algorithm based on the Hilbert transform. In Computers in Cardiology 2000 (pp. 379-382). IEEE.
	In article
[64]	Sufi, F., Fang, Q., & Cosic, I. (2007, August). Ecg rr peak detection on mobile phones. In Engineering in Medicine and Biology Society, 2007. EMBS 2007. 29th Annual International Conference of the IEEE (pp. 3697-3700). IEEE.
	In article
[65]	Zhang, F., & Lian, Y. (2007, November). Novel QRS detection by CWT for ECG sensor. In Biomedical Circuits and Systems Conference, 2007. BIOCAS 2007. IEEE (pp. 211-214). IEEE.
	In article
[66]	Friesen, G. M., Jannett, T. C., Jadallah, M. A., Yates, S. L., Quint, S. R., & Nagle, H. T. (1990). A comparison of the noise sensitivity of nine QRS detection algorithms. Biomedical Engineering, IEEE Transactions on, 37(1), 85-98..
	In article		CrossRef
[67]	Christov, I. I. (2004). Real time electrocardiogram QRS detection using combined adaptive threshold. BioMedical Engineering OnLine, 3(1), 28.
	In article		CrossRef
[68]	Kohler, B. U., Hennig, C., & Orglmeister, R. (2002). The principles of software QRS detection. Engineering in Medicine and Biology Magazine, IEEE, 21(1), 42-57..
	In article		CrossRef
[69]	Vidal, C., Charnay, P., & Arce, P. (2009, January). Enhancement of a QRS detection algorithm based on the first derivative, using techniques of a QRS detector algorithm based on non-linear transformations. In 4th European Conference of the International Federation for Medical and Biological Engineering (pp. 393-396). Springer Berlin Heidelberg.
	In article		CrossRef
[70]	Ahlstrom, M. L., & Tompkins, W. J. (1983). Automated high-speed analysis of Holter tapes with microcomputers. Biomedical Engineering, IEEE Transactions on, (10), 651-657.
	In article		CrossRef
[71]	Macfarlane, P. W., Devine, B., & Clark, E. (2005, September). The university of Glasgow (Uni-G) ECG analysis program. In Computers in Cardiology, 2005 (pp. 451-454). IEEE.
	In article
[72]	Olvera, F. E. (2006). Electrocardiogram waveform feature extraction using the matched filter. ECE SIO Stat Proc.
	In article
[73]	Kay, S. M. (1998). Fundamentals of Statistical signal processing, Volume 2: Detection theory (pp. 345-349). Prentice Hall PTR.
	In article
[74]	Hamilton, P. S., & Tompkins, W. J. (1988, November). Adaptive matched filtering for QRS detection. In Engineering in Medicine and Biology Society, 1988. Proceedings of the Annual International Conference of the IEEE (pp. 147-148). IEEE.
	In article
[75]	Dobbs, S. E., Schmitt, N. M., & Ozemek, H. S. (1984). QRS detection by template matching using real-time correlation on a microcomputer. Journal of clinical engineering, 9(3), 197-212.
	In article		CrossRef
[76]	Ebenezer, D., & Krishnamurthy, V. (1993). Wave digital matched filter for electrocardiogram preprocessing. Journal of biomedical engineering, 15(2), 132-134.
	In article		CrossRef
[77]	Kaplan, D. T. (1990, September). Simultaneous QRS detection and feature extraction using simple matched filter basis functions. In Computers in Cardiology 1990, Proceedings. (pp. 503-506). IEEE.
	In article
[78]	Ruha, A., Sallinen, S., & Nissila, S. (1997). A real-time microprocessor QRS detector system with a 1-ms timing accuracy for the measurement of ambulatory HRV. Biomedical Engineering, IEEE Transactions on, 44(3), 159-167.
	In article		CrossRef
[79]	Rangayyan, R. M. (2002). Biomedical signal analysis (pp. 55-91). New York: IEEE press.
	In article
[80]	Xue, Q., Hu, Y. H., & Tompkins, W. J. (1992). Neural-network-based adaptive matched filtering for QRS detection. Biomedical Engineering, IEEE Transactions on, 39(4), 317-329.
	In article		CrossRef
[81]	Lu, Y., Xian, Y., Chen, J., & Zheng, Z. (2008, May). A comparative study to extract the diaphragmatic electromyogram signal. In BioMedical Engineering and Informatics, 2008. BMEI 2008. International Conference on (Vol. 2, pp. 315-319). IEEE.
	In article
[82]	Gupta, R., Chatterjee, H. K., & Mitra, M. (2012). An online ECG QRS Detection Technique. Int. J. on Recent Trends in Engineering and Technology, 7(2).
	In article
[83]	Mehta, S. S., & Lingayat, N. S. (2007). Comparative study of QRS detection in single lead and 12-lead ecg based on entropy and combined entropy criteria using support vector machine. Journal of Theoretical and Applied Information Technology, 3(2), 8-18.
	In article
[84]	Mehta, S. S., & Lingayat, N. S. (2009). Identification of QRS complexes in 12-lead electrocardiogram. Expert Systems with Applications, 36(1), 820-828.
	In article		CrossRef
[85]	Nouira, I., Abdallah, A. B., Bedoui, M. H., & Dogui, M. A Robust R Peak Detection Algorithm Using Wavelet Transform for Heart Rate Variability Studies, International Journal on Electrical Engineering and Informatics - Volume 5,Number 3, September 2013.
	In article
[86]	Hamilton, P. S., & Tompkins, W. J. (1986). Quantitative investigation of QRS detection rules using the MIT/BIH arrhythmia database. Biomedical Engineering, IEEE Transactions on, (12), 1157-1165.
	In article		CrossRef
[87]	Hargittai, S. (2005, September). Savitzky-Golay least-squares polynomial filters in ECG signal processing. In Computers in Cardiology, 2005 (pp. 763-766). IEEE.
	In article
[88]	Luo, J., Ying, K., & Bai, J. (2005). Savitzky–Golay smoothing and differentiation filter for even number data. Signal Processing, 85(7), 1429-1434.
	In article		CrossRef
[89]	Das, S., & Chakraborty, M. (2012). QRS Detection Algorithm Using Savitzky-Golay Filter. Aceee International Journal on signal & Image processing, 3(1).
	In article
[90]	Fira, C. M., & Goras, L. (2008). An ECG signals compression method and its validation using NNs. Biomedical Engineering, IEEE Transactions on, 55(4), 1319-1326.
	In article		CrossRef
[91]	Suppappola, S. E. T. H., & Sun, Y. (1994). Nonlinear transforms of ECG signals for digital QRS detection: a quantitative analysis. Biomedical Engineering, IEEE Transactions on, 41(4), 397-400.
	In article		CrossRef
[92]	Lin, C. C., Hu, W. C., Chen, C. M., & Weng, C. H. (2008, September). Heart rate detection in highly noisy handgrip electrocardiogram. In Computers in Cardiology, 2008 (pp. 477-480). IEEE.
	In article
[93]	Ulusar, U. D., Govindan, R. B., Wilson, J. D., Lowery, C. L., Preissl, H., & Eswaran, H. (2009, September). Adaptive rule based fetal QRS complex detection using Hilbert transform. In Engineering in Medicine and Biology Society, 2009. EMBC 2009. Annual International Conference of the IEEE (pp. 4666-4669). IEEE.
	In article
[94]	Clifford, G. D., & Azuaje, F. (2006). Advanced methods and tools for ECG data analysis (pp. 55-57). London: Artech house.
	In article
[95]	Cost, A. A., & Cano, G. G. (1989, November). QRS detection based on hidden Markov modeling. In Engineering in Medicine and Biology Society, 1989. Images of the Twenty-First Century., Proceedings of the Annual International Conference of the IEEE Engineering in (pp. 34-35). IEEE.
	In article
[96]	Petrutiu, S., Ng, J., Nijm, G. M., Al-Angari, H., Swiryn, S., & Sahakian, A. V. (2006). Atrial fibrillation and waveform characterization. Engineering in Medicine and Biology Magazine, IEEE, 25(6), 24-30.
	In article		CrossRef
[97]	Coast, D. A., Stern, R. M., Cano, G. G., & Briller, S. A. (1990). An approach to cardiac arrhythmia analysis using hidden Markov models. Biomedical Engineering, IEEE Transactions on, 37(9), 826-836.
	In article		CrossRef
[98]	Krimi, S., Ouni, K., & Ellouze, N. (2008, April). An Approach Combining Wavelet Transform and Hidden Markov Models for ECG Segmentation. InInformation and Communication Technologies: From Theory to Applications, 2008. ICTTA 2008. 3rd International Conference on (pp. 1-6). IEEE.
	In article
[99]	Oppenheim, A. V., Schafer, R. W., & Buck, J. R. (1999). Discrete-time signal processing (Vol. 5). Upper Saddle River: Prentice Hall.
	In article
[100]	Li, C., Zheng, C., & Tai, C. (1995). Detection of ECG characteristic points using wavelet transforms. Biomedical Engineering, IEEE Transactions on, 42(1), 21-28.
	In article		CrossRef
[101]	Taubin, G., Zhang, T., & Golub, G. (1996). Optimal surface smoothing as filter design (pp. 283-292). Springer Berlin Heidelberg.
	In article
[102]	Chu, C. H., & Delp, E. J. (1989). Impulsive noise suppression and background normalization of electrocardiogram signals using morphological operators.Biomedical Engineering, IEEE Transactions on, 36(2), 262-273.
	In article		CrossRef
[103]	Chen, Y., & Duan, H. (2006, January). A QRS complex detection algorithm based on mathematical morphology and envelope. In Engineering in Medicine and Biology Society, 2005. IEEE-EMBS 2005. 27th Annual International Conference of the (pp. 4654-4657). IEEE.
	In article
[104]	Zhang, F., & Lian, Y. (2007, August). Electrocardiogram QRS detection using multiscale filtering based on mathematical morphology. In Engineering in Medicine and Biology Society, 2007. EMBS 2007. 29th Annual International Conference of the IEEE  (pp. 3196-3199). IEEE.
	In article
[105]	Zhang, F., & Lian, Y. (2011). QRS detection based on morphological filter and energy envelope for applications in body sensor networks. Journal of Signal Processing Systems, 64(2), 187-194.
	In article		CrossRef
[106]	Zhang, C. F., & Bae, T. W. (2012). VLSI friendly ECG QRS complex detector for body sensor networks. Emerging and Selected Topics in Circuits and Systems, IEEE Journal on, 2(1), 52-59.
	In article		CrossRef
[107]	Bracewell, R. N., & Bracewell, R. N. (1986). The Fourier transform and its applications (Vol. 31999). New York: McGraw-Hill.
	In article
[108]	Henning, C. (2002). The Principles of Software QRS Detection. IEEE Engineering in Medicine and Biology, 21, 42-57.
	In article		CrossRef
[109]	Pantelopoulos, A., & Bourbakis, N. G. (2010). Prognosis-a wearable health-monitoring system for people at risk: methodology and modeling. IEEE transactions on information technology in biomedicine: a publication of the IEEE Engineering in Medicine and Biology Society, 14(3), 613-621.
	In article		CrossRef
[110]	Addison, P. S. (2010). The illustrated wavelet transform handbook: introductory theory and applications in science, engineering, medicine and finance. CRC Press.
	In article
[111]	Coifman, R. R., & Donoho, D. L. (1995). Translation-invariant de-noising (pp. 125-150). Springer New York.
	In article
[112]	Mukhopadhyay, S., Biswas, S., Roy, A. B., & Dey, N. (2012). Wavelet Based QRS Complex Detection of ECG Signal. arXiv preprint arXiv:1209.1563.
	In article
[113]	Martínez, J. P., Almeida, R., Olmos, S., Rocha, A. P., & Laguna, P. (2004). A wavelet-based ECG delineator: evaluation on standard databases. Biomedical Engineering, IEEE Transactions on, 51(4), 570-581.
	In article		CrossRef
[114]	Baas, T., Gravenhorst, F., Fischer, R., Khawaja, A., & Dossel, O. (2010, September). Comparison of three t-wave delineation algorithms based on wavelet filterbank, correlation and pca. In Computing in Cardiology, 2010 (pp. 361-364). IEEE.
	In article
[115]	Tamil, E. B. M., Kamarudin, N. H., Salleh, R., & Tamil, A. M. (2008, January). A Review on Feature Extraction & Classification Techniques for Biosignal Processing (Part I: Electrocardiogram). In 4th Kuala Lumpur International Conference on Biomedical Engineering 2008 (pp. 107-112). Springer Berlin Heidelberg.
	In article
[116]	Jaswal, G., Parmar, R., & Kaul, A. (2012). QRS Detection Using Wavelet Transform. International Journal of Engineering and Advanced Technology (IJEAT), 1.
	In article
[117]	Ghaffari, A., Golbayani, H., & Ghasemi, M. (2008). A new mathematical based QRS detector using continuous wavelet transform. Computers & Electrical Engineering, 34(2), 81-91.
	In article		CrossRef
[118]	Reddy, G. U., Muralidhar, M., & Varadarajan, S. (2009). ECG De-Noising using improved thresholding based on Wavelet transforms. IJCSNS, 9(9), 221.
	In article
[119]	Chapron, B., & Bliven, L. (1989, July). Wavelet analysis introduction and application to radar scattering from water waves. In Geoscience and Remote Sensing Symposium, 1989. IGARSS'89. 12th Canadian Symposium on Remote Sensing., 1989 International (Vol. 3, pp. 1474-1477). IEEE.
	In article
[120]	Addison, P. S. (2005). Wavelet transforms and the ECG: a review. Physiological measurement, 26(5), R155.
	In article		CrossRef
[121]	Clifford, G. D., & Azuaje, F. (2006). Advanced methods and tools for ECG data analysis (pp. 55-57). London: Artech house.
	In article
[122]	Mallat, S. (1999). A wavelet tour of signal processing. Access Online via Elsevier.
	In article
[123]	Hongyan, X., & Minsong, H. (2008, May). A new QRS detection algorithm based on empirical mode decomposition. In Bioinformatics and Biomedical Engineering, 2008. ICBBE 2008. The 2nd International Conference on (pp. 693-696). IEEE.
	In article
[124]	Tang, J. T., Yang, X. L., Xu, J. C., Tang, Y., Zou, Q., & Zhang, X. K. (2008, October). The Algorithm of R peak detection in ECG based on empirical Mode Decomposition. In Natural Computation, 2008. ICNC'08. Fourth International Conference on (Vol. 5, pp. 624-627). IEEE.
	In article
[125]	Arafat, A., & Hasan, K. (2009, April). Automatic detection of ECG wave boundaries using empirical mode decomposition. In Acoustics, Speech and Signal Processing, 2009. ICASSP 2009. IEEE International Conference on (pp. 461-464). IEEE.
	In article
[126]	Oweis, R. J., & Abdulhay, E. W. (2011). Seizure classification in EEG signals utilizing Hilbert-Huang transform. Biomedical engineering online, 10(1), 38.
	In article		CrossRef
[127]	Zhang, Q. (2010). Cuff-free blood pressure estimation using signal processing techniques (Doctoral dissertation, University of Saskatchewan).
	In article
[128]	Hadj Slimane, Z. E., & Naït-Ali, A. (2010). QRS complex detection using Empirical Mode Decomposition. Digital Signal Processing, 20(4), 1221-1228.
	In article		CrossRef
[129]	Kim, D. H., KIM, J., & Youn, C. H. (2009). Poincare Geometry-Characterized Arrhythmia Identification Scheme in Grid. International journal of engineering science and technology, 1(3).
	In article
[130]	Salahuddin, L., & Kim, D. (2006). Detection of acute stress by heart rate variability using a prototype mobile ECG sensor. Hybrid Information Technology, ICHIT, 6, 453-459.
	In article
[131]	Sovilj, S., Jeras, M., & Magjarevic, R. (2004, May). Real time P-wave detector based on wavelet analysis. In Electrotechnical Conference, 2004. MELECON 2004. Proceedings of the 12th IEEE Mediterranean (Vol. 1, pp. 403-406). IEEE.
	In article
[132]	Miranda, A. A., Le Borgne, Y. A., & Bontempi, G. (2008). New routes from minimal approximation error to principal components. Neural Processing Letters, 27(3), 197-207.
	In article		CrossRef
[133]	Abdi, H., & Williams, L. J. (2010). Principal component analysis. Wiley Interdisciplinary Reviews: Computational Statistics, 2(4), 433-459.
	In article		CrossRef
[134]	Sörnmo, L., & Laguna, P. (2005). Bioelectrical signal processing in cardiac and neurological applications [electronic resource]. Academic Press.
	In article
[135]	Khawaja, A. (2006). Automatic ECG analysis using principal component analysis and wavelet transformation. Univ.-Verlag Karlsruhe.
	In article
[136]	Teodorescu, H. N., & Bonciu, C. (1996, August). Feedforward neural filter with learning in features space. Preliminary results. In Neuro-Fuzzy Systems, 1996. AT'96., International Symposium on (pp. 17-24). IEEE.
	In article
[137]	Costa, E. V., & Moraes, J. C. T. B. (2000). QRS feature discrimination capability: quantitative and qualitative analysis. In Computers in Cardiology 2000 (pp. 399-402). IEEE.
	In article
[138]	Stamkopoulos, T., Diamantaras, K., Maglaveras, N., & Strintzis, M. (1998). ECG analysis using nonlinear PCA neural networks for ischemia detection.Signal Processing, IEEE Transactions on, 46(11), 3058-3067.
	In article		CrossRef
[139]	Vargas, F., Lettnin, D., de Castro, M. C. F., & Macarthy, M. (2002). Electrocardiogram pattern recognition by means of MLP network and PCA: A case study on equal amount of input signal types. In Neural Networks, 2002. SBRN 2002. Proceedings. VII Brazilian Symposium on (pp. 200-205). IEEE.
	In article
[140]	Fatemian, S. Z., & Hatzinakos, D. (2009, July). A new ECG feature extractor for biometric recognition. In Digital Signal Processing, 2009 16th International Conference on (pp. 1-6). IEEE.
	In article
[141]	Martis, R. J., Chakraborty, C., & Ray, A. K. (2009). A two-stage mechanism for registration and classification of ECG using Gaussian mixture model. Pattern Recognition, 42(11), 2979-2988.
	In article		CrossRef
[142]	Upasani, D. E., & Kharadkar, R. D. Automated ECG Diagnosis. IOSR Journal of Engineering Vol.2(5), pp: 1265-1269, May. 2012.
	In article
[143]	Yeh, Y. C., Wang, W. J., & Chiou, C. W. (2009). Heartbeat case determination using fuzzy logic method on ECG signals. International Journal of Fuzzy Systems, 11(4), 250-261.
	In article
[144]	Qidwai, U., & Shakir, M. (2012). Filter Bank Approach to Critical Cardiac Abnormalities Detection using ECG data under Fuzzy Classification. in International Journal of Computer Information Systems & Industrial Management Applications (July 2012) ISSN, 2150-7988.
	In article
[145]	Hasnain, S. K. U., & Asim, S. M. (1999). Artificial Neural network in Cardiology-ECG Wave Analysis and Diagnosis Using Backpropogation Neural network.
	In article
[146]	Nocedal, J., & Wright, S. (2006). Numerical optimization, series in operations research and financial engineering. Springer, New York.
	In article
[147]	Jing-tian, T., Qing, Z., Yan, T., Bin, L., & Xiao-kai, Z. (2007, July). Hilbert-Huang transform for ECG de-noising. In Bioinformatics and Biomedical Engineering, 2007. ICBBE 2007. The 1st International Conference on (pp. 664-667). IEEE.
	In article
[148]	Levenberg, K. (1944). A method for the solution of certain problems in least squares. Quarterly of applied mathematics, 2, 164-168.
	In article
[149]	Marquardt, D. W. (1963). An algorithm for least-squares estimation of nonlinear parameters. Journal of the Society for Industrial & Applied Mathematics, 11(2), 431-441.
	In article		CrossRef
[150]	Jadhav, S. M., Nalbalwar, S. L., & Ghatol, A. A. (2011). Modular neural network based arrhythmia classification system using ECG signal data. International Journal of Information Technology and Knowledge Management, 4(1), 205-209.
	In article
[151]	George Qi Gao. Computerized detection and classification of five cardiac conditions. Thesis submitted in partial fulfillment of the degree of master of engineering, Auckland university of technology, New Zealand, May, 2003.
	In article
[152]	Gupta, K. O., & Chatur, D. P. (2012). ECG Signal Analysis and Classification using Data Mining and Artificial Neural Networks. Int. J. Emerg. Technol. Adv. Eng, 2, 56-60.
	In article
[153]	Golpayegani, G. N., & Jafari, A. H. (2009). A novel approach in ECG beat recognition using adaptive neural fuzzy filter. Journal of Biomedical Science and Engineering, 2(2), 80-85.
	In article		CrossRef
[154]	Bogdanova Vandergheynst, I., Vallejos, R., Javier, F., & Atienza Alonso, D. A Multi-Lead ECG Classification Based on Random Projection Features. 37th IEEE International Conference on Acoustics, Speech, and Signal Processing (ICASSP 2012), pp. 625-628 New York: IEEE Press.
	In article
[155]	Bailón, R., Laguna, P., Mainardi, L., & Sornmo, L. (2007, August). Analysis of heart rate variability using time-varying frequency bands based on respiratory frequency. In Engineering in Medicine and Biology Society, 2007. EMBS 2007. 29th Annual International Conference of the IEEE (pp. 6674-6677). IEEE.
	In article
[156]	Seyd, P. A., Ahamed, V. T., Jacob, J., & Joseph, P. (2008). Time and frequency domain analysis of heart rate variability and their correlations in diabetes mellitus. International Journal of Biological and Life Sciences, 1.
	In article
[157]	Camm, A. J., Malik, M., Bigger, J. T., Breithardt, G., Cerutti, S., Cohen, R. J., & Singer, D. H. (1996). Heart rate variability: standards of measurement, physiological interpretation and clinical use. Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology. Circulation, 93(5), 1043-1065.
	In article		CrossRef
[158]	Kleiger, R. E., Miller, J. P., Bigger Jr, J. T., & Moss, A. J. (1987). Decreased heart rate variability and its association with increased mortality after acute myocardial infarction. The American journal of cardiology, 59(4), 256-262.
	In article		CrossRef
[159]	Kleiger, R. E., Stein, P. K., Bosner, M. S., & Rottman, J. N. (1992). Time domain measurements of heart rate variability. Cardiology clinics, 10(3), 487.
	In article
[160]	Nolan, J., Batin, P. D., Andrews, R., Lindsay, S. J., Brooksby, P., Mullen, M., & Fox, K. A. (1998). Prospective study of heart rate variability and mortality in chronic heart failure results of the United Kingdom heart failure evaluation and assessment of risk trial (UK-Heart). Circulation, 98(15), 1510-1516.
	In article		CrossRef
[161]	Corrales, M. M., de la Cruz Torres, B., Esquivel, A. G., Salazar, M. A. G., & Orellana, J. N. (2012). Normal values of heart rate variability at rest in a young, healthy and active Mexican population. Health, 4(7), 377-385.
	In article		CrossRef
[162]	Uehara, A., Kurata, C., Sugi, T., Mikami, T., & Shouda, S. (1999). Diabetic cardiac autonomic dysfunction: parasympathetic versus sympathetic. Annals of nuclear medicine, 13(2), 95-100.
	In article		CrossRef
[163]	Cugini, P., Bernardini, F., Cammarota, C., Cipriani, D., Curione, M., De Laurentis, T., & Napoli, A. (2001). Is a Reduced Entropy in Heart Rate Variability an Early Finding of Silent Cardiac Neurovegetative Dysautonomia in Type 2 Diabetes Mellitus? Journal of Clinical and Basic Cardiology, 4(4), 289-294.
	In article
[164]	Jelinek, H., Flynn, A., & Warner, P. (2004). Automated assessment of cardiovascular disease associated with diabetes in rural and remote health care practice. In The national SARRAH conference (pp. 1-7).
	In article
[165]	Tian, L., & Tompkins, W. J. (1997, November). Time domain based algorithm for detection of ventricular fibrillation. In Engineering in Medicine and Biology society, 1997. Proceedings of the 19th Annual International Conference of the IEEE (Vol. 1, pp. 374-377). IEEE.
	In article
[166]	Najarian, K., & Splinter, R. (2012). Biomedical signal and image processing. CRC press.
	In article
[167]	Manikandan, M. S., & Soman, K. P. (2012). A novel method for detecting R-peaks in electrocardiogram (ECG) signal. Biomedical Signal Processing and Control, 7(2), 118-128.
	In article		CrossRef
[168]	Farrell, T. G., Bashir, Y., Cripps, T., Malik, M., Poloniecki, J., Bennett, E. D., & Camm, A. J. (1991). Risk stratification for arrhythmic events in postinfarction patients based on heart rate variability, ambulatory electrocardiographic variables and the signal-averaged electrocardiogram.Journal of the American College of Cardiology, 18(3), 687-697.
	In article		CrossRef
[169]	Malik, M., Farrell, T., Cripps, T., & Camm, A. J. (1989). Heart rate variability in relation to prognosis after myocardial infarction: selection of optimal processing techniques. European heart journal, 10(12), 1060-1074.
	In article
[170]	Marwan, N., Wessel, N., Meyerfeldt, U., Schirdewan, A., & Kurths, J. (2002). Recurrence-plot-based measures of complexity and their application to heart-rate-variability data. Physical Review E, 66(2), 026702.
	In article		CrossRef
[171]	Majercak, I. (2002). The use of heart rate variability in cardiology. Bratislavske lekarske listy, 103(10), 368-377.
	In article
[172]	EK, F. S., & jakovljević, m. (2002). Heart rate variability-a shape analysis of Lorenz plots. cell. mol. biol. lett, 7(1), 160.
	In article
[173]	Sevšek, F., & Gomiščcek, G. (2004). Shape determination of attached fluctuating phospholipid vesicles. Computer methods and programs in biomedicine, 73(3), 189-194.
	In article		CrossRef
[174]	Woo, M. A., Stevenson, W. G., Moser, D. K., Trelease, R. B., & Harper, R. M. (1992). Patterns of beat-to-beat heart rate variability in advanced heart failure.American heart journal, 123(3), 704-710.
	In article		CrossRef
[175]	Kamen, P. W., Krum, H., & Tonkin, A. M. (1996). Poincare plot of heart rate variability allows quantitative display of parasympathetic nervous activity in humans. Clinical science, 91(Pt 2), 201-208.
	In article
[176]	Tulppo, M. P., Makikallio, T. H., Takala, T. E., Seppanen, T. H. H. V., & Huikuri, H. V. (1996). Quantitative beat-to-beat analysis of heart rate dynamics during exercise. American Journal of Physiology-Heart and Circulatory Physiology, 271(1), H244-H252.
	In article
[177]	Bogaert, C., Beckers, F., Ramaekers, D., & Aubert, A. E. (2001). Analysis of heart rate variability with correlation dimension method in a normal population and in heart transplant patients. Autonomic Neuroscience, 90(1), 142-147.
	In article		CrossRef
[178]	Theiler, J., Eubank, S., Longtin, A., Galdrikian, B., & Doyne Farmer, J. (1992). Testing for nonlinearity in time series: the method of surrogate data. Physica D: Nonlinear Phenomena, 58(1), 77-94.
	In article		CrossRef
[179]	Small, M., Yu, D., & Harrison, R. G. (2001). Surrogate test for pseudoperiodic time series data. Physical Review Letters, 87(18), 188101.
	In article		CrossRef
[180]	Cohen, M. E., Hudson, D. L., & Deedwania, P. C. (1996). Applying continuous chaotic modeling to cardiac signal analysis. Engineering in Medicine and Biology Magazine, IEEE, 15(5), 97-102.
	In article		CrossRef
[181]	Grassberger, P., & Procaccia, I. (1983). Measuring the strangeness of strange attractors. Physica D: Nonlinear Phenomena, 9(1), 189-208.
	In article		CrossRef
[182]	Judd, K. (1992). An improved estimator of dimension and some comments on providing confidence intervals. Physica D: Nonlinear Phenomena, 56(2), 216-228.
	In article		CrossRef
[183]	Rosenstein, M. T., Collins, J. J., & De Luca, C. J. (1993). A practical method for calculating largest Lyapunov exponents from small data sets. Physica D: Nonlinear Phenomena, 65(1), 117-134.
	In article		CrossRef
[184]	Uzun, I. S., Asyali, M. H., Celebi, G., & Pehlivan, M. (2001). Nonlinear analysis of heart rate variability. In Engineering in Medicine and Biology Society, 2001. Proceedings of the 23rd Annual International Conference of the IEEE (Vol. 2, pp. 1581-1584). IEEE.
	In article		CrossRef
[185]	Sun, Y., Chan, K. L., & Krishnan, S. M. (2000). Arrhythmia detection and recognition in ECG signals using nonlinear techniques. Ann. Biomed. Eng, 28, 32-37.
	In article
[186]	Eckmann, J. P., & World Scientific Series on Nonlinear Science Series A, 16, 365-404.
	In article
[187]	Gomes, M. E. D., Souza, A. V. P., Guimaraes, H. N., & Aguirre, L. A. (2000). Investigation of determinism in heart rate variability. Chaos: An Interdisciplinary Journal of Nonlinear Science, 10(2), 398-410.
	In article		CrossRef
[188]	Grassberger, P., & Procaccia, I. (1983). Characterization of strange attractors.Physical review letters, 50(5), 346-349.
	In article		CrossRef
[189]	Mandelbrot, B. B. (1983). The fractal geometry of nature/Revised and enlarged edition. New York, WH Freeman and Co., 1983, 495 p., 1.
	In article
[190]	Backes, A. R., & Bruno, O. M. (2008). Fractal and Multi-Scale Fractal Dimension analysis: a comparative study of Bouligand-Minkowski method.image, 11, 8.
	In article
[191]	Karim, N., Hasan, J. A., & Ali, S. S. (2011). Heart rate variability–a review. J. Basic Appl. Sci, 7, 71-77.
	In article
[192]	Schroeder, M. R. (2012). Fractals, chaos, power laws: Minutes from an infinite paradise. Courier Dover Publications.
	In article
[193]	C Tricot, C. (1995). Curves and fractal dimension. Springer.
	In article
[194]	Costa, L. D. F. D., & Cesar Jr, R. M. (2000). Shape analysis and classification: theory and practice. CRC Press, Inc..
	In article		CrossRef
[195]	Higuchi, T. (1988). Approach to an irregular time series on the basis of the fractal theory. Physica D: Nonlinear Phenomena, 31(2), 277-283.
	In article		CrossRef
[196]	Katz, M. J. (1988). Fractals and the analysis of waveforms. Computers in biology and medicine, 18(3), 145-156.
	In article		CrossRef
[197]	Kraft, R. (1995). Fractals and dimensions. Munich University of Technology.
	In article
[198]	Krstacic, G., Krstacic, A., Martinis, M., Vargovic, E., Knezevic, A., Smalcelj, A., & Smuc, T. (2002, September). Non-linear analysis of heart rate variability in patients with coronary heart disease. In Computers in Cardiology, 2002 (pp. 673-675). IEEE.
	In article
[199]	Acharya, R., Kannathal, N., & Krishnan, S. M. (2004). Comprehensive analysis of cardiac health using heart rate signals. Physiological measurement, 25(5), 1139.
	In article		CrossRef
[200]	Acharya, U. R., Kannathal, N., Sing, O. W., Ping, L. Y., & Chua, T. (2004). Heart rate analysis in normal subjects of various age groups. Biomed Eng Online, 3(1), 24.
	In article		CrossRef
[201]	Pincus, S. M., & Viscarello, R. R. (1992). Approximate entropy: a regularity measure for fetal heart rate analysis. Obstetrics & Gynecology, 79(2), 249-255.
	In article
[202]	Pincus, S. M. (1991). Approximate entropy as a measure of system complexity. Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
	In article		CrossRef
[203]	Goldberger, A. L., Mietus, J. E., Rigney, D. R., Wood, M. L., & Fortney, S. M. (1994). Effects of head-down bed rest on complex heart rate variability: response to LBNP testing. Journal of Applied Physiology, 77(6), 2863-2869.
	In article
[204]	Richman, J. S., & Moorman, J. R. (2000). Physiological time-series analysis using approximate entropy and sample entropy. American Journal of Physiology-Heart and Circulatory Physiology, 278(6), H2039-H2049.
	In article
[205]	Richman, J. S., & Moorman, J. R. (2000). Physiological time-series analysis using approximate entropy and sample entropy. American Journal of Physiology-Heart and Circulatory Physiology, 278(6), H2039-H2049.
	In article
[206]	Huikuri, H. V., Mäkikallio, T. H., Peng, C. K., Goldberger, A. L., Hintze, U., & Møller, M. (2000). Fractal correlation properties of RR interval dynamics and mortality in patients with depressed left ventricular function after an acute myocardial infarction. Circulation, 101(1), 47-53.
	In article		CrossRef
[207]	Peng, C. K., Havlin, S., Hausdorff, J. M., Mietus, J. E., Stanley, H. E., & Goldberger, A. L. (1995). Fractal mechanisms and heart rate dynamics: long-range correlations and their breakdown with disease. Journal of electrocardiology, 28, 59-65.
	In article		CrossRef
[208]	Mäkikallio, T. H., Høiber, S., Køber, L., Torp-Pedersen, C., Peng, C. K., Goldberger, A. L., & Huikuri, H. V. (1999). Fractal analysis of heart rate dynamics as a predictor of mortality in patients with depressed left ventricular function after acute myocardial infarction. The American journal of cardiology,83(6), 836-839.
	In article		CrossRef
[209]	Eckmann, J. P., Kamphorst, S. O., & Ruelle, D. (1987). Recurrence plots of dynamical systems. Europhys. Lett, 4(9), 973-977.
	In article		CrossRef
[210]	Chua, K. C., Chandran, V., Acharya, U. R., & Lim, C. M. (2008). Computer-based analysis of cardiac state using entropies, recurrence plots and Poincare geometry. Journal of medical engineering & technology, 32(4), 263-272.
	In article		CrossRef
[211]	Kobayashi, M., & Musha, T. (1982). 1/f fluctuation of heartbeat period.Biomedical Engineering, IEEE Transactions on, (6), 456-457.
	In article		CrossRef
[212]	Saul, J. P., Albrecht, P., Berger, R. D., & Cohen, R. J. (1987). Analysis of long term heart rate variability: methods, 1/f scaling and implications. Computers in cardiology, 14, 419-422.
	In article
[213]	Lombardi, F. (2000). Chaos theory, heart rate variability, and arrhythmic mortality. Circulation, 101(1), 8-10.
	In article		CrossRef
[214]	Bigger, J. T., Steinman, R. C., Rolnitzky, L. M., Fleiss, J. L., Albrecht, P., & Cohen, R. J. (1996). Power law behavior of RR-interval variability in healthy middle-aged persons, patients with recent acute myocardial infarction, and patients with heart transplants. Circulation, 93(12), 2142-2151.
	In article		CrossRef