Association of cardiovascular risk using non-linear heart rate variability measures with the framingham risk score in a rural population

Cardiovascular risk can be calculated using the Framingham cardiovascular disease (CVD) risk score and provides a risk stratification from mild to very high CVD risk percentage over 10 years. This equation represents a complex interaction between age, gender, cholesterol status, blood pressure, diabetes status, and smoking. Heart rate variability (HRV) is a measure of how the autonomic nervous system (ANS) modulates the heart rate. HRV measures are sensitive to age, gender, disease status such as diabetes and hypertension and processes leading to atherosclerosis. We investigated whether HRV measures are a suitable, simple, noninvasive alternative to differentiate between the four main Framingham associated CVD risk categories. In this study we applied the tone-entropy (T-E) algorithm and complex correlation measure (CCM) for analysis of HRV obtained from 20 min. ECG recordings and correlated the HRV score with the stratification results using the Framingham risk equation. Both entropy and CCM had significant analysis of variance (ANOVA) results [F(172, 3) = 9.51; <0.0001]. Bonferroni post hoc analysis indicated a significant difference between mild, high and very high cardiac risk groups applying tone-entropy (p < 0.01). CCM detected a difference in temporal dynamics of the RR intervals between the mild and very high CVD risk groups (p < 0.01). Our results indicate a good agreement between the T-E and CCM algorithm and the Framingham CVD risk score, suggesting that this algorithm may be of use for initial screening of cardiovascular risk as it is noninvasive, economical and easy to use in clinical practice.


INTRODUCTION
Identification of the risk of a cardiovascular disease (CVD) is an important attribute of preventative health care.

THE FRAMINGHAM RISK EQUATION AND CARDIOVASCULAR RISK
Multifactorial factors contribute to the increased risk of CVD. Previous health practice guidelines recommend that patients be treated with respect to their underlying coronary heart disease risk (Sheridan et al., 2003). Therefore accurate estimates are required to provide information on treatment strategies and timing of commencement of treatment. The Framingham risk equation has been validated in general populations. (D'Agostino et al., 2001) The Framingham risk equation has also been modified by various countries such as the Australian and New Zealand Cardiovascular Society (Jackson, 2000) as well as alternative risk equations proposed including the Coronary Risk Evaluation (SCORE) for fatal coronary heart disease or CVD and the Diabetes Epidemiology: Collaborative Analysis of Diagnostic Criteria in Europe (DECODE), which incorporates glucose tolerance status and fasting plasma glucose as well. (Conroy et al., 2003), (Balkau et al., 2004) History of CVD, physical inactivity, obesity, and left ventricular hypertrophy diagnosed by echocardiography are not included in current models of cardiovascular risk assessment (Jackson, 2000) as their individual predictive value is unclear, although obesity and left ventricular hypertrophy (LVH) are associated with increased cardiac pathology. There are suggestions that additional clinical measures may be of benefit such as biomarkers including C-reactive protein or D-dimer. (Charo et al., 1998;De Lemos, 2006;Nwose et al., 2007) However, these measures are invasive as they require blood samples to be sent to analytical laboratories and time consuming.

THE FRAMINGHAM SCORE AND AUTONOMIC NERVOUS SYSTEM
Assessment of the autonomic nervous system (ANS) function can provide information on the risk and presence of CVD. As such blood pressure, one component of the Framingham risk score (FRS) is a function of heart rate and peripheral www.frontiersin.org July 2013 | Volume 4 | Article 186 | 1 vascular resistance. In turn heart rate is a function of metabolism, which is also in part regulated by the ANS. Peripheral vascular resistance is equally regulated by the ANS, albeit only by the sympathetic branch. Increased cholesterol levels combined with oxidative stress due to for instance, increased free radicals leads to pathophysiological changes in the vascular endothelial cells and atherosclerosis, which also affects blood pressure. However increased blood pressure can also lead to atherosclerosis (Alexander, 1995, p. 9301). The metabolic changes due to oxidative stress are particularly a factor in diabetes disease progression, where an increase in reactive oxygen species activity associated with an increase in blood glucose levels leads to not only lipid peroxidation and atherosclerosis but also to damage to the kidneys, eye, gastrointestinal tract and the vascular system. Similar pathological changes to the peripheral circulation are due to smoking. The final parameter contained in the Framingham equation is age. Aging also affects heart rate variability (HRV), with increasing age associated with a decrease in HRV. (Pikkujämsä et al., 1999) These factors are the standard clinical parameters incorporated into the Framingham risk equation. Therefore we hypothesized that parameters that indicate characteristics of ANS function such as ANS influence on heart rhythm may be good indicators for CVD risk. The Framingham-based cardiovascular risk categories can be determined on line at http://www.mdcalc.com/framinghamcoronary-heart-disease-risk-score-si-units/ and other websites. The Framingham risk equation uses data on gender, age, cholesterol, blood pressure, diabetes status and smoking to determine the 10 year risk of belonging to one of four risk categories: mild, moderate, high and very high (Wilson et al., 1998).

HRV INDICES AND FRAMINGHAM RISK SCORE
HRV is an important physiological factor that is regulated by the sinoatrial node located in the right atrium of the heart and connected to the remainder of the heart by an elaborate network of cardiac cells that conduct the electrical impulse throughout the heart muscle. The endocrine system via a diverse set of hormones and the ANS also modulate HRV (Naschitz et al., 2005). HRV is also influenced by the same factors as incorporated in the Framingham risk equation including age, gender, cholesterol, blood pressure, diabetes status and smoking.
Standards of measurement and interpretation of HRV have been recommended by the Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology (Tfesc/Naspe, 1996). However the majority of studies investigating cardiovascular risk and HRV have used clinical data such as from the Framingham cohort to assess correlations between parameters measured as part of the FRS and HRV. Thus in the Framingham Heart Study, frequency and time domain analysis of HRV were found to be associated with a higher CVD mortality. (Tsuji et al., 1996) Studies by our laboratory and others have investigated the correlation between HRV and hypertension, cholesterol, gender, age, hyperglycemia and diabetes, smoking, body mass index, oxidative stress and heart disease. (Da Silva et al., 2004;Lampert et al., 2008;McLachlan et al., 2010;Tacoy et al., 2010;Thayer et al., 2010;Johnson et al., 2011;Khandoker et al., 2011;Kotecha et al., 2011;Fakhrzadeh et al., 2012;Matthews et al., 2012;Thiyagarajan et al., 2013) However, these studies did not stratify the cohorts according to the FRS but rather evaluated the magnitude of the HRV associated with presence or absence of the pathology being investigated. Hillebrand et al. (2013) recently reported that HRV is a useful marker for incidence of CVD in a population without any signs of CVD based on lowered time domain HRV values (Hillebrand et al., 2013).
Our search of the literature found one previous paper that investigated the correlation between HRV and FRS in a group of healthy adults (Yoo et al., 2011). The results of the study was based on dichotomizing the FRS into two groups, one with a CVD risk lower than 10% and the other a CVD risk higher than 10%. The results indicated that both time and frequency domain analyses were useful parameters for differentiating between high and low CVD risk.
Time and frequency domain analyses are sensitive to the length and non-stationarity inherent in the R-R series obtained from the ECG recordings. Several HRV analysis methods have been proposed in the last 10 years, which address some of the shortcomings associated with time and frequency domain analysis. The tone-entropy (T-E) algorithm is a method that is robust against data length, non-stationarity of the signal and also against respiratory influence and has been validated as a tool for analyzing HRV. The T-E algorithm has also been shown to correlate with experimentally induced HRV changes by either blocking the parasympathetic component using atropine or increasing the sympathetic influence on heart rate by head-up tilt (Oida et al., 1999). A clinical study analyzing HRV in a cohort of people with diabetes and cardiac autonomic neuropathy indicated a good correlation between cardiac autonomic progression and T-E as suggested by the work of Oida and colleagues (Khandoker et al., 2010). T-E is a combination of tone, which indicates the sympathovagal balance and entropy which indicates the overall activity of the ANS. Higher tone and lower entropy values indicate worsening ANS function and reduced HRV, which by its nature has a certain amount of variability and activity level (Khandoker et al., 2010).
Current HRV parameters used in HRV analysis are essentially static measures of a time signal. To obtain a better understanding of the relationship between autonomic control of the heart rate and the correlation to CVD risk a more dynamic measure of heart rate is required (Peng et al., 1995). The Poincaré Plot is a method proposed and utilized by Tulppo et al. (1996) for analysis of heart rate signals. The Poincaré Plot allows determining linear components of the inter-beat variability associated with short and long term correlations of the signal (Tulppo et al., 1998;Karmakar et al., 2009a). An extension of the Poincaré Plot to determine and measure the temporal dynamics over the recording interval was proposed by Karmakar et al. (2011). The complex correlation method (CCM) relies on computing the point-to-point variation of the signal with increasing lag (beat intervals) rather than the global description the Poincare plot provides based on n, and n + 1 beats.
This study investigated whether T-E and CCM are able to differentiate between ECG recordings obtained from four groups of patients categorized according to their FRS into mild, moderate, high and very CVD risk.

FRAMINGHAM RISK ASSESSMENT AND ECG SIGNALS
Data for the study was obtained from the Diabetes Screening Complications Initiative (DiScRi) clinic at Charles Sturt University. The study was approved by the Charles Sturt University Human Ethics Committee and undertaken between 2004 and 2008. Patient data was only included if it was complete for demographic as well as experimental variables. Under these criteria, 319 records and ECG traces were available. The Framingham cardiac risk score was determined using the protocol outlined by Jackson (2000) and determined automatically from data entered into the DiScRi ACCESS database (Pecoul and Jelinek, 2008). The variables used for determination of the cardiac risk score were age, sex, blood pressure, cholesterol [total and high density lipoprotein (HDL)], and smoking with diabetes status being a categorical variable. Cardiac risk as defined by the Framingham risk score is divided into very high (>20%), high (15-20%), moderate (10-15%) and mild (less than 10%) in accordance with the Australian and New Zealand Guidelines. The CVD risk score for each patient was used to categorize them into one of the four groups.
ECG signals were recorded and edited using the MLS310 HRV module (version 1.0, ADInstruments, Australia) included in the Chart software package. High frequency noise was removed with a 45 Hz low-pass filter and a 0.5 Hz high pass filter adjusted for wandering baseline. Ectopic beats were selected visually and deleted manually. Linear interpolation was used to replace ectopic beats that occur immediately before and after the ectopic interval. Intervals between successive R waves of the QRS complex (i.e., R-R intervals in seconds) were calculated using the algorithm developed by Pan and Tompkins. (Tompkins, 1993) The HRV analysis described in the following sections was performed on 1000 RR intervals.

TIME AND FREQUENCY DOMAIN ANALYSIS
We quantified several time domain HRV parameters: mean RR, standard deviation of normal RR data (SDNN) and the square root of the mean squared difference of the successive RR data (RMSSD). Spectral analysis was performed on linearly resampled (1 Hz) time series using Welch's method (Welch, 1967). The 256point fast Fourier transform was repeatedly computed with 50% overlap between adjacent segments. Then the spectral power of each segment was computed and averaged. Hanning window was applied to avoid spectral leakage. Subsequently, spectral powers in the low frequency (LF) band (0.04-0.15 Hz) and high frequency (HF) band (0.15-0.40 Hz) were obtained by integration (TFESC/NASPE). The normalized LF and HF powers were calculated by LF/(Total Power-VLF) and HF/(Total Power-VLF) respectively as per Task Force recommendation. (Tfesc/Naspe, 1996) The power in the very lower frequency (VLF) band was set at ≤ 0.04 Hz.

TONE-ENTROPY DETERMINATION
The methodology was described in detail in previous reports (Oida et al., 1999;Amano et al., 2005). In brief, acquired heart periods (RR intervals) are transformed into percentage index (PI) time series: where H(n) is a heart period time series, and n a serial number of heart beats shown. PI(n) therefore represents the number of instances each specific heart rate interval (RR) occurs in the time series. The tone is defined as a first order moment (arithmetic average) of this PI time series as: where N is a total number of PI terms. The tone, balance between accelerations (PI > 0) and inhibitions (PI < 0) of the heart, represents the sympatho-vagal balance faithfully as appreciated in all the previous studies (2002). The entropy is defined on PI probability distribution by using Shannon's formula: where [p(i)] is a probability that PI(n) has a value in the range, i < PI (n) < i + 1, i an integer. The entropy evaluates total acceleration-inhibition activities, or total heart period variations, in a familiar unit of bit. Entropy represents therefore the autonomic regulatory activity and tone the sympatho-vagal balance (Khandoker et al., 2010). Tone as a feature has its origin in the investigations in the last century of Rosenbluth and Simeone (Rosenblueth and Simeone, 1984). These authors investigated autonomic control of heart rate as an antagonistic interactive operation between acceleration and inhibition. Entropy evaluates HRV almost the same way as conventional second-order moments, for example, as standard deviation and is based on the Shannon entropy.

COMPLEX CORRELATION MEASURE
CCM measures the point-to-point variation of the signal rather than gross description of the Poincaré plot. It is computed in a windowed manner, which embeds the temporal information of the signal. A moving window of three consecutive points from the Poincaré plot are considered and the temporal variation of the points are measured. If three points are aligned on a line then the value of the variation is zero, which represents the linear alignment of the points. If the Poincaré plot is composed of N points then the temporal variation of the plot, the CCM, is composed of all overlapping three point windows and can be calculated as: where m, the lag, represents the number of consecutive points used from the Poincaré plot, A(i) represents area of the i-th triangle and C n is the normalizing constant which is defined as, C n = π * SD1 * SD2, represents the area of the fitted ellipse over Poincaré plot at lag-m. The length of major and minor axis of the www.frontiersin.org July 2013 | Volume 4 | Article 186 | 3 ellipse are 2SD1, 2SD2, where SD1, SD2 are the dispersion perpendicular to the line of identity (minor axis) and along the line of identity (major axis) respectively as proposed by Tulppo et al. (1996) The detail mathematical formulation of CCM is reported in our previous study Karmakar et al. (2009b).

STATISTICAL ANALYSIS
Results were expressed as means (±SD). One-way analysis of variance (ANOVA) and Bonferroni post hoc examination were carried out for comparisons among the four groups (mild, moderate, high, and very high) to evaluate whether statistically significant differences exist among the groups. In this study, the Lilliefors test was applied to test if the HRV features (tone, entropy, and CCM) comes from a distribution in the normal family. Results from the Lilliefors tests confirmed that HRV features' distributions were normal. Therefore, we decided to perform the statistical analysis with ANOVA. A value of p < 0.05 was considered significant for all examinations.

RESULTS
Three hundred and nineteen patients were recruited into the study with complete results from the DiScRi clinic including data required for determining the FRS. After exclusion of results due to noise in ECG signals or ectopic beats, 170 patients were included in the final analysis. Of these 85 participants were identified with mild (57%), 37 with moderate (18%), 25 with high (13%), and 23 (12%) cardiovascular risk. Demographic and cardiac risk factors according to the Framingham model are shown in Table 1 (Lloyd-Jones, 2010). Both entropy and CCM had significant ANOVA results [F (172, 3) = 9.51; <0.0001]. Bonferroni post hoc analysis indicated a significant difference between mild, high and very high cardiac risk groups with the entropy feature ( Table 2). CCM detected a difference in temporal dynamics of the RR intervals between the mild and very high CVD risk groups.
Post hoc analysis of SDNN for the CVD risk groups was also significant at 0.05 between mild and very high CVD risk. Figure 1 indicates that entropy steadily decreases with increased CVD risk and tone increases with increasing CVD risk (Figure 1).
CCM similarly decreases with advance of CVD risk (Figure 2). The increase in tone indicates a reduction in sympathovagal balance and the decrease in CCM indicates that this reduction in sympathovagal balance may be due to a reduction in parasympathetic influence on HRV. Entropy decreases with CVD risk progression, indicating a loss of ANS activity.
The results of the T-E analysis were then transformed into 2D space shown in Figure 3. Rectangles show the mean and standard errors of Tone and Entropy values of each group. The symbols   show the tone-entropy results for each participant with respect to CVD risk category. The important finding is that the distribution of results on the T-E plot roughly follow a curvilinear trend from top left to bottom right, showing as outlined in the text that there is an increase in tone and a decrease in entropy as one goes from the mild CVD risk group to the very high risk group (see Table 1).
Excitation-inhibition characterizes HRV and is the basis for calculating tone in the T-E algorithm. Figure 4 indicates the histograms based on the probability of excitations and inhibitions occurring in any heart rate recording. The histograms indicate that the relationship between excitations and inhibitions are clearly different between the four CVD risk groups (Figure 4).
The Mild CVD risk histogram has a wide spread of both acceleration and inhibitions in heart rate with inhibitions on the left being slightly more common. The moderate CVD risk group (B) indicates that there is a loss of both inhibitions and excitations with the remaining beat intervals occurring at a higher percentage. With high CVD risk (C) the variance (spread of the histogram) of inhibition and excitations are still reduced but in addition the percentage of the different beat lengths occurring is also reduced. The histogram for very high CVD risk (D) is characterized by a further reduction in variance or spread of the histogram and inhibitions and excitations all occur less often but with inhibitions dominating the heart rate variance.

DISCUSSION
CVD risk as determined by the FRS requires data for gender and age, diabetes status, blood pressure, cholesterol level as well as smoking status. From a clinical perspective determining the final Framingham risk score is based on analysis of blood samples, which is an invasive procedure and requires samples to be sent to a testing laboratory. Simpler methods that are non-invasive but are correlated with the FRS may provide an alternative where the FRS cannot be calculated or provide additional clinical information. Our study focused on whether the magnitude of CVD risk as determined by FRS is correlated with HRV measured by the T-E and CCM algorithms.
Current measures of the autonomic balance of heart rate regulation include time and frequency domain as well as nonlinear measures (Goldberger and West, 1987;Pincus, 1991;Alam et al., 2009;Karmakar et al., 2009a;Voss et al., 2009;Cysarz et al., 2011). Time and frequency domain measures provide a global characterization of the HRV over the recording interval, which can range between 2 min to 24 h (Dekker et al., 2000;Rennie et al., 2003;Kiviniemi et al., 2011). To obtain a better understanding of the optimal ANS control of the heart rate, more dynamic measures of heart rate such as nonlinear methods are required. Quantitative Poincaré plot analysis was used to assess the changes in CCM of HRV signals during parasympathetic blockade was discussed in our previous study (Karmakar et al., 2011). The lowest value of CCM has also been found during atropine infusion which reduced the parasympathetic activity and reduces instantaneous changes in HRV signal. On the other hand, it was also found to be increased with increase in parasympathetic activity during administration of low-dose scopolamine. Variability (increasing or decreasing) in the temporal structure of the Poincaré plot (measured as CCM) reflects the change in parasympathetic activity harmoniously. An advantage of nonlinear methods is that they are not sensitive to the non-stationarity of ECG signals (Peng et al., 1995). T-E has a number of advantages over the conventional methods of HRV analysis in that it is robust against non-stationarity and respiratory influence and not sensitive to differences in length of recordings. Similarly CCM can be applied to short recordings and provides temporal/dynamic information about HRV not obtained by other methods (Karmakar et al., 2011). This latter point makes tone, entropy and CCM especially useful for clinical investigation where short time recordings are more likely to be used. T-E provides also an alternative visualization of the functionality of the ANS with respect to CVD risk as shown in Figure 2 and CCM is a dynamic rather than a static measure of HRV. CCM therefore provides information on the magnitude of the beat-to-beat influence from a lag of one beat to ten beats, which has been shown to be important in HRV and is subject to parasympathetic influence (Lerma et al., 2003), necessitating a means of exploring this influence. CCM, rather than analyzing beat-to-beat variation in fact analyses the relationship between an initial beat and any beat downstream from n to n + 1. In this work we have used n = 3.
Our results show that mild risk of less than 10%, moderate (10-15%), high (15-20%) and very high risk (>20%), in accordance with the Australian and New Zealand Guidelines, is associated with decreasing entropy and CCM and an increasing tone. The lack of a significant result between mild and moderate may be in part due to the minor difference in CVD risk % between the current groups. Previous work by Yoo et al. (2011) divided the cohort into those with less or greater than 10% CVD risk and showed a significant relationship between the FRS and time and frequency domain parameters (Yoo et al., 2011). Our analysis using the four main CVD risk groups of the FRS showed only a significant result for SDNN, the global time domain measure. No frequency domain measures were significant. We propose that entropy and CCM are more suitable measures to identify CVD risk in terms of the FRS due to their robust nature against non-stationarity and measuring total activity and the dynamic nature of the heart rate interval changes over the time of the recording.
Analysis of the distribution of accelerations and inhibitions as shown in Figure 4 are related to the tone-entropy results. The mild CVD risk group has the widest histogram and the lowest peak, which indicates the highest entropy and lowest tone and therefore a higher HRV. Mild CVD risk has the most active sympatho-vagal modulation of heart rate and optimum sympatho-vagal balance. Narrower histograms as seen in Figure 4 with higher peaks indicate loss of ANS activity by a decrease in entropy and loss of sympatho-vagal balance indicated by the increase in tone. The trend of increasing tone and decreasing entropy continues as CVD risk increases. Moderate, high and very high CVD risk histograms are narrower with higher tone values and therefore a lower sympatho-vagal balance. This is especially prominent in the moderate CVD risk group. T-E analysis thus indicates subtle physiological changes in heart rate modulation not previously identified. These results correlate with previous physiological and pharmacological experiments (see Figures 2, 8;Oida et al., 1997).
The results for T-E and CCM analysis correlate with the CVD risk as determined by the FRS. HRV changes in response to various factors including ANS activity. In turn ANS activity is dependent on age, smoking, diabetes and blood pressure, all being variables of the Framingham CVD risk equation. TE and CCM are able to differentiate the mild to high and very high CVD risk groups as defined by the Framingham CVD risk score. As such HRV parameters may be seen as a single measure that is non-invasive and provides an alternative for identifying CVD risk.

Frontiers in Physiology | Computational Physiology and Medicine
July 2013 | Volume 4 | Article 186 | 6 However, we are not suggesting a new or better classification of risk compared to the FRS but rather that T-E and CCM correlate well with the accepted Framingham CVD risk score and the associated categories and provide a better stratification compared to the time and frequency domain parameters as can be seen in Table 2. In essence we are proposing that the entropy and CCM measure are able to stratify CVD risk groups as defined by FRS to mild, high and very high CVR risk groups. This has the advantage that a noninvasive, easy to perform method can be used in the first instance to identify individuals with increased risk of CVD. This becomes important when blood samples for instance cannot be taken for ethical, ethnic or other reasons. In addition HRV has been shown to be a sensitive marker for sudden cardiac death and also for inflammation as an atherosclerosis process as well as risk of arrhythmic events.
A future study of ours is to determine tone, entropy and CCM in a larger, younger cohort as most people below 40 are unlikely to have a very high risk of CVD and therefore may be missed in traditional CVD risk assessment despite other factors included in the Framingham risk equation being abnormal. In addition we are investigating the individual influence of the Framingham CVD risk factors on HRV either separately or in combination. This will provide a basis for reviewing younger patients with raised blood pressure and raised cholesterol levels but in a low CVD risk category due to age.