The Fractal-like Complexity of Heart Rate Variability beyond Neurotransmitters and Autonomic Receptors: Signaling Intrinsic to Sinoatrial Node Pacemaker Cells

The heart rate and rhythm are controlled by complex chaotic neural, chemical and hormonal networks which are not strictly regular, but exhibit fluctuations across multiple time scales. A careful assessment of the Heart Rate Variability (HRV) offers clues to this complexity. A reduction in HRV, specifically in advanced age, is associated with increase in morbidity and mortality. Mechanisms that induce this decrease, however, have not been fully elucidated. The classical literature characterizes changes in HRV as a result of changes in the balance of competing influences of the sympathetic and parasympathetic autonomic impulses delivered to the heart. It has now become clear, however, that the heart rate and HRV are also determined by intrinsic properties of the pacemaker cells that comprise sinoatrial node, and that these properties respond to autonomic receptor stimulation in a non-linear mode. That HRV is determined by both the intrinsic properties of pacemaker cells in the sinoatrial node and the competing influences of the two branches of the autonomic neural input to the cells requires an expansion of our perspective about mechanisms that govern HRV in the normal heart, and how HRV changes with aging in health and in heart diseases.


Introduction
The healthy heart rapidly changes its rate in response to a broad range of integrated input signals that arise from diverse stimuli. Changes in heart rate are governed by robust control mechanisms involving regulatory networks of neural, chemical and hormonal pathways interacting through coupled feedback mechanisms that operate over a wide range of time periods, i.e., ultradian (from milliseconds, minutes to hours) to circadian (24 hours) and infradian (beyond days) [1]. Because these complex networks are chaotic, the periods of rhythms they regulate are never in a steady state, but continuously shift from one period to another. This chaotic complexity is inherent to the rate and rhythm of the healthy heart, and is the reason why the heart rhythm is not strictly regular, but exhibits fluctuations. Careful assessment of the Heart Rate Variability (HRV) is required to probe its complexity. Specifically, decoding the electrogram (ECG), permits access to complex harmonics imbedded in the heart rhythm. Loss of this complexity is manifest as a reduction in HRV, and is related to increases in both morbidity and mortality [2].
Changes in HRV have been mainly interpreted as the balance between the competing influences of the two branches of the autonomic nervous system: the sympathetic tone that accelerates the heart rate, and the parasympathetic (vagal) stimulation that slows the heart rate ( Figure 1A). Based on this logic, the balance of flux between the two arms of the autonomic system determines both the heart rate and HRV. An imbalance between these two afferent systems has been thought to be the major mechanism for reduction in HRV associated with aging and cardiac diseases. Studies near the turn of this century, including recent discoveries [3][4][5][6][7][8], however, call attention to the importance of intrinsic properties of pacemaker cells within the sinoatrial node, the heart's primary pacemaker, and their responsiveness to autonomic signaling, as an important determinant of heart rate and HRV. Novel perspectives have emerged regarding the intrinsic properties of SinoAtrial Node Cells (SANC) that determine the normal automaticity of spontaneous Action Potentials (APs) that emanate from the sinoatrial node ( Figure  1B), and how signaling through autonomic receptors links to these intrinsic signaling pathways to alter normal automaticity ( Figure 1C).

Decoding the Complexity of HRV
If we place our finger on the radial artery to determine our heart rate, we are struck by the apparent regulatory of intervals at which the heart beat occur. In reality, however, a plot of the healthy heart beatto-beat intervals determined from the ECG demonstrates substantial variation, even at rest in the absence of external stimuli (e.g. change in the nervous system input) [9]. Time-series analysis is the simplest means of evaluating HRV and to identify variation in beat-to-beat intervals over time. Common methods utilized are: standard deviation (SDNN) to measure global variation, the square root of mean squared differences of consecutive intervals (RMSSD) to measure short-time variation, and standard deviation of 5-min averages (SDANN) to evaluate long-term variation. Identical means and standard deviations in the time-domain measures, however, can be generated by different C a rd iovascu la r P h ar macol o g y : O pe n Acc es s underling rhythms in the time series of heart-rate intervals. Moreover, differences in the time series of heart-rate intervals lie not only in the distribution of the interbeat variations around the mean heart rate, but also in their time ordering: different dynamic patterns can occur at different time scales that correspond to the actual time scales of the underling physiological process [10]. A nonlinear form and function has been introduced to examine these dynamic patterns. Fractal analysis of the heart beats is one of the method employed to decode the secrets hidden in HRV [11][12][13][14].
The idea of examining fractal-like behavior of a system to demonstrate its complexity was introduced by Mandelbort [15]. A system exhibits fractal-like temporal behavior if any time-based measurement segment within the system resembles the system behavior over the full time scale [15]. Fractal analyses of HRV have indicated that its control mechanisms are not linear, but typically exhibit chaotic behavior over a wide range of time scales that are slower than the heart rate per second. Further, application of fractal analysis to the physiology of HRV clearly demonstrates the presence of a high degree of chaos and entropy within the normal heart rhythm [2]. Moreover, the HRV components occur over a wide range of frequencies and exhibit inverse power-law scaling (1/f-like distribution) suggesting that fractal-like behavior regulates the HRV and therefore the heart rate [11]. Consequently, fractal-like behavior analyses can be used as a tool to monitor heart pathophysiological changes in HRV that occur in a variety of settings.
Two quantitative methods have been employed to examine the fractal-like behavior of HRV. The first approach includes detection of the frequency harmonics embedded within the heart's rhythm by performing Fourier analysis on the inter heart beat intervals. For example, the frequency regimes embedded within the human heart rate include: 1) a high frequency band, between 0.15 and 0.4 Hz (between 2 to 6 heart cycles); 2) a low frequency band, between 0.04 and 0.15 Hz (between 6 to 25 heart cycles); and 3) a very low frequency band, below 0.04 Hz (more than 25 heart cycles) [16]. The slope of the line relating log frequency to log of the power spectrum density (derived from Fourier analysis of heart rate intervals) is the fractal scaling exponent, β ( Figure 2A). The second method is the detrended fluctuation analysis (DFA), in which the degree of correlation among time scales embedded within the heart beat intervals is characterized [17]. In this method, the self-similarity of frequency regimes buried within the heart beat intervals assessed from ECG analysis is assumed to be bi-fractal, and is described by short-and long-term exponents, α 1 and α 2 , respectively ( Figure 2B). Note that if β or α is less than 0.5 the system dynamic is described as white noise (i.e., a random signal with a constant power spectral), and if β or α equals 1.5, the system dynamic is described as Figure 1: The crosstalk between the autonomic nervous system and the sinoatrial node system. (A) In vivo, neural impulses via parasympathetic and sympathetic nerve and epinephrine (not shown) influence the heart rate and rhythm. (B) Intrinsic spontaneous AP firing by sinoatrial node cells, even in the absence of autonomic impulses, is driven by basal Ca 2+ -calmodulin activation of adenylyl cyclases (AC), which produce cAMP-PKA-dependent phosphorylation signaling. Ca 2+ -calmodulin activates, in parallel to AC, calmodulin-dependent kinase II (CaMKII) phosphorylation signaling. cAMP positively shifts the f-channel activation curve. Phosphodiesterases (PDE) degrade cAMP production, while protein phosphatase (PPT) degrade phosphorylation activity. PKA and CaMKII phosphorylate sarcoplasmic reticulum (SR) Ca 2+ cycling proteins (Ryanodine receptor (RyR), Phospholamban (PLB), which binds to and inhibit the sarcoplasmic reticulum Ca 2+ -ATPase (SERCA)), surface membrane ion channel proteins, ATP production mechanisms in the mitochondria and sarcomere. (C) The neurotransmitters noradrenaline (NE) and acetylcholine (ACh) released from sympathetic or parasympathetic nerve terminals bind to β-adrenergic receptors (β-AR) or cholinergic receptors (Ch-R), respectively. Autonomic receptor signaling couples to G-proteins and leads to modulation of the same coupled-clock molecules that drive basal automaticity of SANC.
Brownian noise (i.e., an integral of white noise that has more energy at lower frequencies). That a reduced fractal slope can predict cardiac death [11] indicates that fractal correlation is a powerful tool to probe subtle alterations in heart rate dynamics and complexity.
Fractal scaling of HRV components is a universal organizing principle that unites diverse heart rates among hearts of diverse sizes in diverse members of the mammalian kingdom. When the log of power spectrum density of components that comprise the HRV for different mammals (dog, rat, hamster and human) is plotted against the log of frequency, a unique harmonically modulated inverse powerlaw dependence is observed [18,19]. Such log-log plots are represented by power spectra having a broad band of frequencies with a long low-amplitude, high-frequency tail and fractal scaling exponent, β approaching 1 (Figure 2, Table 1).
Fractal analyses have also been employed to decode changes in the heart's rhythm in response to pharmacological perturbations of β adrenergic receptor or cholinergic receptor signaling. Inhibition of the sympathetic or parasympathetic tone (vagal) increases or decreases the fractal-like behavior slope, respectively (Table 1). Therefore, parasympathetic receptor stimulation shifts the HRV fractal dynamics toward a white noise while the cardiac sympathetic receptor stimulation shifts the fractal dynamics toward a Brownian noise component. The conclusion that a shift in the ratio of sympathetic vs. parasympathetic tone toward sympathetic decreases HRV is supported by the observation that HRV is nearly suppressed in transgenic mice in which the β 1 -adrenergic receptors are over expressed in the atria [20].
We contend that, fractal-like complexity of HRV goes beyond neuro-transmitters and autonomic receptors, and must include both the intrinsic properties within SANC and how these cell properties become organized within sinoatrial node tissue.

Intrinsic Sinoatrial Node Cell Properties: A New Vista to Reveal the Mystery Embedded within the HRV
New experimental and numerical data support the hypothesis that, when isolated from the autonomic nervous system, the SANC AP firing rate is regulated by integrated functions within a system of two coupled clocks [21] (Figure 3): the sarcoplasmic reticulum is a Ca 2+ clock that spontaneously generates roughly periodic local diastolic Ca 2+ releases (LCRs) that activate an inward Na + -Ca 2+ exchanger current, which accelerates the rate of diastolic membrane depolarization. Na + -Ca 2+ exchange current and the f-channel current, another member of the ensemble of sarcolemmal electrogenic molecules ("membrane clock"), concurrently drive the diastolic membrane depolarization to ignite the next AP. The two clocks are tightly coupled by common chemical axes: Ca 2+ -calmodulin activated adenylyl cyclases (neuronal types 1 and 8) generate the second messenger cAMP, which in turn activates Protein Kinase A (PKA), and Ca 2+ /calmodulin-dependent protein kinase II (CaMKII) ( Figure 1B). Both kinases phosphorylate proteins of the "membrane clock" and of the "Ca 2+ clock"; cAMP positively shifts the f-channel activation curve [22]. Moreover, cAMP/PKA and CaMKII-dependent phosphorylation modulate mitochondrial ATP production [23]. The rate at which this intrinsic coupled-clock system ticks becomes further accelerated, or reduced, respectively, by adrenergic or cholinergic neurotransmitter receptor stimulation. Specifically, neurotransmitter-activated receptors signal via many of the same molecules that regulate basic automaticity. Autonomic receptors are coupled to G-proteins and to adenylyl cyclases (likely type 5 or 6) or to guanylyl cyclase, and this leads to activation or suppression of protein kinases that modulate the extent of protein phosphorylation of intracellular Ca 2+ and surface membrane ion channel proteins that drives normal pacemaker clock activity ( Figure 1C). Each one of these biochemical feedback mechanisms can operate over a wide range of time scales, and this heterogeneity, therefore, can induce fluctuations within different frequency regimes embedded within the heart rhythm. For example, in heart failure patients, ivabradine, a pharmacological I f blocker ( Figure 1C), induces a heart rate reduction and an increase in indices of heart rate variability [24]. But the resultant ivabraine-induced changes in heart rate and rhythm embody contributions of changes not only in I f , but also of changes in Ca 2+ cycling within SANC, because the membrane and Ca 2+ clocks crosstalk (in response to any disturbance signal that directly perturbs either clock entrains a change in function of the other clock) [21].
The potential contribution of the intrinsic, coupled-clock signaling within SANC to the complexity of the heart rhythm mechanisms can be probed by deletion of autonomic input by complete cardiac denervation. Such denervation of the heart in situ dramatically reduces the HRV and therefore increases the fractal complexity of the heart rate in dogs [4]. However, even in the presence of this apparent complete denervation, the fractal-like behavior is still present [3], although the fractal-like slope (β) is reduced. The continued presence of fractal-like behavior in the denervated heart suggests that the fractal complexity of the HRV is modulated not only by the autonomic nervous system, but also by intrinsic properties of the sinoatrial node and the cells within it (and also likely by epinephrine released from adrenal cortex). Interestingly, in isolated hearts (i.e., when the heart is completely detached from hormonal and neural input), there is a non-linear relationship between the HRV and heart beating rate in numerous species [8]. These findings are further supported by experiments in single SANC that has been isolated from the sinoatrial node after its removal from the heart. In isolated rabbit SANC, the AP firing rate exhibits a non-linear dependence on cholinergic receptor stimulation (by graded concentrations of acetylcholine) [5,6]. Moreover, there is a tight, but non-linear, relationship between the spontaneous AP firing rate and AP firing rate variation, indicating that changes in spontaneous AP rate and rhythm in single isolated SANC are not simply surrogates of autonomic receptor stimulation [5,6]. In other terms, a non-linear transformation of autonomic receptor signaling occurs within SANC, and this transformation is a determinant of changes in AP firing rate and rhythm initiated by receptor stimulation. Moreover, in isolated Figure 2: Analytic methods of fractal-like behavior. Schematic illustration of (A) the linear function derived from a log-log plot (slope β) of the power spectral analysis (PSD, calculated from Fourier analysis of heart rate intervals) vs. frequency, (B) the short and long-term linear relationship (slope α 1 and α 2 , respectively) of the amplitude of detrended fluctuations (F(n)), calculated from heart rate intervals vs. the block size n (in beats), on a log-log scale. The red dots illustrate experimental data. rabbit SANC, beat-to-beat variations in the spontaneous AP firing rate are directly correlated with average variations in the period of the roughly periodic local Ca 2+ releases that emanate from the sarcoplasmic reticulum via ryanodine receptors [7].
That HRV is determined by both the intrinsic properties of cells in the sinoatrial node and the competing influences of the two branches of the autonomic neural input requires a broader interpretation regarding mechanisms that govern HRV in the normal heart, and how it changes with aging in health and in pathological conditions. Genetic manipulation of protein phosphorylation of intracellular Ca 2+ and surface membrane ion channel proteins of the coupledclock system ( Figure 1B), i.e., I f [25], ryanodine receptors [26], and α 1 subunit of Cav 1.3 [27] induce change in both heart rate and rhythm. In   addition, pharmacological interventions that impact on intrinsic SANC mechanisms (Figure 1), e.g., intracellular Ca 2+ [28], can cause the AP firing rate to become irregular. It is also known that denervation of the autonomic nervous system decreases HRV [4]. These results indicate that under normal conditions the intrinsic properties of the sinoatrial node, per second; contribute Brownian noise-like components to the HRV fractal dynamics. In congestive heart failure syndrome, the HRV is reduced and exhibits a reduced fractal-like slope [29]. That congestive heart failure is associated with increased resting levels of sympathetic tone and with decreased parasympathetic tone can only partially explain the reduced HRV and fractal-like slope [29] because congestive heart failure is also associated with increased heart rate [30]. In other words, as demonstrated by Zaza and Lombardi [6], when the β adrenergic receptor is stimulated, the spontaneous AP firing rate in isolated SANC increases (i.e., the autonomic nervous system is not intact) and the HRV become reduced. Moreover, congestive heart failure is associated with intracellular structural remodeling, change in compartmentation of cAMP activity and alteration in neurotransmitter-activated receptor signaling via phosphorylation mechanisms [31]. Therefore, changes in the properties of cells that comprise the sinoatrial node are also likely implicated in the reduction of HRV and fractal-like slope observed in heart failure syndromes. Similarly, a reduction in HRV and fractal-like slope in advanced age [29] may not only be associated with increased resting levels of sympathetic tone and with decreased parasympathetic tone, and an exaggerated increase in plasma catecholamines in response to stress in the context of a reduced adrenergic receptor stimulation [32], but also to changes in intrinsic intracellular coupled-clock mechanisms that accompany advanced age. In fact, aging is associated with deficient intrinsic cAMP-PKA-Ca 2+ signaling [33]. Therefore, in this case, changes that occur in the properties of the cells in the sinoatrial node likely contribute to age-associated changes in both HRV and fractal-like behavior.

Future Directions
In summary, the heart exhibits fractal-like behavior that is determined by intrinsic properties of sinoatrial node pacemaker cells in addition to the competing influences of the two branches of the autonomic neural input to the heart. Quantification of the fractal-like behavior at different complexity levels (intact heart, denervated heart, isolated sinoatrial node, phosphorylation networks within isolated SANC, and organization of SANC within the sinoatrial node) is required for a complete understanding of the mechanisms that regulate HRV and induce changes in HRV during aging in health and in the presence of cardiovascular diseases.