Hypothetical Control of Heart Rate Variability.

In the last three decades, the analysis of heart rate variability by nonlinear methods demonstrated the complexity of cardiovascular regulation. Additionally to the observations of periodic heart rate regulation by the autonomic nervous system, the long-term statistics of the heart rate has been determined to reminisce a tempered Lévy process. A number of heuristic arguments have previously been made to support a tempering conjecture, using exponentially truncated waiting times for the time intervals between heart beats. Herein we use the fractional probability calculus to frame our arguments and to parameterize the control process that tempers the Lévy process through a collective-induced potential. We also determine that the hypothesis of a self-induced nonlinear potential control resulting in such a tempered Lévy process is consistent with the hypothesis of disease being the loss of physiologic complexity made over 25 years ago.


Data
All data used were from the MIT-BIH RR Interval Database (Goldberger et al., 2000).The Normal Sinus Rhythm RR Interval Database was used as the healthy control group. This database included 54 long-term RR interval recordings of subjects in normal sinus rhythm. The Congestive Heart Failure RR Interval Database (NYHA class I, II and III) and the BIDMC Congestive Heart Failure Database (NYHA class III and IV) were used as the pathological groups. Those databases included 29 and 15 long-term RR interval recordings, respectively. The original ECG signals in all cases were digitized at 128 Hz and the beat annotations were obtained by automated analysis with manual review and correction. Inclusion of NYHA patients of class III in both CHF groups is a consequence of analyzing publicly available datasets and not an intentional choice. Available RR time series are anonymized, the information about the NYHA class of individual patients is not available and more precise separation is not possible. As a result both CHF groups include patients of NYHA class III.

Linear analysis
The mean heart rate, standard deviation of N-N intervals (SDNN) and root mean square of successive differences of N-N intervals (rMSSD) were used as time-domain measures of HRV (Shaffer and Ginsberg, 2017). The power spectrum densities were estimated by Welch's averaged periodogram method applied to R-R time series resampled at 2 Hz with the use of a cubic spine interpolation. Very low frequency power (VLF, 0.0033 to 0.04 Hz), low frequency (LF, 0.05-0.15 Hz) and high frequency power (HF, 0.15-0.40 Hz) were calculated from the entire 24-hour segments.

Nonlinear analysis
Detrended Fluctuation Analysis (DFA) quantifies fractal-like correlations in time series (Peng et al., 1995b). The scaling exponent DFA α indicates the scaling of root-mean-square fluctuations of the integrated and detrended HRV data as a function of observation window of increasing size (Peng et al., 1995a). Parameter α 1 captures brief fluctuations (4-11 beats), while value of α 2 describes long-term scaling (>11 beats).

Non-Gaussianity index
The analysis of non-Gaussianity of HRV starts from integrating and detrending (third order polynomial for the trend was used herein) R-R interval time series (Kiyono et al., 2004(Kiyono et al., , 2006. After the detrending procedure heart rate increments at scale s are defined as is a deviation of the HRV from the polynomial trend. Then the PDF of differences ∆ s B(t) normalized by the SD is evaluated. Finally, the non-Gaussianity index λ s is estimated as where |∆ s B| q denotes an estimated value of the q-th order absolute moment of [∆ s B]. Following (Hayano et al., 2011) we calculate λ s based on the 1/4th order moment (q = 0.25) and we set the scale s at 25 s. The parameter λ s was originally introduced as a parameter of log-normal cascade model describing non-Gaussian distributions in the study of intermittency of hydrodynamic turbulence (Castaing et al., 1990): The value of λ s close to zero corresponds to a PDF close to a Gaussian distribution, while a larger value of λ s means that the observed PDF has fatter tails and a sharper peak in comparison with the Gaussian distribution.

Results
The values (mean ± SD) of listed HRV indices for the group of healthy subjects and two groups of patients with CHF are presented in Table 5.1. All indices overlap with values reported in literature for comparable datasets (Shaffer and Ginsberg, 2017;Beckers et al., 2006;Hayano et al., 2011). In line with previous studies (Hayano et al., 2011), congestive heart failure leads a decrease in values of numerous HRV indices. We observe that the power of the HRV in LF and HF frequency bands, LF/HF ratio or scaling coefficients estimated through DFA are affected by the pathological condition. This corresponds to a change in sympatho-vagal balance, where the CHF leads to a dominance of parasympathetic branch of the ANS. Reduced DFA scaling coefficients reflect diminished correlations present in HRV time series, reflecting reduced complexity of cardiovascular control mechanisms.

FFPE
The conditional PDF P (I, t) ≡ P (I, t |I 0 , t = 0) is the average of the phase space density function ρ (I, t) over an ensemble of realizations of ξ (t) : The phase space density function is defined as where I(t) is the solution to the Langevin equation for a particular realization of dξ (t), in this case given by Eq.(6). The formal expression for the PDF is then given by The evolution of phase space density function is determined by the Liouville equation, or continuity equation The function · I(t) in the continuity equation is defined by differentiating the solution Eq.(6) with respect to time: where for notational ease we introduce the function to replace the integral term. The Fourier transform of ρ (ι, t) is so that the continuity equation can be written as To solve Eq.(30) we introduce the deterministic operator the stochastic operator and rewrite Eq.(30) as ∂φ (ω, t) ∂t The solution to this equation is given by and the mathematical details are given by West and Seshadri (West and Seshadri, 1982).
The function φ (ω, t) can be averaged over an ensemble of realizations of dξ (t) to yield the characteristic function P (ω, t) : which is the Fourier transform of the PDF. For an arbitrary analytic function g(τ ), Doob (Doob, 1942) showed that for a differential Lévy process dB(τ ) of strength σ 2 , with an index α and which is delta correlated in time: Consequently, if we assume the stochastic driver dξ (t) is such a Lévy process, after further analysis (West and Seshadri, 1982), the equation of motion for the characteristic function is determined to be ∂ P (ω, t) ∂t which is solved (West and Seshadri, 1982) to obtain where the time-dependent 'width' of the distribution is The PDF generated by the linear Langevin equation is therefore a symmetric Lévy PDF.

Steady-state solution
Chechkin et al. (Chechkin et al., 2005) argue that the behavior of the steady-state PDF can be determined from an analysis of the contributions to the integral in Eq.(19) in the vicinity of the pole. In particular, the dominant contribution to the integral arises from the region to the left of the pole I = I and consequently the operator may be approximated by λ 2n I 2n+1 P ss (I) − C α d dI The integral term may be expressed to the lowest-order approximation