Voxel-wise seed based correlation analyses (SCA) were performed using the Conn toolbox (15.h). For each subject and each seed region, z-transformed Pearson correlation maps were brought to second level group analysis. We defined a priori seeds in brain regions that in previous literature have been influenced by aerobic exercise. These included (a) left and right hippocampus obtained from the SPM Wake Forest University (WFU) Pickatlas toolbox (Maldjian et al., 2003), (b) left and right parahippocampus (WFU pickatlas), and (c) PCC, based on a cluster reported by Boraxbekk et al. (2016) (center of gravity MNI coordinates -6, -32, 27; 67 voxels).

Variations in the preprocessed and cleaned resting state BOLD time series were quantified using FSL maths¹, rendering one voxel-wise whole brain map of the standard deviations (BOLD_STD) for each subject and each session. As a corroborating measure, we also calculated Amplitude of Low Frequency Fluctuations (ALFF) which quantifies the total power within the lower frequency band. ALFF was calculated for the same preprocessed and cleaned volumes as for the BOLD_STD calculations, using the MatLab toolbox DPARSFA² (Chao-Gan and Yu-Feng, 2010). Thus, for the examined ALFF calculations, the examined frequency window was essentially the same as the bandpass filter used in the resting state data, i.e., 0.008–0.09 Hz.

For the network analyses, we created subject specific weighted graphs based on 264 nodes defined as spherical (4 mm radius) ROIs centered around functionally relevant coordinates (Power et al., 2011; Cole et al., 2013) (see Supplementary Table S1). Edges were defined as Pearson correlation values between each time series of each ROI–ROI pair, where negative correlations were nulled. For measures of hippocampal betweenness centrality, we added two nodes corresponding to the left and right hippocampus (see Seed Based Correlation Analysis) to the original (264 nodes) graph. A high degree of betweenness centrality imply that nodes have central positions, i.e., act as hubs in the given network. Global efficiency is a measure of how efficiently information is exchanged across a network. All graph theoretical measures (i.e., Betweenness Centrality for left and right hippocampus, and global efficiency) were calculated using the brain connectivity toolbox (Rubinov and Sporns, 2010) as implemented in GraphVar1.0 (Kruschwitz et al., 2015).

Independent component analysis (ICA) was performed using the conn toolbox (v16, Whitfield-Gabrieli and Nieto-Castanon, 2012). Comparable ICA’s across subjects were obtained through a dual regression procedure, where we first calculated 20 group level spatial independent components (default value), and subsequently used the time series pertaining to each component to obtain subject-specific spatial components. Longitudinal change in five different networks were analyzed with regard to baseline and change in aerobic capacity, respectively. The ICA components of interest were chosen in an attempt to conceptually replicate the analysis of resting state networks investigated in Voss et al. (2016). Specifically, these included the default mode network (DMN), right and left frontoparietal network (FPN), dorsal attention network (DAN), and Sensorimotor Network (S1M1) (see Supplementary Figure S1).

We also conducted NBS analysis (Zalesky et al., 2010), on the graph consisting of the same 264 nodes as used for graph analysis (see Graph Analyses). NBS allows for identification of group-effects among ROI–ROI connections or subnetworks, while controlling for multiple comparisons. NBS is the conceptual equivalent to cluster statistics for voxel wise whole brain analyses. Analogous to cluster statistics, NBS yields stronger statistical sensitivity for distributed effects, which in turn is traded for spatial resolution which is limited to the size of the identified network (c.f. cluster), rather than individual ROI–ROI connections (c.f. voxels). NBS was used to identify ROI–ROI pairs or subnetworks that revealed group by time interactions, or where connectivity related to aerobic capacity (either at baseline, or in change-change regressions). For this, we employed a typically used network defining ROI–ROI (“primary”) connectivity threshold (Zalesky et al., 2010); p = 0.001 uncorrected. In order to detect weaker but spatially more distributed sub-graphs for 2nd level contrasts where the original ROI–ROI connectivity threshold did not show any effect, we also investigated fitness-network associations using ROI–ROI connectivity thresholds of p = 0.01 and p = 0.05 uncorrected (two sided). For NBS (i.e., “network component intensity”) threshold, we used the significance level of family-wise error correction (FWE) p = 0.05 (1000 permutations).

Additionally, we performed a principal component MVPA, or so called “connectome-MVPA,” to detect patterns of whole brain connectivity that correlated with aerobic capacity. The MVPA analysis complements the SCA and ICA since the investigated connectivity is not restricted to pre-selected seed regions or independent components, respectively. Instead, it provides a regionally unbiased mapping of brain areas with whole brain connectivity patterns that are predicted by aerobic capacity or longitudinal change thereof. In detail, the MVPA measure was obtained by dimension reduction of the whole brain connectivity matrix of each voxel. The connectivity matrix of each subject was reshaped into a row vector and subsequently concatenated over all participants into a matrix N × V, where N was the number of subjects and V is the number of voxels within the brain mask. The dimensionality of the N × V group correlation matrix was reduced by principal component analysis (PCA). This yielded an N × C matrix, where C is the number of maintained principal components. We maintained the first seven principal components that explained the most of the variance of the connectivity matrix (C = 7) [according to the rule of thumb to maintain an approximate 1:7 ratio between number of components and subjects, as proposed by the metric implementer Nieto-Castenon (2015)³]. In other words, the resulting seven component score volumes best represented the whole brain connectivity pattern for each subject. These volumes were included in an F-test at the 2nd level analysis. Thus, we tested for clusters that were predicted by aerobic capacity with regard to whole brain connectivity.

The reported activation maps (of contrasted regression parameter estimates) in group analyses were considered significant at a cluster-level significance of p < 0.05, false-discovery rate (FDR) corrected, and cluster defining voxel threshold of p < 0.001, two-sided. Due to unreliable estimations of smoothness of the MVPA rendered PCA maps, we used non-parametric permutation statistics (1000 permutations) as advocated by the implementer of the metric (Alfonso Nieto-Castanon, personal communication, September, 23, 2016). Connectome-MVPA analyses were thresholded using a cluster based significance threshold of p = 0.05 (FDR). Un-thresholded statistical maps for all reported results (Figures 1–3) are available online at http://neurovault.org/collections/2415/.

For all significant findings of group by time interactions, and associations between intrinsic brain activity either with regard to baseline aerobic capacity or change in aerobic capacity, we further examined the relationships to cognition. More specifically, we defined volumes of interest based on the brain regions that were significantly predicted by aerobic capacity. These volumes were subsequently used as inclusive masks (i.e., small volume correction) within which we investigated effects of either baseline score or longitudinal change in cognition, respectively.

Variation of the BOLD signal time series (BOLD_STD) was negatively correlated with aerobic capacity in three clusters located in midbrain regions at baseline (Figure 1 and Table 2). The largest cluster (973 voxels) extended over right posterior hippocampus and thalamus. A second cluster (298 voxels) covered parts of left pallidum, and to a smaller extent also left thalamus and amygdala. The third cluster (232 voxels) contained right anterior hippocampus and right amygdala. Clusters were anatomically labeled using the WFU pick atlas. BOLD_STD were not positively correlated with aerobic capacity in any brain region.

Amplitude of Low Frequency Fluctuations (ALFF), resulted in virtually identical results as BOLD_STD (see Supplementary Figure S2), and thus ALFF results are not presented further.

Follow-up analyses that tested for the association between cognitive performance and BOLD_STD were restricted to the search volume of regions where aerobic capacity predicted BOLD_STD (displayed in Figure 1). Given that particularly hippocampal and mid-temporal brain regions displayed fitness related BOLD_STD, we tested for EM. In addition, a meta-analysis (Colcombe and Kramer, 2003) found the largest association between improved fitness and cognitive performance to be related to EFs. Furthermore, in our previous study, improved fitness was primarily related to improved general cognitive performance (i.e., CS) (Jonasson et al., 2017). Therefore, we also tested the association between EF, CS and BOLD_STD. Neither EM, EF nor CS were significantly related to BOLD_STD. (All cognitive outcome measures are presented in detail in Jonasson et al., 2017).

None of the a priori defined seeds displayed significant VO₂-peak predicted connectivity at baseline. However, connectivity of the post hoc seed in the right medial temporal cortex (see BOLD-Signal Fluctuations), to frontal-, parietal-, and occipital regions displayed a positive correlation to aerobic capacity. We also detected a negative association between aerobic capacity and the connectivity between the same seed and right thalamus and right occipital cortex (Figure 2A and Table 3).

At baseline, among the five investigated ICA-derived networks we observed a negative relationship between VO₂-peak and the connectivity of a sensorimotor network and a brain region extending over calcarine sulcus and precuneus (see Figure 3A). No associations between VO₂-peak and connectivity of the other four networks were seen.

Several measures of intrinsic brain activity turned out not to be significantly related to aerobic capacity. For graph theoretical measures, neither global efficiency [t(42) = -0.06, p = 0.95] nor left or right hippocampal centrality [L: t(42) = -0.94, p = 0.34, R: t(42) = -1.23, p = 0.22] reached significance.

In the NBS analysis, no ROI–ROI connectivities or subnetworks were correlated with aerobic capacity, neither corrected nor uncorrected. Furthermore, the connectome-MVPA did not reveal any significant clusters predicted by baseline VO₂-peak.

Finally, analysis of ASL perfusion data did not reveal any associations between aerobic capacity and CBF from neither whole brain nor small volume voxel based analysis.

Longitudinal gain in aerobic capacity was positively associated with changes in functional connectivity between the BOLD_STD-derived seed region in right medial temporal lobe, and right angular gyrus and left inferior frontal gyrus (Figure 2B and Table 3). However, change in aerobic capacity was negatively associated with change in connectivity between a priori defined seed regions in the left hippocampus and the right precentral gyrus (Figure 2C and Table 3). None of the other a priori defined seeds related to change in VO₂-peak.

Behavioral follow-up analyses did not reveal any significant relationships between changes in the fitness related SCA results and EM, EF, or CS.

Independent component analysis identified two significant associations between changes in VO₂-peak and connectivity. Firstly, there was a positive relationship between gains in VO₂-peak and connectivity between DMN and left frontal pole (see Figure 3B and Table 3). Secondly, a negative relationship between change in VO₂-peak and the connectivity between the S1M1-network and right hypothalamus was observed (Figure 3C). There was no significant relationship to changes in cognition (EM, EF, and CS).

The NBS revealed a subnetwork constituted by connected nodes in the DMN, the frontal parietal-, visual- and ventral attention network (Figure 4) in which gains in fitness was related to changes in intra-network connectivity. Since the subgraph detected by NBS does not allow for inferences on individual edges (i.e., ROI–ROI connections), we also correlated change in VO₂-peak to connectivity between all individual ROI–ROI pairs (i.e., 34584 individual connections). In this analysis, no connections were significant (Bonferroni correcting for the edges in the full adjacency matrix).

For the graph analysis, neither global efficiency [F(1,42) = 0.70, p = 0.39], nor betweenness centrality of hippocampus [left: F(1,43) = 1.35, p = 0.25; right: F(1,43) = 0.45, p = 0.50] was significantly related to change in fitness.

None of the other measures of intrinsic brain activity (MVPA, BOLD_STD,) or ASL yielded significant change-change relationships to fitness.

In the post hoc analyses we observed several associations between aerobic fitness and intrinsic brain activity, both at baseline and over time. Longitudinal change in aerobic capacity predicted change in hippocampal connectivity. Whereas left hippocampus displayed a decreased connectivity to precentral gyrus with increasing fitness, the right hippocampus was associated with increased connectivity to frontal and parietal areas. The right hippocampus, similar to baseline, was associated with increased frontal and parietal connectivity. Previously, Voss et al. (2013) reported enhanced connectivity between parahippocampus and bilateral temporal cortex, as well as to occipital and parietal areas following a 12 months aerobic (walking) intervention. Similar to the current findings, they observed reduced connectivity of left hippocampus, although to the right prefrontal cortex rather than to the right precentral gyrus as seen in here. Interestingly, the authors showed that increased connectivity between bilateral parahippocampal gyrus and middle temporal gyri were positively predicted by increases in neurotrophic growth factors, although no direct test of the behavioral significance of any of the connectivity changes was presented.

We detected relationships between BOLD_STD and aerobic capacity at baseline, whereas gains in aerobic capacity did not relate to changes in BOLD_STD over time. Contrary to what previously has been reported (Burzynska et al., 2015; Gauthier et al., 2015), the relationships between aerobic fitness and BOLD_STD were negative and primarily observed in mid-temporal rather than cortical areas. Burzynska et al. (2015) reported cross-sectional associations between BOLD_STD and physical activity, suggesting that BOLD_STD could reflect a long-term cardiovascular trait, which coheres with the absence of any intervention related changes in current study. Notably, they did not observe relationship between aerobic capacity and BOLD_STD, but only a positive relationship between BOLD_STD and physical activity. Our findings suggest that high aerobic capacity is associated with a more stable BOLD-signal within the temporal medial lobe. Speculative, the high fit subjects would experience relatively lower levels of physical demand, which could reconcile current findings with those reported by Burzynska et al. (2015). However, no on-line monitoring of degree of physical work during scanning was performed. Thus, further studies are required to establish the relationship between gains in aerobic capacity and BOLD_STD.

Longitudinal changes in intrinsic brain activity that were predicted by change in aerobic fitness did not relate to changes in any of the investigated cognitive domains (EM, EF, and CS). A plausible reason for this is the absent relationships between the cognitive variables and VO₂-peak, both at baseline and in change-change scores (Supplementary Table S2). Unfortunately, this limits the conclusions of the functional meaning of the observed fitness related brain changes. Thus, any functional interpretations of these brain changes are unknown and would have to rely on reverse inference, for which we lack unique cognition-connectivity associations.

Taken together, the findings above support the notion that aerobic fitness primarily affects mid temporal brain regions. In addition to the exercise intervention studies investigating changes in intrinsic brain activity (Voss et al., 2013; Burzynska et al., 2015; Maass et al., 2015; Bär et al., 2016; Tozzi et al., 2016), corroborating evidence of hippocampal involvement is also provided by studies investigating structural changes (e.g., Erickson et al., 2011, for a review, see Erickson et al., 2014; Jonasson et al., 2017). The fact that hippocampus display a high capacity for plastic change in relation to environmental factors (for a review, see Leuner and Gould, 2010), and also undergoes age related resting state connectivity changes (Salami et al., 2016), provides further motivation to specifically probe the intrinsic brain activity in the mid-temporal lobe. By standardizing research protocols and outcome measures, the mechanisms linking hippocampal brain activity, brain structure, vascularity, molecular growth factors and behavior will hopefully be revealed (Duzel et al., 2016; Stillman et al., 2016).

Network based statistics analysis revealed a graph sub-component that was predicted by the longitudinal change in VO₂-peak. Whereas previous cross-sectional reports (Voss et al., 2016) showed a higher within network integrity (particularly in DMN) with higher fitness, our longitudinal findings suggest that increased fitness also modulate between network connectivity. Whether the explanations for the discrepant finding are found in differences in methodological approaches (e.g., graph defining set of nodes) or lack of statistical power remains to be determined.

Resting state networks derived by ICA revealed relationships between fitness and connectivity of both the S1M1 and the DMN. Decreased connectivity between M1S1 and right thalamus was associated with increased aerobic capacity. M1S1 receives input from parts of thalamus, and the decreased connectivity with gains in fitness seems counterintuitive at first. Speculatively, this could reflect enhanced neuronal efficiency, although such a hypothesis should be investigated using methods that complement the correlational fMRI approaches used here. The connectivity between the same S1M1 network and occipital cortex were negatively related to fitness at baseline. In a small exercise intervention among obese children, Krafft et al. (2014) observed that the ICA derived motor network underwent decreased connectivity to visual areas in precuneus in the exercise group, resembling current baseline finding. Although our fitness predicted S1M1 connectivity targeted different regions at baseline compared to over time, it seems that functional segregation of S1M1 reflects better fitness. However, in an aged cohort more similar to ours, Voss et al. (2010) did not see group by time changes in motor connectivity, why this finding would need to be investigated further.

Interestingly, we observed enhanced connectivity between the DMN component and the left orbitofrontal cortex that correlated with gains in fitness. Even though the orbital cluster is localized more laterally than the typical ventromedial frontal DMN hub, the observation bears similarities to previous literature of fitness related connectivity of the posterior anterior midline core of the DMN (Voss et al., 2010, 2016). Both these studies interpret the fitness-related DMN integrity as a rejuvenated connectivity pattern, typical for younger subjects. A similar interpretation could be given here, but the lack of associations between cognition and these fitness-related connectivity patters prevent any firm conclusion of their behavioral significance.