Waves of Maturation and Senescence in Micro-structural MRI Markers of Human Cortical Myelination over the Lifespan

1The cortex contains myelinated axons critical for efficient neuronal communication. More than a century ago, Paul Flechsig (1847–1929) noted that cortical myelination did not occur simultaneously in all regions, but rather in a succession of waves peaking at different ages. Motor and sensory cortical regions myelinated earlier than association cortical regions. extended this notion of waves of myelinogenesis, and noted how intracortical fibers of the association cortex appeared to have a long myelination cycle, demonstrating an immature “greyish pallor” beyond the third decade of life. Convergently, studies in rodents have shown myelin-forming oligodendrocyte cells are produced in waves. Specifically, three waves of oligodendroglial proliferation have been described, peaking at different stages of forebrain development before and after birth. Later in life, primate histology studies have reported cortical myelin degeneration, and aberrant remyelination.

2Magnetic resonance imaging (MRI) has allowed for the indirect study of myelin using volume and other measures likely related to myelin in development, adulthood, and aging. Such studies, focusing mainly on the white matter (WM) of the cerebrum, seem to largely support Flechsig’s idea of myelination heterochronicity. studied myelination during the first postnatal year using T1- and T2-weighted images, delineating five stages, or waves, of myelination. However, there have been no attempts to test Flechsig’s notion of cortical waves during childhood and adolescence, and whether similar patterns could be found late in life.

3There have been only a few studies encompassing a sufficient age range to investigate putative myelinogenesis, from childhood to old age, comprehensively across the human cortex. Studying the WM, fitted linear models to diffusion tensor imaging data before and after the age at peak fractional anisotropy (a marker of WM tract organization) in participants aged 11–90 years. In nine WM tracts, they reported that tracts with a faster rate of maturation had a faster rate of age-related decline. fitted a piecewise linear model of R1 in 24 WM pathways for participants aged 7–85 years. R1, a measure of the longitudinal relaxation rate of protons excited in the imaging process, is a quantitative marker sensitive to myelin. They found no systematic relationship between the age at which a tract matures and the age at which it begins to decline. studied six WM tracts in 61 adults aged 18–84 using myelin water volume fraction (a quantitative measure sensitive to myelin), and reported quadratic age effects. Importantly, a few studies have looked at cortical myelination. found logarithmic increases in myelin water volume fraction with age in children aged 1–6 years. Using methods sensitive to myelin based on magnetization transfer, showed a linear cortical myelin increase for participants aged 14–24 years. reported stable levels of myelin until around 60 years of age, followed by decrease, in several cortical regions including the Heschl’s gyri, for participants aged 19–75 years.

4The lifespan studies yielded mixed support for the concept that decline occurs in reverse order to growth of myelination. This so-called “last in, first out” hypothesis suggests that regions that myelinate last in the course of development, are the first to demyelinate later in life, i.e., are more affected by age. This heterochronic pattern has been related to association cortices compared with sensory cortices, possibly due to their underlying cytoarchitecture. The association cortex is particularly well connected, specialized for multimodal integration, and underpins most complex behaviors. Therefore, these regions might need to remain plastic for as long as possible to refine their connections and consolidate them by a prolonged phase of myelination. In line with this reasoning, but in contrast to the “last in, first out” hypothesis, these regions could be the last to decline (“last in, last out”). The high connectivity of the association cortex, and the importance of myelin in brain network functioning, also raise the question of how patterns of myelination relate to how the brain is organized as a network. Particularly, whether the most well-connected “hub” regions are later to mature and to decline.

5In this context, we aimed to investigate myelination across the cortex in humans (N = 484, aged 8–85 years), using a micro-structural MRI marker putatively linked to myelin: the ratio of T1- and T2-weighted MRI images (T1w/T2w). Although not a truly quantitative measure of intracortical myelin, and likely influenced by other characteristics within each voxel, variation in T1w/T2w has been shown to match histologically-derived myelin content. Studies comparing cortical T1 and T1w intensity with myelin stains in marmosets, have supported the interpretation of the T1w/T2w ratio as a reasonable estimate of relative intracortical myelin content. Thus, while acknowledging that the T1w/T2w ratio is far from a pure measure of myelin content, we adopt this interpretation here, as in our previous work. T1w/T2w variation has been related to individual differences in electrophysiology, and cognition, and has shown positive relationships with age in development and early adulthood, and negative relationships with age in later life. Here, we move further by testing Flechsig’s theory of waves across the lifespan, and by testing for associations between development and aging. From the T1w/T2w maps, we fitted smooth age-trajectories for 360 cortical regions in each of seven cytoarchitectonic classes. The first derivatives of these cross-sectionally estimated cortical growth curves allowed us to define three key micro-structural milestones, namely age at peak growth, age at onset of stability, and age at onset of decline. We analyzed the distribution of these milestones for deviations from unimodality (i.e., at least bimodality). We related the timing of the milestones both to (i) the estimated maximum rates of maturation and decline in each region; and (ii) the topological properties of each node in the structural covariance network estimated on average across a subset of participants, as well as a functional connectivity network based on functional MRI scans.

6On this basis, we tested three hypotheses: (i) the timing of micro-structural milestones – age at peak growth, age at onset of stability, and age at onset of decline - has at least a bimodal distribution indicative of successive waves of change in myelination, (ii) timing of these milestones is related to rates of micro-structural maturation and decline, and (iii) the variation in milestones and rates of micro-structural growth and decline relate to cytoarchitectonic properties, and network topology, namely connectivity strength or hubness of regions in the anatomical network. We particularly expected the association cortex to develop and decline differently compared with the most heavily myelinated primary motor and sensory cortices.

7Participants were included after screening for conditions assumed to affect CNS function (for details regarding Methods, see Supplementary Material). Based on the existence of a T1w and a T2w scan for each individual, and after stringent scan quality control (see below), we included 484 participants: 263 females (54.3%), mean age (SD) = 38.3 (22.5) years, median age = 34.6 years, age range = 8.2–85.4 years. Supplementary Figure S1 shows the gender and mean full-scale intelligence quotient for females and males across age deciles.

8Scans were acquired using a 12-channel head coil on a 1.5-T Siemens Avanto scanner (Siemens Medical Solutions, Erlangen, Germany) at Oslo University Hospital Rikshospitalet, in one session per participant, using a 3D T1w magnetization-prepared rapid gradient echo (T1w), and a 3D T2w sampling perfection with application optimized contrasts using different flip angle evolutions sequence (T2w). All raw scans were visually inspected to assess quality and identify motion-related artefacts. A total of 18 participants were excluded due to incomplete T2w scans, overfolding, or movement artefacts in T1w, T2w, or both. Of these, six scans were excluded due to signs of motion (all males, four participants aged between 9 and 12 years, and two participants aged 52 and 68 years). In evaluating this number of motion-compromised scans, it is relevant to keep in mind that the scans were also evaluated during scanning. If a scan was deemed to be compromised, the scan was repeated if possible. In addition, participants who experienced problems lying in the scanner, and were more likely to move, were less likely to complete more than the first (T1w) scan, thus not being eligible for the study due to lack of the second (T2w) scan.

9Still, to further assess potential confounding effects of image quality due to motion, we repeated our main analyses in two secondary analyses after excluding the 49 participants that scored below the 10th percentile on two measures of image quality putatively related to head motion: (i) the temporal signal-to-noise ratio in diffusion-weighted imaging data acquired in the same scanning session; and (ii) and the Euler number estimated in the T1w images (see Supplementary Material for details).

10T1w/T2w maps were created using the Human Connectome Project (HCP) processing pipeline (https://github.com/Washington-University/Pipelines), including processing with the Freesurfer 5.3 suite (http://surfer.nmr.mgh.harvard.edu). The T1w volume was divided on the aligned (using bbregister and spline interpolation) T2w volume, creating a T1w/T2w ratio volume. To estimate regional T1w/T2w, we used a multi-modal parcellation to divided each cerebral hemisphere into 180 areas. We sampled T1w/T2w values from the WM/gray matter (GM) boundary, at nine intermediate cortical depths (10% intervals), and at the GM/CSF boundary (Supplementary Fig. S2). Here, we focused on T1w/T2w measures at 70% depth from the pial surface to avoid partial volume effects (Supplementary Fig. S3).

11To fit age-trajectories without assuming a specific shape of the relationship a priori, we used penalized splines. Eight piecewise cubic B-spline basis functions were used to fit a smooth, non-linear curve from the weighted sum of these functions (Supplementary Fig. S4). We calculated the derivative, i.e., the slope at each point along the curve, as the differences between the predicted curve (using the predict function in R (www.r-project.org/)) from (i) one set of age values, and (ii) another set with slightly increased age values, divided by this increment:. To find significant increases or decreases, we estimated the 99% CI of the derivative. The CIs were estimated using predict.gam in R (https://cran.r-project.org), with the Bayesian posterior covariance matrix (see Supplementary Material for further details). Based on the CI, we extracted three curve features, or milestones, namely age at peak growth, age at onset of stability, and age at onset of decline, for each region. We also estimated cross-sectional measures of the rate of growth, and the rate of decline (Supplementary Fig. S5). As the CIs of the derivative curve are critical for the milestone estimates, we also calculated the CI using a bootstrap method to evaluate consistency across methods. Specifically, we drew 2000 random samples with replacements across participants, and calculated the growth curve and its derivative as above for each region. The 99% CIs were taken as the 0.05 and 0.995 percentile, respectively, of the resulting distribution of samples, at each point along the curve and the derivative, respectively. The percentile method was chosen for its simplicity compared with the original approach, although more complex, but potentially more accurate, alternatives have been proposed.

12We explored the relationship between cortical histology and the key milestones. Based on previous work, we assigned each of the 360 regions to one of five cytoarchitectonic types classified according to the scheme of von Economo and Koskinas, reflecting the five structural types of isocortex, or to two additional subtypes: limbic cortex, and the insular cortex.

13To investigate how the milestones, which are local in nature, might relate to global network properties of the brain, we constructed an anatomical network using structural covariance analysis. This approach allowed us to make a direct link, by using the same structural MRI data, between the milestone estimation and the network analysis. It has been shown that regions with high structural covariance or correlation are often involved in the same cognitive function and connected via WM pathways. Specifically, we correlated T1w/T2w for each region with all other regions, across participants between the milestones of maturity onset (34 years) and decline onset (72 years) in the global curve (see Fig. 1B), i.e., during the mature period with relatively stable T1w/T2w levels. This procedure yielded a 360 × 360 connectivity matrix, which was binarized employing a minimum spanning tree approach followed by global thresholding, retaining 10% of the strongest connections or edges (Supplementary Fig. S6A; ). From this model, we assessed degree and modularity, two of the most common and, for degree, interpretable network metrics. These analyses were carried out in Matlab (https://www.mathworks.com) using the Brain Connectivity Toolbox. The community structure was obtained using the Louvain algorithm, and consensus clustering. We empirically defined the resolution parameter by finding a local minimum for nodal versatility of modular affiliation (see Supplementary material for details). A network representation of the structural connectivity (Fig. 4ii) was visualized using NetworkX, version 2.1 (https://networkx.github.io).

14To further assess the relationship between the milestones, and brain network topology, we also created a functional network based on functional MRI scans at rest for a subset of the participants (aged between 20 and 72 years). The aim was to provide complementary information as structural and functional networks are robustly, but complexly related (see for instance ). Details of the sub-sample and the pre-processing can be found in the Supplementary material. From the resulting pre-processed fMRI volumes, regional average time series were extracted for the 360 regions. The time series were then fed to a wavelet decomposition to four wavelet scales using the maximal overlap discrete wavelet transform with a Daubechies 4 wavelet, performed by use of the brainwaver package in R (www.r-project.org/). We chose wavelet scale one-frequency band (0.08–0.17 Hz) for further investigation, but also repeated the analyses at the scale two-frequency band (0.04–0.08 Hz) to assess consistency. The wavelet correlation between every pair of regions yielded one 360 × 360 matrix of regional functional connectivity per participant. We then averaged these matrices across participants resulting in one 360 × 360 matrix of regional functional connectivity. Similar network construction and analyses were performed as for the structural network, with focus on the total degree analysis. Thus, we constructed a binary graph by retaining the 10% strongest connections for this average functional connectivity network (Supplementary Material, Fig. S7A), using the minimum spanning tree approach. This network representation was fed to network analyses of community structure, and total degree centrality.

15We tested for unimodality using Hartigans’ dip test statistic. We used an expectation–maximization (EM) algorithm to fit Gaussian finite mixture models with one and two components, where a maximum of two components would provide proof of concept, and not be liable to over-fitting, and tested for best fit using a bootstrap likelihood ratio test with 10 000 bootstraps. To facilitate the fitting process, nine values above 25 years of age, and four values around 40 years were removed from M1 and M3, respectively, before estimation and plotting. Relationships between the regional measures were tested using Spearman’s rank correlation, except between rate of maturation and decline, where a linear regression was performed. To assess differences in rates of maturation and decline, respectively, between early and late waves across the three milestones, we calculated 99% CIs. In the cases where the CI overlapped, we tested for differences between waves by Wilcoxon rank sum test. The same approach was used to test for differences in the degree measures between waves. To test for the differences across cytoarchitectonic classes, we employed Kruskal–Wallis tests by rank, and accounted for the multiple tests performed by applying false discovery rate correction.

16The average T1w/T2w ratio map across all participants shows that highly myelinated primary motor and sensory areas have higher T1w/T2w values compared with association areas (Fig. 1A). The global trajectory of T1w/T2w, on average over all cortical regions in each participant, is approximately an inverted-U function of age over nearly 80 years of the lifespan (Fig. 1B). There is a maturational increase or growth of T1w/T2w in the first part of life, before the curve gradually levels off and a plateau is reached in middle age, followed by a later period of decline from around 60–65 years.

17To quantify this age-related micro-structural variation more precisely, we estimated the slope of the T1w/T2w curve at each point in time, t, i.e., the first derivative of the trajectory, d(T1w/T2w)/dt (Fig. 1C). We then estimated 99% confidence intervals (CIs) for both the curve and its first derivative, which we used to rigorously define three key milestones (M1, M2, M3). M1, peak growth age, was defined as the age of maximum positive rate of change in T1w/T2w; M2, onset stability age, was defined as the earliest age after M1 at which the 99% interval for the derivative curve included zero, indicating that there was no significant growth or decline in T1w/T2w; and M3, onset decline age, was defined as the earliest age after M2 at which the 99% interval for the derivative curve did not include zero, indicating that there was significant negative rate of change in T1w/T2w (Fig. 1D). For the whole brain trajectory, peak growth age (M1) occurred at 19 years, onset stability age (M2) at 34 years, and onset decline age (M3) at 72 years. From the derivative curves we also estimated two additional metrics of interest (Supplementary Fig. S5): the maximum rate of growth in T1w/T2w (i.e., the derivative at M1), and the maximum rate of decline (i.e., the derivative at mid-point between M3 and the end of the age range).

18We estimated all three milestones for each of 360 cortical regions and plotted the age at each milestone for each region as cortical surface maps. Peak growth age (Fig. 2Ai, and Supplementary Material, Fig. S8) occurred earlier, i.e., before 13 years of age, predominantly in areas with high myelin content such as the primary motor and sensory areas, as well as the retrosplenial cortex and right LO2 (please see for explanation of the nomenclature), and area MT. Regions with a later peak growth age, i.e., after 13 years of age, were found in prefrontal, parietal, and temporal cortices, particularly in the right hemisphere, as well as medial occipital cortex.

19Onset of stability was reached between 29 and 50 years (Fig. 2Aii) and was found to occur later particularly in lateral temporal, medial prefrontal and anterior cingulate cortex bilaterally, as well as in left lateral frontal cortex, with regions in general reaching a plateau later in the left hemisphere.

20Onset of significant decline started between 58 and 77 years of age (Fig. 2Aiii) in 170 regions, while 190 regions, particularly in the left hemisphere, did not show a significant decline. In both hemispheres, onset of decline was observed earlier in motor and sensory cortices, and the retrosplenial complex.

21Inspection of the distributions of each of the milestones provides intuitive support for the existence of waves of change in T1w/T2w. For each milestone, there are two distinct peaks (Fig. 2B), suggesting that the distribution is not unimodal but at least bimodal (see also quantile-quantile plots, Supplementary Material, Fig. S9). We tested this observation formally, finding evidence for significant non-unimodality in all milestone distributions (Hartigan’s dip test, P < 0.001). A two-component mixture of Gaussian distributions was fitted to each of the milestone distributions, providing a significantly better fit than a unimodal model in all cases (bootstrap likelihood ratio test, all P < 0.001). Each milestone distribution could thus be partitioned into an early wave and a late wave (Fig. 2B), and each cortical region could be assigned to the early or late waves of each milestone (Fig. 2C).

22The first wave of peak growth occurred at a mean age of 9.4 years (95% CI [9.0, 9.7]) in relatively few areas of primary cortex, specialized for motor and somatosensory function; the second wave occurred at a mean age of 19.5 years (95% CI [19.2–19.7]) in a larger number of association areas, specialized for semantic and social functions, such as theory of mind (determined through a search of neurosynth.org, Fig. 2D).

23The first wave of onset of stability occurred at a mean age of 33.3 years (95% CI [33.0, 33.6]), and involved areas of parietal cortex specialized for attention. A smaller number of areas, mostly located in the left hemisphere and specialized for semantic language functions, reached stability later with mean age of onset 45.1 years (95% CI [44.4, 45.8]).

24The first wave of onset of decline occurred at a mean age of 68.8 years (95% CI [67.8, 69.8]), and involved mainly right hemisphere frontal and parietal regions specialized for attention and executive functions. The second wave in this case represented a large number of regions, mainly in the left hemisphere and specialized for semantic language functions, which did not show evidence of significant decline within the age range of the sample (up to 85.4 years).

25Surface maps of the regional variation in rate of maturation, and the rate of decline, are shown in Figure 3A and Supplementary Figure S10. We observed higher rates of maturation and decline in highly myelinated areas of primary motor and sensory areas and lower rates of maturation and decline in lightly myelinated frontal and temporal cortical areas. The average T1w/T2w map (Fig. 1A) and the rate of peak growth map (Fig. 3Ai) were strongly positively correlated (Spearman’s ρ = 0.75, P < 0.001). Conversely, the average T1w/T2w map was strongly negatively correlated with the rate of decline map (Fig. 3Aii; Spearman’s ρ = −0.67, P < 0.001). As expected given these results, the regions with high rates of maturation also had high rates of decline (Fig. 3B; Spearman’s ρ = −0.59, P < 0.0001), i.e., a “fastest in, fastest out” pattern.

26These maturation and decay rates of T1w/T2w, respectively, indicative of rates of cortical myelination and demyelination, were significantly different between waves as defined by the three milestones, M1–M3 (Fig. 3C). The rate of maturation was significantly more positive, and the rate of decline was significantly more negative, in the first M1 wave of regions that reached peak growth before puberty than in the second M1 wave of regions that reached peak growth after puberty (Supplementary Table S1A). Likewise, the first M2 wave of regions, which reached onset of stability earlier, had more positive rates of maturation than the second M2 wave of regions, which reached stability later. Finally, a similar pattern was evident for the M3 milestone: regions in the first wave, which reached onset of decline earlier, had more positive rates of maturation. Regions in the second M3 wave had less positive rates of maturation and less negative rates of decline (close to zero).

27The findings presented above for onset stability age (M2) and onset decline age (M3) point to a “last in, last out” pattern. That is, in general the regions which reach stability later (in the second wave of M2) are more likely to also show later onset of age-related decline (in the second wave of M3). This result is partly at odds with the “last in, first out” hypothesis, which is supported by studies using macro-structural measures, such as cortical thickness. Thus, we tested the relationships between age and cortical thickness in the early and late waves of onset decline age. To this end, similarly to previous analysis for the “last in, first out” hypothesis, we calculated the regression slope F value (weighted by direction of the slope as 11 slopes were not negative) of the cortical thickness-age relationships in each region, and interpreted more negative age-relationships as indicative of “first out”. That is, the “first out” manifests as steeper age-related decline.

28The sample used was formed of participants aged 34.5 (onset stability age for the average curve, Fig. 1C) or older (n = 242). Then, we tested whether the weighted F values were different for the two waves of onset of decline age. As hypothesized based on, regions in the second wave showed a stronger decline with age compared with regions which reached maturity earlier (Wilcoxon rank sum test W = 21194, P < 0.0001). This relationship was not altered by taking the T1w/T2w level of the region into account. This pattern is accordance with the “last in, first out” hypothesis. That is, the regions reaching onset of stability age late (“last in”), showed stronger decline in cortical thickness with higher age (“first out”).

29We tested the relationship between cortical histology and milestones of the T1w/T2w trajectory by classifying each of the 360 regions according to its membership of one of seven cytoarchitectonic classes defined by the von Economo and Koskinas atlas (Supplementary Fig. S11A). Different cytoarchitectonic classes and cerebral hemispheres were differently represented in the two waves at the M1 milestone (Supplementary Fig. S11B). The proportion of regions reaching peak growth before puberty was significantly different between cytoarchitectonic classes, with a higher proportion of primary motor and primary sensory cortices in the first wave than in the second wave (χ2 = 65, df = 6, PFDR = 1 × 10−10). The differences in proportion of cytoarchitectonic classes in the first and second waves of the M2 milestone were not significant after multiple comparison correction (χ2 = 14, df = 6, PFDR = 0.126, Puncorrected = 0.035). However, there were significant differences between the first and second waves of the M3 milestone, with a higher proportion of association cortical areas, especially in the left hemisphere, represented in the second wave (χ2 = 37, df = 6, PFDR = 1.3 × 10−5).

30To investigate how the local cortical milestones of age-related differences in T1w/T2w might relate to the topology of the cortical connectome, we estimated the pair-wise correlation of T1w/T2w between each possible pair of regions across participants between the ages of stability and decline onset in the global curve (see Fig. 1B; 37 and 69 years, respectively), that is, a period with relative stable T1w/T2w levels. From this micro-structural network model (Supplementary Fig. S6A), we constructed a binary graph by retaining the 10% strongest connections which could be decomposed into a community structure of modules (Supplementary Fig. S6B). For each regional node, we then estimated its degree centrality, i.e., the number of edges connecting it to the rest of the network; its intra-modular degree, i.e., the number of edges connecting it with nodes in the same module; and its inter-modular degree, i.e., the number of edges connecting it with nodes in other modules. The network had a fat-tailed degree distribution (Supplementary Fig. S6C), indicative of so-called hub regions. These high-degree regional nodes were mainly found in frontal, parietal, and temporal association cortices (Fig. 4Ai–ii).

31There were significant differences in nodal topology between the two waves of the M1 milestone (Fig. 4B, and Supplementary Fig. S12 and Table S1B). Regions which reached peak growth before puberty had significantly reduced total degree, intra-modular degree, and inter-modular degree, compared with regions which reached peak growth after puberty. Likewise, there were significant differences in nodal topology between the two waves of the M3 milestone. Regions in the second wave, which did not show onset of significant decline, had significantly higher total degree, and inter-modular degree, than regions in the first wave which had onset of decline in the 60-s. There were no significant differences in nodal topology between regions in the first and second waves of the M2 milestone.

32The functional network could also be decomposed into a community structure of modules (Supplementary Fig. S7B), and had a fat-tailed degree distribution (Supplementary Fig. S7C). The high-degree regional nodes were mainly found in primary visual and motor cortices, as previously reported (Fig. 5Ai–ii).

33As for the structural network, there were significant differences in functional network nodal topology between the two waves of the milestones (Fig. 5B, and Supplementary Table S1B). However, as functional network hubs were mainly found in primary visual and motor cortices, the patterns were different, with early wave regions showing higher total degree. For M1, regions which reached peak growth before puberty had significantly higher total degree compared with regions which reached peak growth after puberty. Likewise, there were significant differences in nodal topology between the two waves of the M3 milestone: regions which had onset of decline in the 60 s had significantly higher total degree than regions which did not show significant age-related decline. In contrast to the structural network, there was also a significant difference in nodal topology between regions in the first and second waves of the M2 milestone. Again, regions in the first wave, which reached stability earlier, had significantly higher total degree than regions in the second wave which reached onset of stability later. Repeating the analysis using the wavelet scale two-frequency band (0.04–0.08 Hz) produced nearly identical results (Supplementary Table S1B).

34We repeated parts of our analysis by calculating the CIs based on the bootstrap method. As the CIs based on the bootstrap method were highly similar to the original CIs (Supplementary Fig. S14Ai–ii), the milestones based on the bootstrapped CIs were also highly similar (Supplementary Fig. S14Aiii–v).

35We also performed two additional analyses based on sub-samples of 435 participants with 49 participants excluded on the grounds that they were in the lowest 10% of the sample for image quality metrics. For the analysis where 49 participants were excluded based on the tSNR, the fitted spline curve was highly similar to the curve estimated on the full sample (Pearson’s r = 0.98, Supplementary Fig. S14Bi–ii). For the analysis where participants were excluded based on the Euler number, the fitted curve in the younger half of the sample was again very similar to the full sample (r = 0.95, Supplementary Fig. S14Ci–ii). For the latter part of the age range, there was less evidence of a decline in T1w/T2w ratio (for 24 regions compared with 190 regions) compared with the results of analysis on the full sample. This reduction in regions showing decline when data were quality controlled by the Euler number criterion also made the distribution for M3 less clearly non-unimodal than it was in the whole sample (Hartigan’s dip test, P = 0.02 versus P < 0.001, Supplementary Fig. S14Cv). For the other milestones, the histograms of both quality-controlled sub-samples showed a non-unimodal distribution similar to the distributions identified by analysis of the whole sample (Supplementary Fig. S14Biii–v, and S14Ciii–iv).

36Based on a cross-sectional sample covering nearly eight decades of the lifespan, we delineated the age-related trajectory of maturation and senescence of a micro-structural MRI marker linked to myelin. We found indications of maturation and decay occurring in stages, or waves, at three key milestones: age at peak growth, age at onset of stability, and age at onset of decline. Importantly, we also observed that regions in the first wave of developmental milestones had faster rates of maturation and faster rates of decline. This result contradicted the “last in, first out” hypothesis, according to which early maturing regions should be more protected from age-related decline than late maturing regions. We found the opposite. By linking the timing of milestones to the cytoarchitectonic atlas of von Economo and Koskinas, we found that early maturing regions were histologically classified mainly as heavily myelinated primary motor and sensory cortices, whereas later maturing regions were mainly association, limbic and insular cortices. We also found that the high-degree hubs of the anatomical connectome had significantly later age at peak growth and later age at onset of decline than less well-connected nodes in the structural covariance network.

37We observed a bimodal distribution of all milestones. In accordance with histological results based mostly on early childhood brains, we interpret this finding as suggestive of waves of myelinogenesis. This concept has not been studied within a lifespan perspective, and the observations could be conceptualized as “echoes” in later life of earlier waves of neuro- and glio-genesis. Another interesting hypothesis relates post-pubertal myelination to sexual maturity, and the changes in executive and social cognition associated with adolescence. This suggestion relates to the timing of peak myelination which occurred in one wave before and one wave after puberty. The post-pubertal wave regions were more concentrated in association, insular, and limbic cortices, and this accords with the suggested link between post-pubertal myelination and the onset of psychiatric diseases during adolescence and early adulthood.

38The present micro-structural results indicated a “last in, last out” pattern. This result is in line with a recent study in marmosets using a similar approach to the one taken here, reporting that GM volume in primary motor and sensory regions generally declined earlier than associations regions. In contrast, the macro-structural cortical thickness results were in line with the “last in, first out” hypothesis. MRI lifespan studies using GM macro-structural measures have found support for this hypothesis (e.g., ), while others using micro-structural measures of the cortex (e.g., ) have not. One possible explanation for this discrepancy between the lifespan trajectories for macro- versus micro-structural aspects of the cortex may be that loss of myelin in primary motor and sensory areas could offset cortical thinning. As both the cortex and the subjacent WM become “grayer” due to demyelination, the WM/GM boundary recedes, making thinning appear less extensive. This process is essentially an inversion of the suggested explanation of thinning in childhood and adolescence, where increased myelination makes the deeper layers of cortex appear “whiter”. The more negative macro-structural relationships with age observed in the lightly myelinated late areas might be driven more by processes other than demyelination, for instance reduction in dendritic spine density of pyramidal neurons. Another possible explanation is that demyelination with age is offset by other types of myelin degeneration, such as the formation of redundant myelin and increasing thickness of myelin sheaths, as reported in older aged monkeys.

39Although we did not find a clear relationship between age at onset of stability and age at onset of decline, we did find that regions showing higher rates of peak growth, and doing so relatively early, also showed an early onset of decline and a faster rate of decline. The mainly primary motor and sensory regions showed earlier and faster rates of cortical myelination and earlier and faster rates of cortical demyelination with high age. Thus, in addition to a “last in, last out” pattern, our results fit partly with a “gain-predicts-loss” hypothesis. This hypothesis predicts that the rate of change in development will be mirrored by the rate of decline in aging. We found that the aging process is not exactly the reverse of development, in contrast to R1 trajectories of WM tracts. Instead, the T1w/T2w trajectories increased more steeply in development compared with less rapid rates of decline later in life, as in the diffusivity trajectories of the same tracts. Thus, there is a similarity between the models of micro-structural development over the lifespan in WM and cortex. Whether these effects are governed by similar or different mechanisms warrants closer scrutiny in future studies. Interestingly, the differences in rate of development and decline in primary motor and sensory regions compared with association, limbic and insular regions, fits well with the processing-speed theory of adult age differences in cognition. This theory proposes that processing speed reduction leads to the age-related decline in higher-order cognition, such as memory and reasoning. In development, a similar idea, but in an opposite direction with increased speed, has been proposed. Accordingly in this context, processing speed relates more strongly to sensorimotor areas than association cortices, while association cortices relate more strongly to higher-order cognitive functioning.

40Another feature of our findings was the relationship between waves of putative (de-)myelination and high-degree network regions of the anatomical connectome, so-called hubs. We found that structural hubs had significantly later age at peak growth and later age at onset of decline than less well-connected nodes in the structural covariance network. observed that association regions were the focus of a relatively late wave of myelinogenesis. In our data and previous studies, structural hub regions are disproportionately located in association cortex. Thus, the observation of different micro-structural trajectories for hubs compared with primary motor and sensory areas aligns with Flechsig’s idea. Here, we extend this notion by showing that structural hub regions seem to be the focus of a relatively late wave of myelinogenesis also later in life. Although the relationships were not strong, the result followed from our earlier findings that association cortices (and limbic and insular cortices) show later and lower rates of micro-structural maturation and decline compared with motor and sensory cortex. As discussed in, cortical myelination in development might be topologically targeted to optimize structural hub performance. Similar optimization processes might be at play in aging, which might break down in disorders such as Alzheimer’s disease.

41We also found relationships between waves of putative (de-)myelination and the hubs of the functional network. However, the pattern here was different, as the functional hubs were more posteriorly located, particularly in primary visual and motor cortices, which mainly comprised the first wave of earlier maturing regions at each milestone. Compared with structural network hubs, the interpretation of degree-based hubs in functional networks hubs is less straightforward, as functional connections reflect indirect communication occurring within the underlying structural network. These findings should thus be interpreted with caution.

42First, the results are based on cross-sectional analyses, and longitudinal data are needed to estimate measures of individual change. Cross-sectional studies may suffer from cohort effects, potentially exacerbated by differences in recruitment bias between age groups. However, few longitudinal studies cover the lifespan, impeding measurements of change rates as a function of age, and such studies may also be influenced by selective attrition. Longitudinal lifespan studies, measuring change over relatively few years (<10 years) for each participant, may also suffer cohort effects. Hence, given the relative paucity of studies including most of the lifespan (i.e., from children to older adults), both cross-sectionally and longitudinally, we believe the current approach has some merits.

43Second, the T1w/T2w ratio measure used here stems from T1- and T2-weighted scans, at a conventional resolution. Using a quantitative measure at high resolution likely would have been preferable for sensitivity and specificity. Moreover, the link between the T1w/T2w ratio and myelin content has recently been questioned. Indeed, it seems unlikely that the T1w/T2w ratio, or other MRI indices related to myelin such as R1, should depend solely on myelin given the complex biophysical process underlying MRI relaxation rates. More specifically, reported positive correlations between the T1w/T2 ratio and myelin-associated genes, supporting its use as a myelin marker, but concluded that molecule size, transcription marking oligodendrocytes, axon caliber, and pH were more strongly associated with T1w/T2w., in 20 participants aged between 24 and 70 years, observed positive correlations between T1w/T2w and myelin water fraction within six WM tracts, but negative correlations between different tracts. The authors suggest that the T1w/T2w ratio relates more to variation in caliber and packing density of the axons than myelin content. In this regard, it is interesting to note that R1, in both ex vivo and in vivo imaging of six and four WM regions, respectively, correlated with histology-derived axon size, but not with myelin content in rats. caution against interpreting T1 differences between regions or participants as necessarily reflecting differences in myelin content, even when using potentially more promising T1 measurements than used here to quantitatively map myelin content. Accordingly, as previously noted, interpretation of signal intensity differences across the lifespan should be made cautiously without histology data, as a host of neuronal and non-neuronal processes likely play out and potentially interact at various intensities and stages across development and aging.

44Interestingly, combined T1w/T2w imaging and histology, in nine late stage secondary progressive multiple sclerosis patients aged 51–86 years. They correlated T1w/T2w ratio values in five normal-appearing cortical regions with histology-derived measures of myelin, axon and dendrites. The result showed that T1w/T2w ratio related significantly to dendrite, but not myelin (or axonal), density. Although, as acknowledged by the authors, the nature of the sample (late stage multiple sclerosis brain donors) might impede the generalizability to healthy participants, the finding is potentially important, and highlights the importance, and difficulty, of histological validation of in vivo MRI indices. Taken together, it might be more appropriate here, to use the term “apparent myelin”, similar to “apparent cortical thickness” and “apparent diffusion coefficient”, to highlight the aforementioned challenges in in vivo MRI.

45Third, many regions (particularly in the left hemisphere) did not exhibit age-related decline, which contrasts with other MRI markers such as cortical thickness. Hemispheric asymmetry in patterns of age-related cortical volume shrinkage in aging have been inconsistent. Still, verbal abilities are relatively well preserved in aging, and a potential relationship between intracortical micro-structure and preservation of verbal abilities, which are known to depend mainly on the left hemisphere, constitutes an interesting hypothesis for future research. Similarly, assessment of relationships between individual differences in the timing of these milestones and behavioral or neuropsychological changes associated with development would be an important task for future research.

46Finally, recent studies have shown that even subtle motion might bias anatomical estimates in studies of brain structure and age (see for instance ). Here, all images were thoroughly visually inspected for motion, both during and after scanning. Still, we further probed our data by excluding 10% of the sample with the lowest data quality defined by two criteria of image quality: tSNR from diffusion-weighted scans and the Euler number from the T1w scans. Repeating the main analyses on these sub-samples, the results in the younger half of the age range were virtually unaltered. This consistency of results was also apparent in the latter part of the age range when using tSNR as the quality index. However, when using the Euler number, the number of regions showing decline in older age was reduced to 24 regions compared with 170 regions in the full sample. However, it is unknown if this reduction in regions showing age-related decline in T1w/T2w ratio is due to head movement. For instance, it is possible that the Euler number could be increased in older adults with negligible motion during scanning because it is a measure of the topological complexity of the reconstructed cortical surface and age-related change in cortical structure could adversely impact on the performance of the cortical segmentation algorithm in FreeSurfer.

47The consistency of results in the younger half of the sampled age range, where head motion has been identified as pervasive (see for instance ), and in aging when using tSNR as quality index, supports the notion that our principal findings are not unduly influenced by movement during scanning. That said, these analyses were all post-hoc. Future lifespan studies of age-related structural brain effects should aim for direct and prospective estimates of head motion during scanning.

48As indicated by the present cross-sectional age differences, peak growth, stability, and decline of an MRI measure related to intracortical myelination proceeds in stages, or waves, in different regions. Primary motor and sensory regions reached all milestones earlier than association, limbic, and insular regions. The results contrasted with the so-called “last in, first out” hypothesis. Instead, the findings suggested a “last in, last out” pattern, where the regions less likely to show decline in aging were slowly maturing, association cortical hubs. In addition, the results fit partly with a gain-predicts-loss hypothesis, where regions which mature fastest, decline fastest. Future studies might benefit from looking more closely into how timing of waves relate to the onset of neuropsychiatric disorders, as well as the changes in cognitive and social skills, in development and aging.