Figures
Abstract
The proposed approach evaluates complexity of the cardiovascular control and causality among cardiovascular regulatory mechanisms from spontaneous variability of heart period (HP), systolic arterial pressure (SAP) and respiration (RESP). It relies on construction of a multivariate embedding space, optimization of the embedding dimension and a procedure allowing the selection of the components most suitable to form the multivariate embedding space. Moreover, it allows the comparison between linear model-based (MB) and nonlinear model-free (MF) techniques and between MF approaches exploiting local predictability (LP) and conditional entropy (CE). The framework was applied to study age-related modifications of complexity and causality in healthy humans in supine resting (REST) and during standing (STAND). We found that: 1) MF approaches are more efficient than the MB method when nonlinear components are present, while the reverse situation holds in presence of high dimensional embedding spaces; 2) the CE method is the least powerful in detecting age-related trends; 3) the association of HP complexity on age suggests an impairment of cardiac regulation and response to STAND; 4) the relation of SAP complexity on age indicates a gradual increase of sympathetic activity and a reduced responsiveness of vasomotor control to STAND; 5) the association from SAP to HP on age during STAND reveals a progressive inefficiency of baroreflex; 6) the reduced connection from HP to SAP with age might be linked to the progressive exploitation of Frank-Starling mechanism at REST and to the progressive increase of peripheral resistances during STAND; 7) at REST the diminished association from RESP to HP with age suggests a vagal withdrawal and a gradual uncoupling between respiratory activity and heart; 8) the weakened connection from RESP to SAP with age might be related to the progressive increase of left ventricular thickness and vascular stiffness and to the gradual decrease of respiratory sinus arrhythmia.
Citation: Porta A, Faes L, Bari V, Marchi A, Bassani T, Nollo G, et al. (2014) Effect of Age on Complexity and Causality of the Cardiovascular Control: Comparison between Model-Based and Model-Free Approaches. PLoS ONE 9(2): e89463. https://doi.org/10.1371/journal.pone.0089463
Editor: Mathias Baumert, University of Adelaide, Australia
Received: November 11, 2013; Accepted: January 20, 2014; Published: February 24, 2014
Copyright: © 2014 Porta et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Funding: This study was supported by the São Paulo Foundation for Research Support/Brazil (FAPESP, process number 2010/52070-4), and CAPES (AUXPE-CSF-PVES-2619/2013 process number 23028.007721/2013-41). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The spontaneous fluctuations of heart period (HP) about its mean value observable in five minutes’ recordings are the apparent manifestation of the short-term cardiovascular control [1], [2]. Short-term cardiovascular regulation is carried out by a set of interacting neural and non neural components simultaneously operating over a range of frequencies from 0.04 to 0.5 Hz in humans [3]. Since these regulatory mechanisms work according to similar but not coincident temporal scales and they are coordinated by the autonomic nervous system but maintain a certain degree of autonomy to accomplish specific local tasks (e.g. the maintenance of the peripheral vasomotion at the district level in presence of vasoconstriction), the dynamics of HP changes cannot be fully described by a finite number of strictly periodic, fully predictable, oscillations. Complexity analysis quantifies the departure of a given signal from a fully predictable course [4]–[11]: the smaller the predictability, the higher the complexity. The improvement of predictability of an assigned effect signal when a presumed cause is introduced in the multivariate data set has been suggested to be a measure of the strength of the causal relation from the cause to the effect [12]: the larger the predictability improvement, the strong the intensity of the cause-effect link.
It is well known that aging influences the complexity of the cardiovascular control, as assessed from the analysis of HP variability, by reducing the number of temporal scales involved into the regulatory process, especially in the high frequency band (i.e. above 0.15 Hz) [4]–[9]. This information is clinically relevant because it was suggested that complexity analysis of HP variability can provide noninvasive indexes for monitoring the aging process and the susceptibility of individuals to injury and illness [10]. Nonetheless, two main issues deserve elucidation. The first issue is related to the traditional approach to assess the influence of age on the cardiovascular control: it is almost exclusively based on the analysis of HP variability. However, recently it has been pointed out that complexity analysis of systolic arterial pressure (SAP) variability can provide additional information [11]. We hypothesize that tracking the course of complexity of SAP variability and respiration (RESP) with age can provide information closely related to senescence of vascular and respiratory systems, thus complementing the traditional view exploring solely the senescence of cardiac control according to the analysis of HP variability. The second issue is linked to the possibility provided by causality tools [12]–[16] in interpreting changes of complexity of a designated variable in terms of modifications of the strength of the relation between the variable and its determinants [17]. For example, it is well-known that SAP variability and RESP contribute to HP oscillations respectively through the cardiac baroreflex [18]–[20] and the coupling between respiratory activity and vagal outflow [18], [21]–[23]. In a more complete universe of knowledge including SAP and RESP variability in addition to the HP one, causality analysis might explain the variation of the complexity of HP variability observed, for example, during STAND [24] as the result of the variation of the strength of the causal relation from SAP to HP and/or from RESP to HP. We hypothesize that complementing traditional univariate complexity analysis with the assessment of the strength of the causal relations via causality analysis might provide a more insightful description of the evolution of the cardiovascular control with age by linking the observed change of complexity to physiological mechanisms and their impairment.
According to these hypotheses the aim of this study is to assess the influence of aging on the complexity of HP, SAP and RESP variabilities and on the strength of the causal relations among them. The study exploits a multivariate compact framework devised to favor the interpretation of changes of complexity in terms of modifications of the strength of the physiological links. The main features of this approach are the construction of a multivariate embedding space accounting for the interactions among HP, SAP and RESP variabilities, the optimization of the number of samples necessary to describe cardiovascular variability interactions (i.e. the embedding dimension), the estimation of the interaction delays between different variability series and between samples of the same series, the evaluation of the performance of two alternative classes of methods for the assessment of complexity and causality, i.e. the linear model-based (MB) and nonlinear model-free (MF) class, and the comparison between the two techniques most frequently utilized in the MF class, i.e. local predictability (LP) and conditional entropy (CE). The effect of senescence was assessed at baseline in supine resting condition (REST) and during sympathetic activation induced by changing posture from REST to active standing (STAND) in a cohort of 100 healthy subjects from 21 to 70 years subdivided into five bins of age each including 20 subjects.
Methods
General Definitions
Given M series y1 = {y1(n), n = 1,…,N}, …, yM = {yM(n), n = 1,…,N} where N is the series length and n is the progressive counter, the series are first normalized to have zero mean and unit variance. We define Ω = {y1,…,yi,…,yj,…,yM} as the universe of the knowledge about the system under study. After labeling as Yji(n) = [yj(n-τji),…,yj(n-pi)] with 1≤i,j≤M, the embedding vector formed by pji = pi-τji+1 components of yj, where τji and pi represent the minimal and the maximal delay of influence of yj on yi, the multivariate embedding vector accounting for all signals present in Ω is Zi(n) = [Y1i(n),…,Yii(n),…,Yji(n),…,YMi(n)] and has dimension(1)where 1≤τii≤pi and 0≤τmi≤pi with m≠i. The possibility that τmi = 0 is introduced to account for instantaneous (i.e., non-delayed) effects from ym to yi. All the multivariate embedding vectors form a set Zi = {Zi(n), n = pi+1,…,N}. These definitions can be easily extended to the universe Ω after the exclusion of yj, i.e. Ω\yj = {y1,…,yi,…,yM}. In this case the multivariate embedding vector obtained from Zi(n) after excluding Yji(n) is Zi(n)\Yji(n) = [Y1i(n),…,Yii(n),…,YMi(n)] and the set of all Zi(n)\Yji(n) is indicated as Zi\Yji.
Linear MB Approach for the Assessment of Complexity
This technique computes complexity as the degree of unpredictability of yi in Ω according to a MB approach [25]. The multivariate linear regression model was exploited to describe the interactions of yi with all the signals present in Ω [25]. More specifically,(2)
where Ai is the qix1 vector containing the coefficients of the model, Ai = [A1i,…,Aii,…,Aji,…,AMi] with Aji = [aj(τji),…,aj(pi)], pi is the order of the model, the symbol T is the transpose operator, and εi is the white noise with zero mean and variance λi2. Equation (2) describes yi(n) as a linear combination of yj(n-k) with τji≤k≤pi and 1≤j≤M weighted by the constant coefficient aj(k) plus εi(n) representing the part of yi(n) that cannot be predicted in Ω. We follow the criterion leading to the minimization of λi2 in Ω to identify Ai [25]. Given the estimate of Ai, , the prediction of yi(n) is
(3)
Defined the prediction error, ei(n), as the difference between yi(n) and its prediction, , a measure of the unpredictability of yi in Ω is the mean square prediction error indicated as
in the following.
is computed as a function of the model order, pi, according to a procedure leading to add one sample for each signal present in Ω at every increment of pi (i.e. M samples at a time). Since at every increment of pi the in-sample ability of the model in fitting the data improves,
progressively decreases with pi, thus making practically useless the monitoring of
with pi to decide the optimal model order [25]. Therefore, instead of tracking
it is a common practice to monitor a figure of merit, here the Akaike figure of merit [26]. The Akaike figure of merit adds a term gradually increasing with pi to a function depending on
, thus forcing the creation of a minimum. The model order at the minimum of the Akaike figure of merit will be indicated as pio.
at pio is taken as a normalized complexity index (NCI) of yi in Ω based on the MB approach, and labeled as
. Unless necessary to stress the universe of knowledge, NCIiMB will be adopted instead of
and, in this case, we will assume that the universe of knowledge is Ω. Since all the series are normalized to have unit variance, NCIiMB ranges between 0 and 1, where 0 indicates null complexity of yi and its full predictability and 1 indicates maximal complexity of yi and its full unpredictability.
Linear MB Approach for the Assessment of Causality
Causality from yj to yi based on MB approach is evaluated according to the causality ratio (CR)(4)
where and
represent NCIiMB assessed respectively in Ω and Ω\yj. CRj→iMB quantifies the strength of the causal relation from yj to yi as the unpredictability decrement of yi due to the inclusion of yj in Ω\yj.
was computed by estimating the coefficients of the polynomials in Ω\yj while maintaining pio optimized in Ω. CRj→iMB<0 indicates that yj carries unique information about the future evolution of yi that cannot be derived from any signal in Ω\yj, and, according to the concept of Granger causality [12], it can be stated that yj Granger-causes yi in Ω.
Nonlinear MF Approach for the Assessment of Complexity Based On LP
This technique computes complexity as the degree of unpredictability of yi in Ω [11]. The method hypothesizes that there is a function f(.) linking Zi(n) to yi(n). yi(n) is usually indicated as the image of Zi(n) through f(.) and, vice versa, Zi(n) is said to be the vector associated to yi(n) through f(.). The k-nearest-neighbor approach provides a local approximation of f(.) under smooth conditions and without making any specific assumption on the dynamical relationship linking Zi(n) to yi(n) [27], thus allowing the prediction of yi(n), , according to a MF approach. More specifically,
is computed in Ω as a combination over a subset of values of yi whose associated embedding vectors belong to the set of the k points closest to the embedding vector, Zi(n), associated to the value to be predicted, yi(n). Usually, the set of the closest points does not include vectors with time indexes close to n to reduce the correlation among the selected values of yi solely due to the temporal closeness of the associated embedding vectors [28]. To combine the k selected values of yi, the weighted mean, where the weights are the inverse of the distance from Zi(n) to each neighbor vector is a natural choice (i.e. zero-order predictor) [29]. The distance between vectors is calculated according to a predefined norm (here we use the maximum norm, i.e. the absolute value of the maximal difference between corresponding components) [30]. After calculating the prediction of yi,
, the square correlation coefficient between yi and
, ri2, can be assessed [31]. ri2 is bounded between 0 and 1, respectively indicating null and perfect predictability of yi. Predictability depends on the construction of the multivariate embedding space. Here we exploit a strategy growing up progressively the embedding space by adding one lagged component at a time selected according to a given criterion [14], [32]. Let us refer to as candidate the lagged component being possibly selected as a new coordinate of the embedding space. Initially (i.e. when qi = 0), the set of candidates tested to predict yi(n) is {y1(n-τ1i),…,y1(n-pi),…,yi(n-τii),…,yi(n-pi),…,yj(n-τji),…,yj(n-pi),…,yM(n-τMi),…,yM(n-pi)}. All candidates are tested by assessing ri2 and the new component is selected as the one maximizing ri2. Then, the selected component is retained, qi is increased by 1 and the set of candidates is reduced by excluding samples of a signal with time indexes more recent or equal to any component of the same signal already exploited to form the multivariate embedding space, thus avoiding duplicate selections and speeding up reconstruction. The process continues until the set of candidates is empty. ri2 varies with qi and, provided that Ω is helpful to predict yi, its course is the result of two opposite tendencies: i) a larger number of components increases the ability of Zi to predict yi by allowing a more accurate unfolding of the trajectories in the multivariate embedding space and reducing the ambiguities in mapping Zi to future samples, thus leading to an increase of ri2 with qi; ii) in presence of noise and time series of limited length, vectors in Zi progressively take apart with qi in the multivariate embedding space and, consequently, the ability of the k closest points to predict yi worsens, thus leading to a decrease of ri2 with qi. Therefore, ri2 exhibits a maximum over qi,or a minimum if the functional 1-ri2 is considered [31]. Conversely, if yi cannot be predicted in Ω, ri2 remains close to 0 regardless of qi. The largest value of ri2 over qi, max(ri2), is complemented to 1 (i.e. 1-max(ri2)) to obtain the NCI of yi based on LP in Ω labeled as
. Unless necessary to stress the universe of knowledge, from now on NCIiLP will be adopted instead of
and, in this case, we will assume that the universe of knowledge is Ω. NCIiLP ranges from 0, indicating null complexity and perfect predictability of yi, to 1, representing maximal complexity and full unpredictability of yi. qi at NCIiLP, qiLPo, represents the optimal number of components of Zi(n) leading to the best prediction of yi(n). qiLPo is bounded between 0 and qi given in (1).
Nonlinear MF Approach for the Assessment of Causality Based on LP
Causality from yj to yi based on LP is evaluated according to the CR(5)
where and
represent NCIiLP assessed respectively in Ω and Ω\yj. CRj→iLP quantifies the strength of the causal relation from yj to yi as the decrement of uncorrelation between yi and
due to the inclusion of yj in Ω\yj. In the assessment of CRj→iLP,
was computed by excluding the components of yj from the optimal embedding vector leading to the calculation of
. With this expedient, if the optimal embedding vector does not contain any component of yj, CRj→iLP = 0, indicating the absence of causality from yj from yi. Conversely, CRj→iLP<0 indicates that yj carries unique information about the future evolution of yi that cannot be derived from any signal in Ω\yj, and, according to the concept of Granger causality [12], it can be stated that yj Granger-causes yi in Ω.
Nonlinear MF Approach for the Assessment of Complexity Based on CE
The technique computes complexity as the amount of information carried by yi that cannot be derived from any of the series present in Ω through the calculation of the CE [33]. A k-nearest-neighbor approach can be exploited to compute the CE of yi in Ω [34]. CE is computed as the average value of the Shannon entropy of the conditional distribution of yi(n) given Zi(n). The distribution of yi(n) given Zi(n) is built by considering the subset of values of yi whose associated embedding vectors belong to the set of the k points closest to Zi(n). The Shannon entropy of the distribution of yi(n) given Zi(n) is calculated as the negative natural logarithm of the average probability that two samples of the distribution of yi(n) given Zi(n) are at distance closer than ε [34]. In this study the tolerance ε was set to 10% of the difference between the 84th and the 16th percentile of yi and the distance was computed, as in the LP approach, using the maximum norm. CE is bounded between 0 and the Shannon entropy of yi indicating respectively the minimal and the maximal amount of information contained in yi. CE is smaller than the Shannon entropy of yi when the set of conditioning vectors, Zi, is helpful to reduce the uncertainty associated to yi. CE depends on the construction of the multivariate embedding space. Here we exploit a strategy to construct and optimize the embedding space analogous to the one exploited by the LP approach. Only the optimization criterion is different. Indeed, the new component is selected among the set of candidates as the one producing the maximal decrement of the uncertainty of yi (i.e. the minimum of CE) compared to the maximal level of uncertainty of yi quantified by the Shannon entropy of yi [14]. CE varies with qi and, provided that Ω is helpful to reduce the uncertainty of yi, the course of CE is the result of two opposite tendencies: i) a larger number of components increases the ability of Zi to reduce the uncertainty of yi because longer conditioning patterns have better possibilities in fixing future samples, thus leading to a decrease of CE with qi; ii) in presence of noise and time series of limited length, vectors in Zi progressively take apart with qi in the multivariate embedding space and, consequently, the ability of the k closest points to limit uncertainty vanishes, thus leading to the increase of CE with qi. Therefore, CE exhibits a minimum over qi. Conversely, if the information carried by yi cannot be reduced in Ω, CE remains close to the Shannon entropy of yi regardless of qi. The minimum of the CE over qi is normalized by the Shannon entropy of yi to obtain a NCI of yi based on CE in Ω, labeled as . Unless necessary to stress the universe of knowledge, from now on NCIiCE will be adopted instead of
and, in this case, we assume that the universe of knowledge is Ω. NCIiCE ranges from 0, indicating null information and complexity of yi, to 1, representing maximal information and complexity of yi. qi at NCIiCE, qiCEo, represents the optimal number of components of Zi(n) leading to the maximal reduction of uncertainty of yi. qiCEo is bounded between 0 and qi given in (1).
Nonlinear MF Approach for the Assessment of Causality Based on CE
Causality from yj to yi based on CE is computed according to the CR(6)
where and
represent NCIiCE assessed respectively in Ω and Ω\yj. CRj→iCE quantifies the strength of the causal link from yj to yi as the decrement of information carried by yi attributable solely to the inclusion of yj in Ω\yj. We use the same expedient exploited in the calculation of CRj→iLP resulting in CRj→iCE = 0 when no components of yj are present in the optimal embedding vector leading to the computation of
. Conversely, CRj→iCE<0 indicates that yj is capable to reduce the uncertainty of yi to a level that cannot be achieved by exploiting any signal in Ω\yj and, according to the concept of transfer entropy [35], it can be stated that yj causes yi in Ω.
Experimental Protocol and Data Analysis
Ethics Statement
The study was performed according to the Declaration of Helsinki and it was approved by the Human Research Ethics Committee of the Federal University of São Carlos (protocol number 173/2011). A written informed consent was obtained from all subjects.
Experimental Protocol
We studied 100 nonsmoking healthy humans (54 males, age from 21 to 70 years, median = 45 years; weight from 43 to 100 Kg, median = 71 Kg; height from 146 to 197 cm, median = 167 cm; body mass index (BMI) from 17.4 to 33.4 Kg.m−2, median = 25 Kg.m−2). The population was composed by 20 subjects in each of the following bins: from 21 to 30 years (10 males, median age = 26 years, median weight = 71 Kg, median height = 168 cm, median BMI = 23.9 Kg.m−2); from 31 to 40 years (11 males, median age = 34 years, median weight = 69 Kg, median height = 168 cm, median BMI = 24.8 Kg.m−2); from 41 to 50 years (10 males, median age = 45 years, median weight = 70 Kg, median height = 167 cm, median BMI = 25.4 Kg.m−2); from 51 to 60 years (10 males, median age = 55 years, median weight = 71 Kg, median height = 169 cm, median BMI = 25.1 Kg.m−2); from 61 to 70 years (13 males, median age = 65 years, median weight = 72 Kg, median height = 164 cm, median BMI = 26.7 Kg.m−2). The population was balanced in terms of gender to limit the influences of this confounding factor on the analysis [36]. The peak oxygen uptake, peak , was computed during a maximum cardiopulmonary exercise test performed on a treadmill using a ramp protocol. The criteria for exercise interruption were based on [37] and the assessment was made according to [38]. In our population the median (25th–75th percentile) of the peak
in each age bin was 2271 (1574–2894), 2340 (1993–3135), 2384 (1712–3080), 2026 (1363–2484) and 1825 (1436–2025) mL.min−1. All the subjects were apparently healthy, had no history and no clinical evidence of any disease based on clinical and physical examinations, laboratory tests, standard electrocardiogram (ECG) and on a maximum cardiopulmonary exercise test conducted by a physician. They were not taking any medication known to interfere with cardiovascular control. Smokers and habitual drinkers were excluded from this study. All subjects were evaluated in the afternoon. The experiments were carried out in a climatically controlled room (22–23°C) with relative air humidity at 40–60%. Subjects were instructed not to consume caffeinated and alcoholic beverages as well as not to perform strenuous exercises on the day before the recording. They were also instructed to ingest a light meal at least 2 hours prior to the test. On the experimental day, the subjects were interviewed and examined before the test to verify whether they were in good health and they had a regular night sleep. Prior to the recording, the volunteers were made familiar with the equipment and with the experimental procedure. During the entire protocol the subjects breathed spontaneously and they were not allowed to talk.
ECG (modified lead I), continuous plethysmographic arterial pressure (Finometer PRO, Finapress Medical System, The Netherlands) and respiratory movements via thoracic belt (Marazza, Monza, Italy) were digitalized using a commercial device (BioAmp Power Lab, AD Instruments, Australia). Signals were sampled at 400 Hz. The arterial pressure was measured from the middle finger of the left hand being maintained at the level of heart by fixing the subject’s arm to his thorax. All the experimental sessions of the protocol included two periods in the same order: 1) 15 minutes at REST; 2) 15 minutes during STAND. Before REST we allowed 10 minutes for stabilization. The arterial pressure signal was cross-calibrated in each session using a measure provided by a sphygmomanometer at the onset of REST. The auto-calibration procedure of the arterial pressure device was switched off after the first automatic calibration at the onset of the session. Analyses were performed after about 2 minutes from the start of each period.
Extraction of the Beat-to-beat Variability Series
After detecting the R-wave on the ECG and locating its peak using parabolic interpolation, HP was approximated as the temporal distance between two consecutive parabolic apexes. The maximum of arterial pressure inside of the n-th HP, HP(n), was taken as the n-th SAP, SAP(n). The signal of the thoracic movements was down-sampled once per cardiac beat at the occurrence of the first R-wave peak delimiting HP(n), thus obtaining the n-th RESP measure, RESP(n). HP(n), SAP(n) and RESP(n) were expressed in ms, mmHg and arbitrary units (a.u.) respectively. The occurrences of R-wave and SAP peaks were carefully checked to avoid erroneous detections or missed beats. After extracting the series HP = {HP(n), n = 1,…,N}, SAP = {SAP(n), n = 1,…,N} and RESP = {RESP(n), n = 1,…,N}, where n is the progressive cardiac beat counter and N is the total cardiac beat number, sequences of 256 consecutive measures were randomly selected inside REST and STAND periods, thus focusing short-term cardiovascular regulatory mechanisms [3]. If evident nonstationarities, such as very slow drifting of the mean or sudden changes of the variance, were present in spite of the linear detrending, the random selection was carried out again.
Traditional Index Calculation
Traditional time domain parameters such as mean and variance of HP and SAP were calculated and indicated as µHP, µSAP, σ2HP and σ2SAP. They were expressed in ms, mmHg, ms2 and mmHg2 respectively. We assessed cardiac baroreflex sensitivity (BRS) according the baroreflex sequence method [39]. The technique relies on the search for sequences characterized by the contemporaneous increase (positive sequence) or decrease (negative sequence) of HP and SAP. Both positive and negative sequences are referred to as baroreflex sequences. They are identified according to the following prerequisites: 1) the length of the sequences was four beats (three increases or decreases); 2) the lag between HP and SAP values was set to 0; 3) the total HP variation was larger than 5 ms; 4) the total SAP variation was larger than 1 mmHg; v) the correlation coefficient in the plane [SAP(n),HP(n)] was larger than 0.85. When a baroreflex sequence matched with the above mentioned prerequisites the slope of the regression line in the plane [SAP(n),HP(n)] was calculated and averaged over all baroreflex sequences. This average was indicated hereafter as BRS and expressed as ms.mmHg−1. The total number of baroreflex sequences was also retained and divided by the number of meaningful SAP ramps where a SAP ramp was identified as a sequence of three consecutive SAP increases or decreases leading to a total SAP variation larger than 1 mmHg and a correlation coefficient in the plane [n,SAP(n)] larger than 0.85. This index was referred to as cardiac baroreflex effectiveness index (BEI) [40] and expressed as dimensionless units. This index can be considered an index of causality from SAP to HP because it gives the proportion of SAP ramps capable of evoking a meaningful cardiac baroreflex response.
Complexity and Causality Indexes Calculation
NCI was calculated from HP, SAP and RESP series (i.e. NCIHP, NCISAP and NCIRESP) according to MB, LP and CE approaches. CR was computed from any pair of series (i.e. from SAP and RESP to HP, CRSAP→HP and CRRESP→HP, from HP and RESP to SAP, CRHP→SAP and CRRESP→SAP, from HP and SAP to RESP, CRHP→RESP and CRSAP→RESP) according to MB, LP and CE approaches. In the MB approach the optimal model order was chosen as the one minimizing the Akaike figure of merit in the range from 1 to 8. In the case of LP and CE techniques the set of initial candidates was formed by 8 components for each series. The delays τSAPHP and τRESPHP were set to 0 to allow the description of the fast vagal reflex (within the same cardiac beat) capable to modify HP in response to changes of SAP and RESP [41], [42]. The delay τRESPSAP was set to 0 to account for the potential rapid effect of RESP on SAP due to the immediate transfer of an alteration of intrathoracic pressure on SAP value [1]. The delay τHPSAP was set to 1 due to the measurement conventions preventing that the modification of HP(n) can affect SAP(n) [18]. According to [43] actions of HP and SAP on RESP were slower (i.e. they cannot occur in the same beat), thus leading to τHPRESP = 1 and τSAPRESP = 1. The number of nearest neighbors, k, was set to 30 for both LP and CE approaches [31], [34].
Statistical Analysis
Two way repeated measures analysis of variance (one factor repetition, Holm-Sidak test for multiple comparisons) was utilized to test the significance of the differences between the optimal model order assessed according to the MB approach and the number of components utilized to build the optimal multivariate embedding space according to LP and CE methods within the experimental condition (i.e. REST or STAND) while varying the method (i.e. MB, LP and CE) and within the method while varying the experimental condition. The null hypothesis of Normal distribution of all series and of all parameters extracted from them was tested according to Kolmogorov-Smirnov test. Linear regression analysis of µHP, σ2HP, µSAP, σ2SAP, BRS and BEI on age was carried out. If Normality test over µHP, σ2HP, µSAP, σ2SAP, BRS and BEI passed Pearson product moment correlation coefficient was calculated. Otherwise, Spearman rank order correlation coefficient was computed. The same procedure was carried out to check the dependence of NCI and CR on age. Statistical analysis was carried out using a commercial statistical program (Sigmastat, SPSS, ver.3.0.1). A p<0.05 was always considered as significant.
Methodological Results
Comparison of the MB, LP and CE Approaches
Figure 1 shows the grouped bar-graphs of the optimal model order, po, according to the minimum of the Akaike figure of merit derived from the MB approach (white bar) and the number of components utilized to build the optimal embedding space, qo, according to the smallest unpredictability obtained from the LP technique (gray bar), and to minimal amount of information derived from the CE method (black bar). Values of po and qo are plotted as mean plus standard deviation over the entire cohort of subjects regardless of age and as a function of the experimental condition (i.e. REST and STAND). po and qo were assessed while varying the designated effect series, i.e. HP, SAP and RESP in Figs. 1a,b,c respectively. Independently of the assigned effect series methods provided significantly different values of po and qo within the same experimental condition (i.e. REST or STAND). These differences indicated that the optimal description of the interactions among variability series was achieved using the largest and the smallest number of components by the MB and CE approaches respectively, while the LP technique was in between them. It is worth noting that no difference between experimental conditions within the same method was observable except in the case of qHPCEo (Fig. 1a): qHPCEo during STAND was significantly smaller than qHPCEo at REST.
Grouped bar-graphs show the mean (plus standard deviation) of the optimal model order, po, leading to the minimum of the Akaike figure of merit according to the MB approach (white bar) and of the number of components, qo, utilized to build the optimal multivariate embedding space leading to smallest unpredictability according to the LP method (gray bar) and to the minimal amount of information according to the CE technique (black bar). po and qo are shown as a function of the experimental condition (i.e. REST and STAND) and values are pooled together regardless of age. The designated effect series are HP, SAP and RESP in (a), (b) and (c) respectively. The symbols # indicates a significant difference with p<0.05 within the same experimental condition (i.e. REST or STAND) while varying the method. The symbols * indicates a significant difference with p<0.05 within the method while varying the experimental condition.
Experimental Results
Linear Regression Analysis of Traditional Parameters on Age
Correlation analysis was carried out to assess the association of traditional parameters derived from HP and SAP variability on age at REST. While µHP and BEI was unrelated to age, σ2HP, µSAP, σ2SAP and BRS were found significantly correlated with age. The correlation coefficient was negative in the case of σ2HP and BRS, thus indicating that σ2HP and BRS progressively decreased with age, and it was positive in the case of µSAP and σ2SAP, thus evidencing that µSAP and σ2SAP increased with age. Remarkably, the probability of the type I error assessed over σ2HP, σ2SAP and BRS were at least three orders of magnitude smaller than that relevant to µSAP.
Correlation analysis was carried out to assess the association of traditional parameters obtained from HP and SAP variability on age during STAND as well. At difference with REST µHP and BEI was significantly linearly correlated with age, while no linear relation was detected between σ2SAP and age. The sign of the correlation coefficient of µHP and BEI suggests that the orthostatic challenge induced a tachycardic response and an association from SAP to HP becoming less and less important with age. Similarly to REST the progressive decrease of σ2HP and BRS and the gradual increase of µSAP was significant.
Linear Regression Analysis of Complexity Indexes on Age
Table 1 reports the results of the linear regression analysis of complexity parameters on age at REST. Independently of the approach (i.e. MB, LP or CE) NCIHP was significantly linearly correlated with age. The correlation coefficients of NCIHPMB, NCIHPLP and NCIHPCE were negative, thus suggesting a progressive loss of complexity of HP dynamics with age. The analysis of SAP complexity points out differences among the methods. Indeed, while NCISAP derived from MF techniques, NCISAPLP and NCISAPCE, were significantly associated with age, NCISAPMB was unrelated to it. The correlation coefficients of NCISAPLP and NCISAPCE were negative, thus indicating a progressive loss of complexity of the SAP series with age. Regardless of the method complexity of RESP dynamics was unaffected by aging.
Table 2 reports the results of the linear regression analysis of complexity indexes on age during STAND. In this condition solely the complexity of HP series was influenced by senescence. Indeed, NCISAP and NCIRESP were unrelated to age. NCIHP was found significantly related to age when computed according to MB and LP approaches, while NCIHPCE was independent of it. The correlation coefficients of NCIHPMB and NCIHPLP on age were positive, thus suggesting that HP variability during STAND became more and more unpredictable with age.
Linear Regression Analysis of Causality Indexes on Age
Table 3 reports the results of the linear regression analysis of causality parameters on age at REST. Regardless of the method (i.e. MB, LP and CE) CRSAP→HP, CRRESP→SAP, CRHP→RESP and CRSAP→RESP showed the same association with age: while CRSAP→HP, CRHP→RESP and CRSAP→RESP were unaffected by age, CRRESP→SAP was linearly associated to it. The correlation coefficient of CRRESP→SAP was positive, thus suggesting that aging reduced the strength of the causal link from RESP to SAP. The methods were not in full agreement in the case of CRHP→SAP and CRRESP→HP. Indeed, the MB approach was the unique technique detecting a significant linear relation of CRHP→SAP on age. The correlation coefficient of CRHP→SAPMB was positive, thus indicating that the strength of the causal link from HP to SAP diminished with age. In the case of CRRESP→HP, both LP and CE approaches found a significant linear association with age, while the MB one was unable to detect it. The correlation coefficients of CRRESP→HPLP and CRRESP→HPCE were positive, thus suggesting that the influence of RESP on HP became more and more ineffective with age.
Table 4 reports the results of the linear regression analysis of causality indexes on age during STAND. The association of CRRESP→HP, CRHP→RESP and CRSAP→RESP with age did not depend on the method (i.e. MB, LB and CE). Indeed, all the three approaches found that CRRESP→HP, CRHP→RESP and CRSAP→RESP were unrelated to age. Conversely, the association of CRSAP→HP, CRHP→SAP and CRRESP→SAP with age depended on the method. Only the LP approach was able to detect the significant linear relation of CRSAP→HP with age. The correlation coefficient of CRSAP→HPLP was positive, thus indicating that the strength of the causal link from SAP to HP became weaker and weaker with age. The significant linear association of CRHP→SAP with age was detected only by the MB technique. The correlation coefficient of CRHP→SAPMB was positive, thus indicating that a progressive decrease of the intensity of the causal relation from HP to SAP with age. In the case of CRRESP→SAP, only the CE method was unable to detect a significant linear relation of the strength of the causal link from RESP to SAP with age. Indeed, the correlation coefficients of CRRESP→SAPMB and CRRESP→SAPLP were significant and positive, thus suggesting that the senescence reduced the association between RESP and SAP variability in the temporal direction from RESP to SAP.
Discussion on the Methodological Findings
The methodological findings of this study can be summarized as follows: i) the study proposes a compact framework for the assessment of complexity of cardiovascular variability and causality of their interactions in a multivariate embedding space; ii) the common feature of all techniques designed within the proposed framework is the optimization of the embedding dimension; iii) this framework facilitates the comparison between linear MB and nonlinear MF classes and, inside the MF class, between LF and CE techniques; iv) the simultaneous assessment of complexity and causality favors the interpretation of changes of complexity in terms of modifications of the strength of physiological relations.
A Multivariate Compact Framework for the Assessment of Complexity and Causality in Cardiovascular Variability Series
We applied a multivariate compact framework for the assessment of dynamical complexity and causality over an arbitrary set of signals, Ω. In the case of MB and LP approaches the complexity of an assigned effect series was evaluated as its degree of unpredictability in Ω (i.e. the portion of normalized variance of the assigned effect series that cannot be explained based on all signals present in Ω). In the case of the CE method the complexity of the selected effect series was assessed as the amount of information carried by the designated effect series that cannot be derived from the signals present in Ω. In the case of MB and LP approaches the strength of the causal relation from a cause series to an effect one in Ω was computed, according to the notion of Granger causality [12], as the fractional decrement of unpredictability of the assigned effect series resulting from the inclusion of the cause series into the incomplete Ω disregarding the supposed cause. In the case of CE method causality was assessed, according to the notion of information transfer [35], as the fractional decrement of the information carried by the selected effect series resulting from the inclusion of the cause series into the incomplete set of conditioning signals formed by Ω devoid of the presumed cause.
One of the main features of the framework is the assessment of complexity of a designated effect series into a multivariate embedding space built in Ω. This characteristic distinguishes the proposed approach from univariate methods exploiting embedding spaces built over a unique series [29], [30], [44]–[47]. In our application the multivariate embedding space might provide a more efficient description of the overall complexity of the cardiovascular control than the univariate embedding space usually constructed using solely the HP components because it allows a more complete representation of the dynamical behavior of the subsystems contributing to the overall functioning of the cardiovascular system. Indeed, the joint exploitation of variables directly linked to cardiac, vascular and respiratory subsystems (i.e. HP, SAP and RESP) might reveal portions of the cardiovascular system that otherwise might remain unveiled or undervalued using the time course of a single variable. Since the multivariate embedding space is the starting point for causality analysis [12], [48], the proposed framework allows the joint evaluation of complexity and causality indexes. At difference with most of the applications of causality approaches, the proposed framework is grounded on the optimization of the multivariate embedding dimension, thus avoiding the common practice of its arbitrary assignment [13], [49], [50] or its opportunistic setting to a value allowing the best separation among populations and/or experimental conditions [16]. In our study the multivariate embedding dimension was optimized by MB, LP and CE approaches according to, respectively, the classical Akaike information criterion, the maximization of the predictability of the assigned effect series, and the minimization of the information carried by the assigned effect series. The optimization of the multivariate embedding dimension allows the reduction of the number of parameters needed to be set, thus favoring the automatic computation of complexity and causality markers.
Another element of novelty of the proposed framework is the possibility of performing the direct comparison between MB and MF classes of methods. The MB methods interpret the dynamical interactions according to a model, while the MF ones are completely data-driven and do not impose any particular predefined structure to the dynamics. While several studies advocated the use of MF approaches for the computation of causality due to their straight adherence to the data and their ability to account for nonlinearities of any type [16], [48], the presumed superiority of MF tools was never checked over real, short and noisy data via a direct comparison to MB methods. Understanding the additional information that might be derived from MF methods compared to MB ones is a relevant issue because MB approaches, especially if based on multivariate linear regression models, are well-established, robust and computationally efficient [25]. The proposed framework compares MF approaches with a traditional MB method based on multivariate linear regression model [25]. In this MB approach the multivariate embedding space was built sequentially by adding M delayed components at a time, one for any of the signals present in Ω. This strategy, quite common in applications of multivariate linear regression models [13], [51]–[55], might lead to an overparametrization of the model and to possible suboptimality of the description of the interactions among signals because the number of components exploited in any regression is constrained to be the same. Conversely, the strategy followed by both the proposed MF approaches (i.e. LP and CE) leads to a reconstruction of the multivariate embedding space using fewer components. Defined as parsimoniousness of a model its ability to describe interactions among signals using a small number of components, it can be stated that the proposed MF approaches are more parsimonious than the MB one. Indeed, at every increment of the embedding dimension we tested all the candidate components and we selected the one providing the maximal correlation between the original effect series and its prediction in the case of the LP approach, and the maximal reduction of uncertainty about the future evolution of the effect series in the case of the CE method. Parsimoniousness is a quality because limits overfitting even though, when it is excessive, it might become a disadvantage because it leads to underfitting. The optimization of the embedding dimension allows the best unfolding of the system dynamics as well as it solves the issue of setting the delay of the interactions from any cause signal to the effect one and the lag between components belonging to the same signal. The setting of these parameters is a non trivial issue. Indeed, the interaction delay is usually set according to physiological considerations about the latency of reflexes [56], while the delay among components of the same signal is set such a way to minimize correlation among samples [57]. However, although the values of these parameters can be found in literature or can be easily estimated by computing the decorrelation time, their reliability is questionable. Indeed, while the latency of the interactions is usually estimated by opening closed loop circuits using surgical procedures and pharmacological challenges, thus being inadequate when closed loop reflexes are active, decorrelation time is an average parameter that cannot be optimal for any time scale, thus being inappropriate in presence of multiple temporal scales.
In addition, the proposed framework allows the comparison between two different techniques for the evaluation of complexity and causality inside the same class of MF approaches (i.e. LP and CE). The differences between LP and CE methods depend on the exploitation of two different functionals. While the LP approach exploited the correlation between observed and predicted data, the CE method quantified directly the uncertainty in the information domain. The comparison was carried out without altering the logic of coarse graining of the system evolution in the multivariate embedding space. Indeed, both LP and CE methods exploited the k-nearest-neighbor technique [27], [34]. The k vectors closest to the current one, regardless of their actual distance from it, formed the library of patterns utilized to build the conditional distribution of the future values (i.e. the distribution of the images of the k-nearest-neighbor vectors). The conditional distributions are utilized to set the predictor and assess the uncertainty as the variance of the prediction error in the case of the LP approach or as the average value of the Shannon entropy of the conditional distributions in the case of the CE method. Therefore, this framework provides a unique opportunity to compare different MF approaches under similar settings.
The framework proposes the contemporaneous evaluation of complexity and causality indexes. This simultaneous evaluation is important because the assessment of causality might be helpful to interpret variations of complexity with age in terms of physiological mechanisms. For example, since at REST the decrease of the strength of the causal link from RESP to HP combined with the invariable intensity of the causal relation from SAP to HP is incompatible with the decreased complexity of HP variability with age, we suggest that mechanisms other than cardiorespiratory coupling and cardiac baroreflex are responsible for the loss of complexity of the HP variability during senescence. During STAND the increase of the complexity of the HP variability with age was associated to a progressive decrease of the strength of the causal link from SAP to HP in presence of an invariable intensity of the causal relation from RESP to HP. Therefore, during STAND the observed increase of the complexity of the HP variability with age can be attributed to the progressive reduction of the effectiveness of the cardiac baroreflex. This finding suggests that cardiac baroreflex plays an important role in keeping low the complexity of HP variability and its impairment leads to an increase of the complexity of the cardiac control.
Comparison Between MB and MF Approaches
Our findings pointed out differences between MB and MF approaches. These differences might be explained either in terms of the ability of the MF approach to account for potential nonlinearities present at the level of the dynamics of the series and/or at the level of the relations among different series, or in terms of different capabilities of the approaches in limiting the curse of dimensionality (i.e. the gradual decrease of the reliability of complexity and causality indexes with the embedding dimension [58] unavoidable in presence of a limited amount of data).
At REST, at variance of the MB technique, the MF approaches were able to detect the decrease of complexity of the SAP variability with age. The decline of the SAP complexity with age might be associated to the well-known gradual increase of sympathetic activity during senescence [59], [60] leading to the gradual increase of synchronization of sparse vasomotor activities at peripheral vascular levels [61]. Given the inherent nonlinear nature of synchronization it is not surprising to find out that the age-related modifications of SAP dynamics cannot be detected by a linear approach such as the MB one.
Moreover, at REST, again at variance with the MB approach, the MF techniques were able to detect the dependence of the intensity of the causal relation from RESP to HP on age. The presence of a causal relation from RESP to HP is the result of a nonlinear coupling between respiratory activity and vagal outflow determining the respiratory-related HP variations [22], [23], [62]. Since nonlinear mechanisms provides the basis of the coupling between respiratory activity and vagal outflow [22]. it is not surprising to find out that a linear MB approach cannot detect the gradual reduction of the strength of the causal link from RESP to HP with age.
The MB approach was also unable to identify the progressive decrease of the strength of the causal link from SAP to HP during STAND detected by the LP approach. The causal relation from SAP to HP during STAND is due to the activation of cardiac baroreflex in response of the orthostatic challenge [63]. Since it was observed that the likelihood of finding nonlinear interactions between SAP and HP variabilities increased in old subjects [64], again it is not surprising that a nonlinear MF approach such as the LP one is more powerful than a linear MB method in detecting the decreased efficiency of the cardiac baroreflex with age.
Regardless of the experimental condition MF approaches were unable to identify the progressive decrease of the intensity of the causal relation from HP to SAP. The causal link from HP to SAP is the consequence of two opposite tendencies [18]: 1) a negative relation of SAP on HP due to the diastolic runoff leading to a decrease of diastolic arterial pressure while lengthening HP and, thus, to a reduced SAP at the next cardiac beat in presence of an unaltered pulse pressure; 2) a positive relation of SAP on HP due to the Frank-Starling mechanism leading to an increase of SAP at the next cardiac beat due to the larger ventricular filling induced by the HP lengthening. The ability of the MB approach in detecting the effect of age on the causal link from HP to SAP can be explained by hypothesizing that this relation needs high dimensional embedding spaces to be perfectly unfolded. Unfortunately, MF approaches exploited in this study might be inefficient in presence of high dimensional dynamics, thus giving the reason for their incapacity in detecting the relation revealed by the MB approach. Indeed, MF approaches escape to the curse of dimensionality by drastically limiting the embedding dimension at the cost of reducing their efficiency in presence of high dimensional dynamics and imposing the enlargement of the data set to reliably explore higher dimensional spaces. Conversely, as a consequence of the less parsimonious strategy to construct the multivariate embedding space, the MB approach can be more efficient in presence of high dimensional dynamics.
Comparison between LP and CE Approaches
Both MF approaches are parsimonious compared to the MB one, but the CE technique makes use of a smaller number of components than the LP one. Parsimoniousness is certainly helpful to limit the rate of false detection of causality when causality is not present (i.e. false positives), but it might increase the rate of false negatives (i.e. the probability of missing causality when a causal link does exist). This observation can explain the more frequent detection of age-related changes of complexity and causality indexes observed when the LP approach was applied compared to the CE one. It is worth noting that when both MF techniques detected a significant association of the parameters with age the sign of the relation is the same. The greater parsimoniousness of the CE approach might explain its inability to detect the progressive increase of complexity of the HP dynamics with age observed during STAND and to find out the progressive decrease of the strength of the causal link from SAP to HP and RESP to SAP with age during STAND. Both relations were identified by the LP method. However, since LP and CE techniques exploit different functionals it might happen that the differences between the two MF approaches might be related to the quantification of different aspects of dynamics as well.
Discussion on the Experimental Findings
The experimental results about the age dependency of the complexity of the cardiovascular control can be summarized as follows: a) while HP complexity progressively decreases with age at REST, it gradually increases with age during STAND, thus suggesting an age-related impairment of the cardiac control and of the response of the cardiac regulation to stressors; b) SAP complexity gradually decreases with age at REST, thus indicating the progressive increase of sympathetic activity and/or modulation with age; c) no association between SAP complexity and age is observed during STAND, thus being indicative of a reduced responsiveness of the vasomotor control to stressors; d) RESP complexity is unrelated to age, thus demonstrating that respiratory control is preserved in elderly people.
The experimental results about the age dependency of the strength of the causal relations among cardiovascular variabilities can be summarized as follows: 1) while at REST aging preserves the strength of the causal link from SAP to HP along the cardiac baroreflex, during STAND a negative relation with age was found, thus indicating that in elderly individuals the cardiac baroreflex remains a fundamental reflex for homeostasis at REST but its efficiency in response to stressors is progressively lost with age; 2) at REST and during STAND the strength of the causal link from HP to SAP decreases with age due to the progressive exploitation of Frank-Starling mechanism at REST and the progressive increase of peripheral resistances during STAND; 3) regardless of the experimental condition the intensity of the causal relation from RESP to HP decreases with age, thus reflecting the increase of sympathetic tone and the uncoupling between respiratory activity and vagal outflow; 4) regardless of the experimental condition the intensity of the causal relation from RESP to SAP decreases with age as a probable result of the progressive increase of left ventricular thickness and vascular stiffness and of the gradual decrease of the respiratory-related HP variations; 5) the causal link from SAP to RESP is not affected by aging, as a likely result of the negligible role played by this pathway in short-term cardiovascular control; 6) the causal relation from HP to RESP is unmodified by aging due to the preserved nerve conduction velocity and neural processing time with age.
Effect of Age on HP and SAP Traditional Parameters
We confirm that at REST the mean HP is unrelated to age [65], mean SAP progressively increases with age [65], and HP variance gradually decreases [4], [8], [66]. The tendency of SAP variance to increase with age at REST observed in [67] was found to be significant in this study. Several mechanisms have been advocated to explain these relations with age: i) the depressed pacemaker activity of sinoatrial node myocytes [68]; ii) the gradual augmentation of tonic sympathetic activity as measured from post-ganglionic sympathetic nerves directed to skeletal muscles [59], [60]; iii) the progressive increase of norepinephrine concentrations [36], [60], [69]; iv) the continuing decline of vagal modulation as assessed from the amplitude of respiratory sinus arrhythmia in the time or frequency domain [8], [36], [65], [70]; v) the gradual alteration of the adrenoceptor function [71]; vi) the progressive diminution of the responsiveness of the sinus node to sympathetic outflow [36], [65], [72]; vii) the steady decrease of BRS [36], [60], [65], [73]. The association between the BRS decline and the increase of the SAP variance with age observed at REST stresses the buffering role of baroreflex.
During STAND we confirm the positive dependence of HP mean and the negative relation of HP variance on age [36], [66], the positive correlation of SAP mean with age [73] and the lack of a linear relation between SAP variance and age [36], [73]. During STAND, at difference of REST, the BRS decline with age was not associated to a progressive increase of the SAP variance. This finding might be the consequence of the simultaneous decrease of the intensity of the causal relation from HP to SAP with age, thus leading to a situation of progressive HP-SAP uncoupling with age. These results have been explained by the reduced effect of the postural maneuver on the cardiovascular variables due to the diminished responsiveness of the sinus node to neural inputs in response to stressors [36], [65], [72], [74], to the reduced responsiveness of the vasculature to vasodilatator agents [75] and in reaction to stimuli [36], [65], [73], to the increase of peripheral resistances [65], and to the decreased cardiac baroreflex efficiency in response to the postural challenge [65].
Effect of Age on Complexity Indexes
This study confirms the gradual decrease of the complexity of the HP dynamics with age detected at REST using a univariate approach [4], [6]–[9]. Since this finding was obtained by constructing a multivariate embedding space, we conclude that accounting for series more closely related to vascular and respiratory systems (i.e. SAP and RESP) was inessential for the quantification of complexity of the HP variability [11]. This finding appears to be robust because it was detected by all the proposed approaches. As a new finding STAND was associated with a progressive increase of the HP complexity with age. Since the HP complexity during STAND decreased [33], [76] as a result of the sympathetic activation and vagal withdrawal [73], [77]–[79], this finding indicates a reduced ability of the cardiovascular system to cope with the postural challenge with age. This finding confirms, at the level of the complexity of the cardiovascular control, the difficulty of old individuals to deal with sympathetic stressors such as exercise and/or orthostatic challenge [65], [80]. Therefore, we suggest the use of complexity analysis of HP series during orthostatic challenge to quantify the reduced ability of the cardiovascular control of elderly subjects to cope with stressors.
At REST the complexity of the SAP series was found to decrease with age. We speculate that the progressive increase of the sympathetic activity with age [59], [60] determines a decrease of the number of temporal scales observable in the SAP series in the low frequency band (i.e. from 0.04 to 0.15 Hz) by synchronizing mechanisms operating in this band (e.g. peripheral vasomotion) [61]. Conversely, the complexity of the SAP variability was unrelated to age during STAND. Since sympathetic activity progressively increased with age [59], [60] and STAND led to an additional sympathetic overactivation [77], [78], a progressive reduction of the SAP complexity could be expected during STAND. Since the expected progressive reduction was not found in this study, the lack of association between SAP complexity and age during STAND might be the result of the reduced responsiveness of vasculature to vasomotor sympathetic control with age [36], [65], [73].
Regardless of the experimental condition RESP complexity was unrelated to age. Therefore, it can be concluded that senescence does affect respiratory centers and brain stem level and, more generally, respiratory activity. Since respiration is a strong periodical function, even though with a certain degree of irregularity about the dominant frequency [33], it is not surprising to find out that this pattern is maintained during aging independently of the experimental condition (i.e. REST or STAND).
Effect of Age on Causality Indexes
At REST we did not find any linear relation of the strength of the causal link from SAP to HP on age. This finding suggests that the intensity of the causal relation from SAP to HP was preserved in old individuals and, thus, in elderly people the relation from SAP to HP (i.e. the cardiac baroreflex) continues to play an important role in regulating cardiovascular variables. This result was confirmed by BEI. This result might appear surprising at the first sight because we found that BRS gradually fell with age [36], [60], [81]. Nevertheless, it is worth noting that the decrease of the gain of the relation from SAP to HP does not necessary imply a diminished strength of the causal link from SAP to HP because gain and strength measure different aspects of a relation between variables [15]. For example, during graded head-up tilt BRS progressively decreased [77], while the strength of the causal relation from SAP to HP progressively increased [23], [63], [82]. On the contrary, during STAND we observed a gradual reduction of the intensity of the causal relation from SAP to HP with age, thus suggesting a progressive reduction of the efficiency of the cardiac baroreflex control with age that was unveiled by a maneuver challenging cardiac baroreflex regulation (i.e. STAND). This finding was again confirmed by BEI, thus remarking the association between this traditional index [40] and state-of-art causality markers.
The causal relation from HP to SAP, i.e. the so-called mechanical feedforward pathway, forms with the cardiac baroreflex feedback the closed loop regulating HP-SAP interactions. This closed loop control was found active both at REST and during STAND [15], [63]. Since the mechanical feedforward pathway relies on the opposite actions of Frank-Starling mechanism and diastolic runoff, with a prevalence of the diastolic runoff in humans [18], the reduction of the strength of the causal link from HP to SAP with age could be observed in presence of the progressive balancing of the two opposite influences during senescence. At REST we can hypothesize that Frank-Starling mechanism gains importance with age in regulating HP-SAP variability interactions due to the diminished importance of the cardiac baroreflex [65], thus reducing the dominance of the diastolic runoff. During STAND the main mechanism underpinning the decrease of the strength of the relation from HP to SAP with age might be the increase of peripheral resistances progressively reducing the importance of the diastolic runoff during senescence [65]. The decreased strength of the causal link from HP to SAP combined with the diminished intensity of the causal relation in the opposite temporal direction (i.e. from SAP to HP) observed during STAND led to a situation of progressive HP-SAP uncoupling with age. This finding stresses further the gradual loss of the ability to deal with orthostatic challenge with age.
Another relevant finding of this study is the progressive decrease of the strength of the causal link from RESP to HP with age at REST. Since the relation from RESP to HP is the result of the continuous modifications of the membrane potentials of preganglionic vagal motoneurones and their responsiveness to stimulatory inputs imposed by respiratory activity [22], [62], the gradual decrease of the strength of the link from RESP to HP is likely to be the consequence of the gradual increase of tonic sympathetic activity with age [59], [60], [83], of the progressive vagal withdrawal leading to the progressive reduction of the respiratory sinus arrhythmia [8], [36], [65], [70] and of the gradual uncoupling between respiratory activity and the heart. This finding corroborates recent observations based on a transfer entropy approach using a fixed embedding dimension [49] and on a method quantifying phase synchronization between nonlinear models of the cardiac and respiratory oscillators [84]. This trend was not detected during STAND, thus suggesting that the limitation of the respiratory sinus arrhythmia induced by the sympathetic activation evoked by STAND did mask this phenomenon.
A significant negative linear association between age and the strength of the causal link from RESP to SAP was detected both at REST and during STAND. In a previous study we found a significant causal relation from RESP to SAP at REST and we demonstrated that it was maintained during orthostatic challenge [82]. The causal relation from RESP to SAP is the result of the mechanical influences on venous return and on large vessels due to respiratory-related modifications of intrathoracic pressure [85]–[87]. The observed decrease of the strength of the causal link from RESP to SAP might be a direct consequence of a less efficient modification of the stroke volume and a reduced effect of the intrathoracic pressure changes on large vessels resulting from the increase of cardiac and vascular stiffness with age [83], [88] due to the progressive left ventricular wall and large elastic artery intimal media thickening [83], [88]. Since the modification of the stroke volume at the respiratory rate depends also on the amplitude of the respiratory sinus arrhythmia, its well-known decline with age might contribute to the reduced respiratory-related variation of the diastolic filling and, thus, to the decline of the respiratory-related SAP changes with age.
None of the approaches was able to detect a linear association between age and the strength of the causal relation from HP or SAP to RESP both at REST and during STAND. Since previous studies suggested that a causal link from SAP to RESP was extremely unlikely in short-term cardiovascular regulation at REST and during postural change [15], [82], the uncorrelation between age and the strength of the causal relation from SAP to RESP was expected. Conversely, previous studies detected a significant causal link from HP to RESP as a result of the quickness of the neural cardiac influences at the respiratory rate compared to the slowness of the respiratory-related changes assessed according to respiratory inductive plethysmography exploiting thoracic belts [82], [89]. Since fast respiratory-related HP variations decreased with age, the lack of association between age and the strength of the causal link from HP to RESP might be an unexpected result. This surprising finding might be the consequence of the preservation of the nerve conduction velocity and neural processing time with age regardless of the experimental conditions.
Limitations of the Study and Future Developments
The strength of causal relations was estimated according to the concept of Granger causality [12]. The Granger approach to the inference of causality has been chosen among others [90] due to its direct connection with the quantification of complexity based on predictability and information content. We advocate the practical exploitation of alternative definition of causality, including those implying an intervention (i.e. the physical manipulation of the cause to stimulate an effect) to elucidate on real data the dependence of the conclusions on the adopted paradigm and the difference between the assessment of causality based on spontaneous variability and physical interventions. Among possible extensions we encourage the use of techniques based on permutations of ordinal patterns [91], thus performing the analysis directly in the framework of symbolic dynamics. Since the null hypothesis of Normality of HP, SAP and RESP series was rejected in 41%, 54% and 88% at REST and in 46%, 34% and 82% during STAND a portion of the MB complexity might be the result of the relevant percentage of non-Normal series. The significant percentage of non-Normal distribution might have some impact on the MB causality indexes as well. Future efforts should be devoted to complete the automation of the MF analysis by proposing an automatic procedure for setting the few parameters, notably k, that still remain under control of the user. From a more physiological standpoint the description of the relation governing variability interactions should be rendered more complete by accounting for new, contemporaneously recorded, physiological quantities such as diastolic arterial pressure, left ventricular contractility, stroke volume and venous return. Moreover, instead of comparing young sedentary volunteers to elderly sedentary humans, future protocols should contrast the same control group with elderly athletes with preserved peak to clarify the role played by fitness deconditioning and separate it from that of senescence. Further applications should be devoted to make clear the role of diet in preserving complexity of cardiac regulation and the intensity of the causal relations among cardiovascular variabilities by comparing groups of people following different eating habits.
Conclusions
This study proposed a compact framework for the assessment of complexity of cardiovascular variability series and causality of their interactions. The approach provided quantitative indexes that have been demonstrated helpful in elucidating the effect of age on cardiovascular control in humans. For example, results suggest that mechanisms other than cardiorespiratory coupling and baroreflex are responsible for the decrease of complexity of cardiac regulation with age at REST, while the impairment of the baroreflex is responsible for the increase of the complexity of the cardiac control with age during STAND. The proposed markers might be exploited in future clinical applications addressing the issue of monitoring the aging process and assessing the performance of countermeasures to the senescence of the cardiovascular control. The framework is particularly powerful because it allows the direct comparison between linear MB and nonlinear MF approaches and the minimization of the number of parameters needed to be set to perform complexity and causality analyses. Therefore, from a methodological standpoint, it is helpful to understand whether the exploitation of MF methods provides additional insights compared to a traditional simpler linear MB approach, and, from a more applicative standpoint, it favors user-independent applications and improves reproducibility of the results. The framework clearly demonstrates the importance of the MF methods especially in presence of nonlinear dynamics and nonlinearities in the interactions among series, the better reliability of the MB approach when the complexity of the interactions needs high dimensional embedding spaces to be fully unfolded, the greater statistical power of the LP approach compared to the CE one in detecting the age-related modifications of cardiovascular control, and the importance of computing causality indexes together with the more traditional complexity markers in the context of the evaluation of the senescence of cardiovascular regulatory mechanisms.
Author Contributions
Conceived and designed the experiments: AP ACMT AMC. Performed the experiments: NMP JM VM ABS ACMT AMC. Analyzed the data: AP VB AM TB NMP JM VM ABS ACMT AMC. Contributed reagents/materials/analysis tools: AP VB AM TB. Wrote the paper: AP LF GN ACMT AMC.
References
- 1. Cohen MA, Taylor JA (2002) Short-term cardiovascular oscillations in man: measuring and modelling the physiologies. J Physiol 542: 669–683.
- 2.
Koepchen HP (1991) Physiology of rhythms and control systems: an integrative approach. In: Rhythms in physiological systems, edited by Haken H, Koepchen HP. Berlin: Springer-Verlag, 3–20.
- 3. Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology. Standard of measurement, physiological interpretation and clinical use (1996) Circulation. 93: 1043–1065.
- 4. Kaplan DT, Furman MI, Pincus SM, Ryan SM, Lipsitz LA, et al. (1991) Aging and the complexity of cardiovascular dynamics. Biophys J 59: 945–949.
- 5. Lipsitz LA (1999) Age-related changes in the “complexity” of cardiovascular dynamics: a potential marker of vulnerability in disease. Chaos 5: 102–109.
- 6. Pikkujamsa SM, Makikallio TH, Sourander LB, Raiha IJ, Puukka P, et al. (1999) Cardiac interbeat interval dynamics from Childhood to Senescence. Comparison of conventional and new measures based on fractals and chaos theory. Circulation 100: 393–399.
- 7. Takahashi ACM, Porta A, Melo RC, Quiterio RJ, da Silva E, el al (2012) Aging reduces complexity of heart rate variability assessed by conditional entropy and symbolic analysis. Intern Emerg Med 7: 229–235.
- 8. Beckers F, Verheyden B, Aubert AE (2006) Aging and nonlinear heart rate control in a healthy population. Am J Physiol 290: H2560–H2570.
- 9. Viola AU, Tobaldini E, Chellappa SL, Casali KR, Porta A, et al. (2011) Short-term complexity of cardiac autonomic control during sleep: REM as a potential risk factor for cardiovascular system in aging. PLoS ONE 6: e19002.
- 10. Goldberger AL, Amaral LAN, Hausdorff JM, Ivanov PCh, Peng C-K, et al. (2002) Fractal dynamics in physiology: alterations with disease and aging. Proc Natl Acad Sci 99: 2466–2473.
- 11. Porta A, Castiglioni P, di Rienzo M, Bari V, Bassani T, et al. (2012) Short-term complexity indexes of heart period and systolic arterial pressure variabilities provide complementary information. J Appl Physiol 113: 1810–1820.
- 12. Granger CWJ (1980) Testing for causality. A personal viewpoint. J Econ Dyn Control 2: 329–352.
- 13. Chen Y, Rangarajan G, Feng J, Ding M (2004) Analyzing multiple nonlinear time series with extended Granger causality. Phys Lett A 324: 26–35.
- 14. Faes L, Nollo G, Porta A (2011) Information-based detection of nonlinear Granger causality in multivariate processes via a nonuniform embedding technique. Phys Rev E 83: 051112.
- 15. Porta A, Castiglioni P, di Rienzo M, Bassani T, Bari V, et al. (2013) Cardiovascular control and time domain Granger causality: insights from selective autonomic blockade. Phil Trans R Soc A 371: 20120161.
- 16. Vicente R, Wibral M, Lindner M, Pipa G (2011) Transfer entropy - a model-free measure of effective connectivity for the neurosciences. J Comput Neurosci 30: 45–67.
- 17. Porta A, Faes L (2013) Assessing causality in brain dynamics and cardiovascular control. Phil Trans R Soc A 371: 20120517.
- 18. Baselli G, Cerutti S, Badilini F, Biancardi L, Porta A, et al. (1994) Model for the assessment of heart period and arterial pressure variability interactions and respiratory influences. Med Biol Eng Comput 32: 143–152.
- 19. Porta A, Baselli G, Rimoldi O, Malliani A, Pagani M (2000) Assessing baroreflex gain from spontaneous variability in conscious dogs: role of causality and respiration. Am J Physiol 279: H2558–H2567.
- 20. Mullen TJ, Appel ML, Mukkamala R, Mathias JM, Cohen RJ (1997) System identification of closed loop cardiovascular control: effects of posture and autonomic blockade. Am J Physiol 272: H448–H461.
- 21. Triedman JK, Perrott MH, Cohen RJ, Saul JP (1995) Respiratory sinus arrhythmia: time domain characterization using autoregressive moving average analysis. Am J Physiol 268: H2232–H2238.
- 22. Eckberg DL (2003) The human respiratory gate. J Physiol 548: 339–352.
- 23. Porta A, Bassani T, Bari V, Tobaldini E, Takahashi ACM, et al. (2012) Model-based assessment of baroreflex and cardiopulmonary couplings during graded head-up tilt. Comput Biol Med 42: 298–305.
- 24. Porta A, Gnecchi-Ruscone T, Tobaldini E, Guzzetti S, Furlan R, et al. (2007) Progressive decrease of heart period variability entropy-based complexity during graded head-up tilt. J Appl Physiol 103: 1143–1149.
- 25.
Soderstrom T, Stoica P (1988) System identification. Englewood Cliffs: Prentice Hall.
- 26. Akaike H (1974) A new look at the statistical novel identification. IEEE Trans Autom Contr 19: 716–723.
- 27. Abarbanel HDI, Carroll TL, Pecora LM, Sidorowich JJ, Tsimring LS (1994) Predicting physical variables in time-delay embedding. Phys Rev E 49: 1840–1853.
- 28. Theiler J (1986) Spurious dimension from correlation algorithms applied to limited time-series data. Phys Rev A 34: 2427–2432.
- 29. Sugihara G, May RM (1990) Non linear forecasting as a way of distinguishing chaos from measurement error in time series. Nature 344: 734–741.
- 30. Pincus SM (1995) Approximate entropy (ApEn) as a complexity measure. Chaos 5: 110–117.
- 31. Porta A, Guzzetti S, Furlan R, Gnecchi-Ruscone T, Montano N, et al. (2007) Complexity and non linearity in short-term heart period variability: comparison of methods based on local non linear prediction. IEEE Trans Biomed Eng 54: 94–106.
- 32. Vlachos I, Kugiumtzis D (2010) Nonuniform state-space reconstruction and coupling direction. Phys Rev E 82: 016207.
- 33. Porta A, Guzzetti S, Montano N, Pagani M, Somers VK, et al. (2000) Information domain analysis of cardiovascular variability signals: evaluation of regularity, synchronisation and co-ordination. Med Biol Eng Comput 38: 180–188.
- 34. Porta A, Castiglioni P, Bari V, Bassani T, Marchi A, et al. (2013) K-nearest-neighbor conditional entropy approach for the assessment of short-term complexity of cardiovascular control. Physiol Meas 34: 17–33.
- 35. Schreiber T (2000) Measuring information transfer. Phys Rev Lett 85: 461–464.
- 36. Barnett SR, Morin RJ, Kiely DK, Gagnon M, Azhar G, et al. (1999) Effects of age and gender on autonomic control of blood pressure dynamics. Hypertension 33: 1195–1200.
- 37. Balady GM, Arena R, Sietsema K, Myers J, Coke L, et al. (2010) Clinician’s Guide to cardiopulmonary exercise testing in adults: a scientific statement from the American Heart Association. Circulation 122: 191–225.
- 38. Neves LTM, Karsten M, Neves VR, Bertrame T, Borghi-Silva A, et al. (2011) Relationship between inspiratory muscle capacity and peak exercise tolerance in patients post-myocardial infarction. Heart Lung 41: 137–145.
- 39. Bertinieri G, di Rienzo M, Cavallazzi A, Ferrari AU, Pedotti A, et al. (1985) A new approach to analysis of the arterial baroreflex. J Hypertens 3: S79–S81.
- 40. di Rienzo M, Parati G, Castiglioni P, Tordi R, Mancia G, et al. (2001) Baroreflex effectiveness index: an additional measure of baroreflex control of heart rate in daily life. Am J Physiol 280: R744–R751.
- 41. Porta A, Bassani T, Bari V, Pinna GD, Maestri R, et al. (2012) Accounting for respiration is necessary to reliably infer Granger causality from cardiovascular variability series. IEEE Trans Biomed Eng 59: 832–841.
- 42. Eckberg DL (1976) Temporal response patterns of the human sinus node to brief carotid baroreceptor stimuli. J Physiol 258: 769–782.
- 43. Saul JP, Berger RD, Albrecht P, Stein SP, Chen MH, et al. (1991) Transfer function analysis of the circulation: unique insights into cardiovascular regulation. Am J Physiol 261: H1231–H1245.
- 44. Richman JS, Moorman JR (2000) Physiological time-series analysis using approximate entropy and sample entropy. Am J Physiol 278: H2039–H2049.
- 45. Porta A, Baselli G, Liberati D, Montano N, Cogliati C, et al. (1998) Measuring regularity by means of a corrected conditional entropy in sympathetic outflow. Biol Cybern 78: 71–78.
- 46. Farmer JD, Sidorowich JJ (1987) Predicting chaotic time series. Phys Rev Lett 59: 845–848.
- 47. Porta A, Baselli G, Guzzetti S, Pagani M, Malliani A, et al. (2000) Prediction of short cardiovascular variability signals based on conditional distribution. IEEE Trans Biomed Eng 47: 1555–1564.
- 48. Hlavackova-Schindler K, Palus M, Vejmelka M, Bhattacharya J (2007) Causality detection based on information-theoretic approaches in time series analysis. Phys Rep 441: 1–46.
- 49. Nemati S, Edwards BA, Lee J, Pittman-Polletta B, Butler JP, et al. (2013) Respiration and heart rate complexity: effects of age and gender assessed by band-limited transfer entropy. Resp Physiol Neurobi 189: 27–33.
- 50. Uddin Ahmed M, Mandic DP (2012) Multivariate multiscale entropy analysis. IEEE Signal Proc Lett 19: 91–94.
- 51. Xiao X, Mullen TJ, Mukkamala R (2005) System identification: a multi-signal approach for probing neural cardiovascular regulation. Physiol Meas 26: R41–R71.
- 52. Porta A, Aletti F, Vallais F, Baselli G (2009) Multimodal signal processing for the analysis of cardiovascular variability. Phil Trans R Soc A 367: 391–408.
- 53. Blinowska KJ (2011) Review of the methods of determination of directed connectivity from multichannel data. Med Biol Eng Comput 49: 512–529.
- 54. Schelter B, Timmer J, Eichler M (2009) Assessing the strength of directed influences among neural signals using renormalized partial directed coherence. J Neurosci Methods 179: 121–130.
- 55. Takahashi DY, Baccalà L, Sameshima K (2007) Connectivity inference between neural structures via partial directed coherence. J Appl Stat 34: 1259–1273.
- 56. Borst C, Karemaker JM (1983) Time delays in the human baroreceptor reflex. J Auton Nerv Syst 9: 399–409.
- 57. Fraser AM, Swinney HL (1986) Independent coordinates from strange attractors from mutual information. Phys Rev A 33: 1134–1140.
- 58. Runge J, Heitzig J, Petoukhov V, Kurths J (2012) Escaping the curse of dimensionality in estimating multivariate transfer entropy. Phys Rev Lett 108: 258701.
- 59. Seals DR, Esler MD (2000) Human ageing and sympathoadrenal system. J Physiol 528: 407–417.
- 60. Parker Jones P, Christou DD, Jordan J, Seals DR (2003) Baroreflex buffering is reduced with age in healthy men. Circulation 107: 1770–1774.
- 61. Baselli G, Porta A, Pagani M (2006) Coupling arterial windkessel with peripheral vasomotion: modeling the effects on low-frequency oscillations. IEEE Trans Biomed Eng 53: 53–64.
- 62. Gilbey MP, Jordan D, Richter DW, Spyer KM (1984) Synaptic mechanisms involved in the inspiratory modulation of vagal cardio-inhibitory neurones in the cat. J Physiol 356: 65–78.
- 63. Porta A, Catai AM, Takahashi ACM, Magagnin V, Bassani T, el al (2011) Causal relationships between heart period and systolic arterial pressure during graded head-up tilt. Am J Physiol 300: R378–R386.
- 64. Nollo G, Faes L, Porta A, Pellegrini B, Ravelli F, et al. (2002) Evidence of unbalanced regulatory mechanism of heart rate and systolic pressure after acute myocardial infarction. Am J Physiol 283: H1200–H1207.
- 65. Laitinen T, Niskamen L, Geelen G, Lansimies E, Hartikainen J (2004) Age dependency of cardiovascular autonomic responses to head-up tilt in healthy subjects. J Appl Physiol 96: 2333–2340.
- 66. O’Brien IAD, O’Hare P, Corrall RJM (1986) Heart rate variability in healthy subjects: effect of age and the derivation of normal ranges for tests of autonomic function. British Heart Journal 55: 348–354.
- 67. Laitinen T, Hartikainen J, Niskamen L, Geelen G, Lansimies E (1999) Sympathovagal balance is major determinant of short-term blood pressure variability in healthy subjects. Am J Physiol 276: H1245–H1252.
- 68. Larson ED, St ClairJR, Summer WA, Bannister RA, Proenza C (2013) Depressed pacemaker activity of sinoatrial node myocytes contributes to the age-dependent decline in maximum heart rate. Proc Natl Acad Sci 110: 18011–18016.
- 69. Ziegler MG, Lake CR, Kopin IJ (1976) Plasma noradrenaline increases with age. Nature 261: 333–334.
- 70. Hrushesky WJM, Fader D, Schmitt O, Gilberstsen V (1984) The respiratory sinus arrhythmia: a measure of cardiac age. Science 224: 1001–1004.
- 71. Kelly J, O’Malley K (1984) Adrenoceptor function and ageing. Clin Sci 66: 509–515.
- 72. Lakatta EG (1993) Cardiovascular regulatory mechanisms in advanced age. Physiol Rev 73: 413–465.
- 73. Veermann DP, Imholz BPM, Wieling W, Karemaker JM, van Mortfrans GA (1994) Effects of aging on blood pressure variability in resting conditions. Hypertension 24: 120–130.
- 74. Esler MD, Thompson JM, Kaye DM, Turner AG, Jennings GL, et al. (1995) Effects of aging on the responsiveness of the human cardiac sympathetic nerves to stressors. Circulation 91: 351–358.
- 75. Elliott HL, Summer DJ, McLean K, Reid JL (1982) Effect of age on the responsiveness of vascular alpha-adrenoceptors in man. J Cardiovasc Pharmacol 4: 388–392.
- 76. Turianikova Z, Javorka K, Baumert M, Calkovska A, Javorka M (2011) The effect of orthostatic stress on multiscale entropy of heart rate and blood pressure. Physiol Meas 32: 1425–1437.
- 77. Cooke WH, Hoag JB, Crossman AA, Kuusela TA, Tahvanainen KUO, et al. (1999) Human responses to upright tilt: a window on central autonomic integration. J Physiol 517: 617–628.
- 78. Furlan R, Porta A, Costa F, Tank J, Baker L, et al. (2000) Oscillatory patterns in sympathetic neural discharge and cardiovascular variables during orthostatic stimulus. Circulation 101: 886–892.
- 79. Javorka M, Trunkvalterova Z, Tonhajzerova I, Javorkova J, Javorka K, et al. (2008) Short-term heart rate complexity is reduced in patients with type 1 diabetes mellitus. Clin Neurophysiol 119: 1071–1081.
- 80. Rodeheffer RJ, Gerstenblith G, Becker LC, Fleg JL, Weisfeldt ML, et al. (1984) Exercise cardiac output is maintained with advancing age in healthy human subjects: cardiac dilatation and increased stroke volume compensate for a diminished heart rate. Circulation 69: 203–213.
- 81. Laitinen T, Hartikainen J, Vanninen E, Niskamen L, Geelen G, et al. (1998) Age and gender dependency of baroreflex sensitivity in healthy subjects. J Appl Physiol 84: 576–583.
- 82.
Porta A, Takahashi ACM, Catai AM, Montano N (2013) Assessing causal interactions among cardiovascular variability series through a time domain Granger causality approach. In: Methods in brain connectivity inference through multivariate time series analysis, edited by Baccalà L and Sameshima K. Boca Raton Ann Arbor London Tokyo: CRC Press, 223–242.
- 83. Lakatta EG, Levy D (2003) Arterial and cardiac aging: major shareholders in cardiovascular disease enterprises. Part II: the aging heart in health: links to heart disease. Circulation 107: 346–354.
- 84. Iatsenko D, Bernjak A, Stankovski T, Shiogai Y, Owen-Lynch PJ, et al. (2013) Evolution of cardiorespiratory interactions with age. Phil Trans R Soc A 371: 20110622.
- 85. Caiani E, Turiel M, Muzzupappa S, Porta A, Baselli G, et al. (2000) Evaluation of respiratory influences on left ventricular function parameters extracted from echocardiographic acoustic quantification. Physiol Meas 21: 175–186.
- 86. Innes J, De Cort S, Kox W, Guz A (1993) Within-breath modulation of left ventricular function during normal breathing and positive-pressure ventilation in man. J Physiol 460: 487–502.
- 87. Toska K, Eriksen M (1993) Respiration-synchronous fluctuations in stroke volume, heart rate and arterial pressure in humans. J Physiol 472: 501–512.
- 88. Lakatta EG, Levy D (2003) Arterial and cardiac aging: major shareholders in cardiovascular disease enterprises. Part I: aging arteries: a “set up” for vascular disease. Circulation 107: 139–146.
- 89. Yana K, Saul JP, Berger RD, Perrott MH, Cohen RJ (1993) A time domain approach for the fluctuation analysis of heart rate related to instantaneous lung volume. IEEE Trans Biomed Eng 40: 74–81.
- 90. Eichler M (2013) Causal inference with multiple time series: principles and problems. Phil Trans R Soc A 371: 20110613.
- 91. Riedl M, Muller A, Wessel N (2013) Practical considerations of permutation entropy. Eur Physl J Special Topics 222: 249–262.