After approximately 3 years from the baseline assessment, in order to evaluate cognitive change, participants were invited to attend a second neuropsychological evaluation using the same battery of neuropsychological tests employed at baseline. A geriatric assessment was also carried out with the same instruments used to measure functional and global cognitive status at baseline, i.e., the BADL, IADL and MMSE scales.

MCI was diagnosed as previously described. Dementia was diagnosed by means of standard criteria: the National Institute on Aging-Alzheimer’s association (NIA-AA) criteria for all-cause and Alzheimer’s disease (AD) dementia (McKhann et al., 2011) and the American Heart Association/American Stroke Association (AHA/ASA) criteria for vascular dementia (VAD) (Gorelick et al., 2011). Further details on the diagnosis of dementia are given in the Supplementary Materials and Methods.

It should be noted that, although a diagnosis of cognitive status was also made at follow-up, in our study the reference for the diagnosis of cognitive status was taken to be the same time of the HRV assessment, i.e., baseline. Since at baseline subjects with NC and MCI were enrolled in the study while subjects with dementia were excluded from it, our hypotheses were tested only in the two groups of subjects diagnosed with NC and MCI at baseline (also see Limitations).

The study was conducted in accordance with the Declaration of Helsinki and was approved by the ethics committee of the Fondazione IRCCS Ca’ Granda Ospedale Maggiore Policlinico in Milan, Italy. All participants gave written informed consent to participation in the study.

Data are reported as mean (standard deviation) for continuous variables and as number (percentage) for categorical variables.

In the two cognitive groups at baseline (NC and MCI) the neuropsychological test scores were compared at baseline and follow-up by means of the paired t-test or the Wilcoxon sign-ranked test as appropriate.

Linear regression was performed in the NC and MCI groups at baseline to investigate if the autonomic response at baseline predicted cognitive change from baseline to follow-up. The annualized cognitive change score (see later) was taken as dependent variable while the independent variables were the Δ HRV indices (the variables of interest) as well as potential confounders which were selected a priori from the literature according to their being known risk factors for cognitive decline: socio-demographic factors (age, sex and education), physical activity and physical/mental comorbidity (e.g., Roberts and Knopman, 2013). Adjustment was also carried out for baseline cognition and resting HRV. Both simple (unadjusted) and multiple (adjusted) linear regression were conducted.

In accordance with the study hypotheses, two different regression models were fitted for each of the NC and MCI groups: one with the episodic memory change score as the outcome and the sympathetic Δ HRV indices (i.e., active standing Δ LFn and Δ LF/HF) among the predictors (model 1), the other with the executive functioning change score as the outcome and the parasympathetic Δ HRV indices (i.e., paced breathing Δ LFn and Δ LF/HF) among the predictors (model 2).

Composite scores were computed for episodic memory and executive functioning by means of a Z-score transformation. The raw scores of the individual neuropsychological tests were converted to Z-scores based on the mean and standard deviation of published normative data for the Italian population (see Supplementary Table 1) and then averaged to yield a domain-specific score. Scores quantifying response time (i.e., Trail-Making Tests A and B) or number or errors (i.e., Cognitive Estimates total and bizarre) were multiplied by −1 so that lower Z-scores consistently indicated poorer performance. We chose to use such composite cognitive scores because they offer several advantages: they minimize the likelihood of type I error associated with multiple testing, have greater reliability and reduce floor and ceiling effects (Morris et al., 1999; Williams et al., 2019). The cognitive change score for each domain was calculated by subtracting the baseline Z-score from the follow-up Z-score and it was then annualized by dividing it by the number of years of follow-up, so as to account for different follow-up durations among the study participants.

As far as comorbidity was concerned, with the aim of achieving a more parsimonious model, physical and mental comorbidities were combined to create an additive index. The physical burden of illness was measured by the CIRS-m score. Psychological distress was quantified by the GDS-s and STPI-TA scores. For comparability, all three scores were rescaled between 0 and 1 using a simple linear stretch according to the formula X’ (rescaled score) = [X (original score) – X min (minimum score value)]/[X max (maximum score value) – X min], and were then summed. The choice to standardize these measures to a common metric by a minimum-maximum normalization rather than a Z-score transformation, as was the case for the cognitive scores, is explained hereafter. First, normative neuropsychological test data for the Italian population are well-established and are the cornerstone of neuropsychological assessment [see e.g., Spinnler and Tognoni (1987) in Supplementary Table 1]. On the contrary, normative data for the CIRS-m, GDS-s and STPI-TA scores are either lacking (e.g., for the CIRS-m) or limited by representativeness issues [e.g., norms for the GDS-s are available for mentally and physically healthy older adults in New Zealand (Knight et al., 1983); norms for the STPI-TA are available for a different 10-item form of the STPI-TA in a French study (Bergua et al., 2016)]. Second, although the current study sample could have been taken as a reference for norming, we believed this would have been a less appropriate approach given its relatively small size. Nonetheless, if the three scales are aggregated by using a Z-score transformation based on our sample data, the pattern of significance is unaffected (see Supplementary Tables 2, 3).

In order to meet the assumptions of normality and/or homoscedasticity, two outliers (standardized residual > |3|) (Stevens, 2012) were omitted from the regression analyses: one in the NC group for model 2 and one in the MCI group for model 1. None of the omitted (or included) data points were influential (Cooke’s distance < 1) (Stevens, 2012). Indeed, retaining all data and repeating the analyses with non-parametric testing (i.e., Spearman’s simple and partial correlations) does not substantively change the results despite some expected loss of statistical power (Sedgwick, 2015; see Supplementary Tables 4, 5).

The regression models satisfied the assumptions of normality (Shapiro–Wilk’s test, p > 0.05), homoscedasticity (Breusch Pagan test, p > 0.05), independence of errors (Durbin Watson values between 1.5 and 2.5) and no significant multicollinearity [variance inflation factor (VIF) < 5] (Stevens, 2012).

Inflation of type I error due to multiple testing was controlled for by means of the Benjamini–Hochberg procedure (Benjamini and Hochberg, 1995) with the false discovery rate (FDR) set at the conventional level for alpha (5%).

Analyses were performed by means of the statistical packages SPSS version 25.0 (SPSS Inc., Chicago, IL, United States) and R version 3.5.2 (The R Foundation for Statistical Computing, Vienna, Austria) for Windows. A P value ≤ 0.05 was considered statistically significant.

Due to the novelty of the original study no power calculation was performed a priori (Nicolini et al., 2014). However, we conducted a power determination analysis for the current sample size with G*power (Faul et al., 2009) based on the data of Mahinrad et al. (2016) from a similarly aged population. Cohen’s d for differences in cognitive test scores (Stroop, Letter-Digits Coding and Picture-Word Learning total) between the high and low HRV tertiles (d = 2.7–3.1) was converted to Cohen’s f² (f² = 1.8–2.4) according to Cohen (1988) so as to obtain an effect size estimate for linear regression. Both measures indicated large effect sizes and, in order to be conservative, the lowest f² value (f² = 1.8) was used in the power analysis. The achieved sample sizes (n = 36 NC, n = 33 MCI) yielded an almost 100% power for multiple linear regression with a 5% alpha level (two tailed) and eight predictors.

The observation that an increased autonomic response to a parasympathetic challenge predicted a slower decline in executive functioning generally aligns with a large body of cross-sectional studies reporting an association between parasympathetic HRV indices and executive functioning (e.g., Williams et al., 2016, 2019; Forte et al., 2019a). More specifically, it conforms well to findings from the sparse longitudinal studies available in the literature. In particular, the standard deviation of the normal to normal intervals (SDNN), which is largely dependent on the parasympathetic nervous system when recorded under resting conditions (Malik et al., 1996), has been shown to have an inverse relationship with decline in the Letter-Digit Coding Test over a 3-year follow-up (Mahinrad et al., 2016) as well as a direct relationship with performance on the Stroop Test 5 years later (Zeki Al Hazzouri et al., 2017) and the Digit-Symbol Coding Test 10 years later (Schaich et al., 2020). Also, a greater parasympathetic responsivity to cognitive tasks, in terms of HF recovery and reactivity, has been demonstrated to predict an attenuated decline in executive functioning (and, to a lesser degree, in episodic memory) across a 9-year period (Knight et al., 2020). Although this latter study was unable to capture an association between the HF reactivity to an orthostatic challenge and changes in cognitive functioning, two points should be noted. First, the orthostatic stress was less intense since it involved standing from a sitting rather than a supine position. Second, there was no recovery period following the orthostatic challenge and no “true” baseline since HF reactivity was computed as HF during standing minus HF during the recovery epoch of the cognitive challenges. Given that mental stress produces parasympathetic withdrawal (Castaldo et al., 2015), as does orthostatic stress, such methodological choice is bound to reduce the magnitude of the orthostatic HF response and thus the likelihood of detecting any potential association it may have with cognition. Furthermore, even if Britton et al. (2008) concluded for no association between HRV and cognitive impairment, due to the lack of consistency of the observed associations across cognitive domains, they did identify a relationship between lower SDNN and HF and a greater 5-year decrease in the Mill Hill Test, which can be considered a measure of executive functioning (Raven, 1983). Lastly, a very recent study by Costa et al. (2021) has found that greater decline in executive functioning tests like the Digit-Symbol Coding and the Digit Span Forwards and Backwards over 6 years was predicted by greater heart rate fragmentation (HRF). This is a novel HRV metric whose physiological underpinnings are still unresolved, but are likely reflective of a degradation of the parasympathetic nervous system (Costa et al., 2017). The authors do not find an association between cognitive decline and traditional parasympathetic HRV indices [SDNN, the root mean square of successive differences of the normal to normal intervals (RMSSD) and HF], but the use of ECG recordings from in-home overnight polysomnography from an ancillary sleep study is an approach that, albeit convenient, carries some limitations. The presence of sleep-disordered breathing, highly prevalent in the study population (Chen et al., 2015), can have a confounding effect on HRV analysis both because it actually alters autonomic control and because abnormal breathing patterns distort the interpretation of the power spectrum (Tobaldini et al., 2013). The same would be true of leg movements during sleep, which are linked to autonomic activation as well as a potential source of motion artifacts (Guggisberg et al., 2007; Tobaldini et al., 2013). Moreover, despite the fact that the stationarity requirements of spectral analysis warrant the averaging of values from 5-min windows over the entire 12-h period, such averaging obscures detailed information on autonomic modulation (Malik et al., 1996).

The finding that an increased autonomic response to a sympathetic challenge predicted a faster decline in episodic memory was contrary to our expectations. Indeed, sparse cross-sectional studies have shown that better episodic memory is associated with enhanced sympathetic activation, as indexed by increased sympathetic HRV or decreased parasympathetic HRV. In fact, LF/HF has been found to positively correlate with the episodic memory subdomain of the Montreal Cognitive Assessment (MOCA) score in older adults (Frewen et al., 2013), orthostatic Δ LFn and Δ LF/HF have been reported to positively correlate with the Prose Delayed Recall Z-score in subjects with amnestic MCI (Nicolini et al., 2020), and a stress-induced increase in SDNN has been noticed to positively correlate with the California Verbal Learning Task in healthy individuals of different ages (Hilgarter et al., 2021). Also, resting TP, which is mostly parasympathetically-mediated (Malik et al., 1996), and HF, which is a parasympathetic index, have been found to negatively correlate with Word Recall and Word/Picture Recognition in late-life depression (Vasudev et al., 2012) and with the Seoul Verbal Learning Test in subjects with MCI/dementia (Kim et al., 2018). Along the same lines, Dalise et al. (2020) have described a positive adjusted association between LF/HF and the MMSE score in geriatric outpatients representative of a real-life setting, and it is recognized that, especially at older ages, the MMSE primarily captures episodic memory (Soubelet and Salthouse, 2011; Dichgans and Leys, 2017).

There could be different potential explanations for the association between the sympathetic response and the change in episodic memory being in the opposite direction than anticipated.

First, it could be a paradoxical effect stemming from survivor bias. In fact, low HRV has been demonstrated to increase mortality in older adults (Nicolini et al., 2012) and individuals with higher HRV, by living longer, could have a greater opportunity to develop cognitive impairment. However, this interpretation seems unlikely since, in our study, loss to follow-up was small (11%) and well below the 20–40% threshold recommended for cohort studies (Ramirez-Santana, 2018). Also, sympathetic Δ HRV indices were not significantly different between the mortality and no-mortality groups (see Supplementary Table 14).

Second, it could reflect a regression to the mean phenomenon by which, on repeated neuropsychological testing, more extreme scores at baseline tend to converge toward the mean at follow-up, so that there is a negative association between baseline and change cognitive scores. If this is coupled with a positive association between baseline sympathetic Δ HRV indices and baseline episodic memory, the former would then have a negative longitudinal relationship with change in episodic memory, i.e., greater baseline sympathetic Δ HRV indices would be linked to greater episodic memory decline at follow-up. Nonetheless, this does not appear to be plausible because there was no significant cross-sectional association between the Δ HRV indices and cognitive tests at baseline (see Supplementary Tables 15, 16 and the Supplementary Discussion) and because, in any case, regression to the mean was controlled for by adjusting the regression analyses for baseline cognitive scores (Yu and Chen, 2015).

Therefore, it is reasonable to assume that our results do not represent a statistical artifact, but, rather, have a true pathophysiological basis. Support for this view comes from a wealth of functional neuroimaging studies that have primarily focused on individuals with MCI (O’Brien et al., 2010; Sperling et al., 2010; Huijbers et al., 2015), mostly of the amnestic type, but have also spanned the cognitive continuum to include individuals with mild AD dementia (Alsop et al., 2008; Sperling et al., 2010) as well as cognitively normal older adults at genetic risk (Sperling et al., 2010; Rao et al., 2015) or not (Leal et al., 2017) for AD. They have revealed that in these subjects, relatively to controls, there is mostly a hyperactivity of the hippocampal/parahippocampal region (Alsop et al., 2008; O’Brien et al., 2010; Sperling et al., 2010; Huijbers et al., 2015; Rao et al., 2015) and, additionally, that hyperactivity is a longitudinal predictor of greater cognitive decline (O’Brien et al., 2010; Sperling et al., 2010; Huijbers et al., 2015; Rao et al., 2015; Leal et al., 2017), especially in episodic memory (O’Brien et al., 2010; Rao et al., 2015; Leal et al., 2017). These data fit in nicely with the cognitive reserve theory (Stern, 2012) according to which the brain can compensate for damage by engaging neural networks. Within this context, hyperactivity is viewed as a compensatory attempt to maintain cognitive performance in the face of cerebral pathology (Sperling et al., 2010; Stern, 2012). Likewise, hyperactivity is a harbinger of cognitive decline because it signifies that the system is working closer to its maximum capacity, i.e., closer to the threshold beyond which reserve is depleted and compensation is no longer possible (Sperling et al., 2010; Stern, 2012).

Somewhat more recently, similar evidence has emerged for the insula. Hyperactivity of the insular cortex has been demonstrated in subjects with MCI (Scarmeas et al., 2004; Alexopoulos et al., 2012; Alfini et al., 2019), mild AD dementia (Scarmeas et al., 2004) and cognitively normal older adults at genetic risk of AD (Thambisetty et al., 2010; Lin et al., 2017). Also, in subjects with MCI, insular hyperactivity has been associated with greater cognitive decline at follow-up (Devanand et al., 2006).

Although there is little functional (and even structural) neuroimaging work targeting the locus coeruleus, due to the technical challenges of imaging a small brainstem nucleus (Mather, 2021), other lines of research converge toward the notion that it, too, may exhibit an inverse U-shaped activation curve in AD (e.g., Ross et al., 2015). Thus, noradrenergic transmission would be aberrantly increased in the earlier stages and then decline later on (e.g., Ross et al., 2015).

On the whole, it can be conjectured that a greater HRV response to a sympathetic challenge at baseline could be indexing an increased activity of those cerebral regions that generate sympathetic outflow (i.e., the hippocampus/parahippocampus, insula and locus coeruleus) and that such hyperactivation may be signaling the impeding breakdown of compensatory mechanisms. Accordingly, as brain damage progresses, a point will be reached during the follow-up when compensation can no longer occur and there is a rapid deterioration of cognitive performance in the domains underpinned by these brain structures, i.e., episodic memory. It should be noted that the negative (and not positive) longitudinal association between the response to a parasympathetic challenge and executive decline would appear to imply that compensation is not operating here, and a potential reason for this is given the Supplementary Discussion.

The study has strengths and limitations. The study strengths include its novelty as well as the performance of a detailed and ad hoc neuropsychological, autonomic and clinical assessment including extensive neuropsychological testing, the use of two physical maneuvers to differentially challenge the two branches of the ANS and the collection of information on several potential confounders.

A number of limitations must also be acknowledged. Since HRV was not measured at follow-up it was not possible to track the trajectory of autonomic function over time. This would have enabled us to more thoroughly examine the relationship between HRV and cognition by determining how change in HRV associates with change in cognition (Elias and Torres, 2017), an issue which still remains unexplored. Thus, future longitudinal research should include concurrent HRV and cognitive measures.

The sample size was relatively small. Although an a posteriori power analysis with a large expected effect size revealed an almost 100% power for multiple linear regression (see “Statistical analysis” section), only in the MCI group the observed effect size was large according to Acock (2014) (standardized beta coefficient > 0.5). Therefore, the study may have been underpowered to detect smaller effect sizes in the NC subjects and our findings will need to be replicated in larger cohorts.

Only the HRV indices that were the main focus of the study (i.e., Δ LFn and Δ LF/HF) were found to be significant predictors of cognitive change. No meaningful results were observed for other Δ HRV indices (time- and frequency-domain) or for all resting HRV indices. This was probably the case because it is generally recognized that, in short-term ECG recordings performed in controlled experimental conditions, there is a greater sensitivity of the spectral and Δ indices (Malik and Camm, 1990; Montano et al., 2009; Malik et al., 2019). The fact that this was not true of other Δ spectral indices (i.e., Δ HF, Δ TP and Δ LF) is likely due to their statistical properties and/or physiological significance. Although HF is an index of parasympathetic modulation it is also a highly dispersed variable (Nunan et al., 2010). TP and LF are indices of joint sympathetic and parasympathetic modulation, which implies that, since sympathetic activity is related to episodic memory and parasympathetic activity is related to executive functioning, a “mixed” measure is less able to capture these domain-specific associations, because their strength is diluted by the presence of a second autonomic component unrelated to the specific cognitive domain. It cannot be excluded that with a larger sample size and/or PSA techniques based on the Welch method, which reduces HRV variance (Stoica and Moses, 2004), our findings could have been different. Also, it should be noted that, even if transformed HRV indices (i.e., LFn and LF/HF) have received ample support in the literature (e.g., Malliani et al., 1998; Montano et al., 2009; Malik et al., 2019) and have been extensively employed in HRV studies across topics (e.g., Niu et al., 2018; Dalise et al., 2020; Sun et al., 2021), they have also been questioned (e.g., Heathers, 2014). Future works should address other measures of sympathetic activity such as plasma adrenaline/noradrenaline levels (Goldstein et al., 2011) and the pre-ejection period (Wiley et al., 2021).

The nature of the sample (i.e., outpatients referred to a Geriatric Unit) can limit the generalizability of our study to the older community-dwelling population at large, and the exclusion of individuals with dementia at baseline implies that HRV was not investigated as a predictor of cognitive decline in this subset of older adults.

Lastly, although the role of HRV as indicator of brain activity is grounded in much literature (see “Introduction,” Supplementary Introduction, and “Discussion” sections), the study lacked functional neuroimaging which could have provided deeper insight into the complexity of the relationship between HRV, cognition and the activity of specific CAN regions.