1. Introduction
Mapping the human brain development and aging longitudinally from a unique dataset is barely impossible due to our extensive lifespan. Suppose one wants to understand the time evolution of the brain morphology through the whole lifespan. In this case, it is mandatory to combine multiple datasets and devise methods that allow a fair comparison among them. The recent advent of large heterogeneous datasets (to ensure the legit results across several populations) and data curation innovations allowed neuroscience to explore the brain’s changes on an enormous scale, from functional data to structural studies. Here, we propose a novel methodology of combining structural MRI from different acquisition sites and equipment, acknowledging their heterogeneity and providing insights about cortical folding during development and aging, with a practical application of these results in Alzheimer’s Disease, the most common dementia worldwide. This effort is essential to provide experimental information about the evolution span of the cortex morphology and build a cortical folding theory.
In the past years, the study of cortical folding in humans and other mammals was extensively promoted by the application of translational science, open-access databases, and data science tools which allowed the field to grow in multiple and sometimes convergent directions (Madan and Kensinger, 2016; Llinares-Benadero and Borrell, 2019; Mota et al., 2019; Essen, 2020). Due to its biological background, cortical folding measurements have been included in human brain structure analysis, discriminating disease from healthy controls (Cao et al., 2017; Lamballais et al., 2020), describing its correlation with cognition (Núñez et al., 2020), and relation to aging (Hogstrom et al., 2013; Madan, 2021). Despite the lack of consensus defining a unique theory that explains cortical folding on every scale. Mota and Herculano-Houzel (2015) proposed a cortical folding model respected by more than 50 mammals’ brain hemispheres. It predicts a power-law relationship between cortical thickness (T), exposed (AE), and total area (AT) (Equation 1)
where α is the self-similarity index, a universal constant, with a theoretical value of 1.25 and calculated as 1.305 (Mota and Herculano-Houzel, 2015) when considering a heterogeneous dataset of different species. It is an index of complexity, reflecting how much detail in a real fractal pattern changes with the scale at which it is measured. The constant k is dimensionless and associated with conserved viscoelastic properties of cortical matter. Following this equation, the usual geometrical variables used to describe the cortex (T, AT, and AE) are not independent in this model.
Wang et al. continued this study, validating the theory for different groups of humans (Wang et al., 2016) in specific regions of interest (Wang et al., 2019b), obtaining new evidence that the cortex is indeed a self-similar structure. Wang et al. (2021) showed that from this model one could derive a more natural set of nearly independent morphological variables, K, S, and I, that could be used to improve disease discrimination from regular and typical structural changes of the brain, such as aging.
In this framework, the morphology of each cortex is expressed as a point in a three-dimensional abstract space, with each coordinate component corresponding to the log-value of an independent morphometric variable (Wang et al., 2021). The use of log-values guarantees that linear combinations of the basis vectors correspond to power-law relations. Conversely, as long as any new variables are expressible as power-laws, different but equivalent sets of variables can, thus, be related to each other by a change of Cartesian coordinate bases. As a starting point, we use the log-values of the commonly morphometric variables, the total area log10AT, exposed area log10AE and average thickness . We then derive our new variables thusly: K = log10k (Equation 2) is a near-invariant quantity obtained by isolating k in Equation (1). S (Equation 3), also dimensionless, encapsulates the aspect of brain shape that changes more significantly across cortices, and I (Equation 4), represents brain isometric volume and carries the information about overall cortical size. Both are derived from planes perpendicular to K, noting that K remains unchanged when applying an isometric transformation on the brain (a direction in which all areas scale equally) resulting in I, and the perpendicular dimensionless vector S resultant from K × I . Those three variables combined may help distinguish pathological events similar to age effects, such as AD.
Baseline values (age-specific norm values) of those cortical morphological variables and their inherent uncertainties over the human lifespan were never defined due to methodological limitations of combining multicenter studies of brain structural images (Gronenschild et al., 2012; Fortin et al., 2018). However, both are crucial for applied studies. Uncertainties determine measuring limitations and allow a proper comparison between different cohorts. Baseline values allow clinical applications and, when necessary, simulations of healthy control data. Multi-center studies have become more common in the last few decades due to improved data sharing and the open science trend. However, combining MRI data from multiple experiments is not a trivial task. Part of these uncertainties is explained by random errors, a direct consequence of the natural variance within the species individuals, and the impossibility of replicating the exact same condition over multiple data acquisitions. Theoretically, these random errors are the same if one considers multiple experiments.
Confounding components could be added when gathering MRI images with differences in acquisition parameters, acquisition equipment, and versions of the post-processing software (Dickerson et al., 2008; Gronenschild et al., 2012). Together, they add uncertainty to the data due to a systematic effect, ultimately reflecting a systematic difference in morphological variables calculated from different databases. In practice, the natural random fluctuations are convoluted with the uncertainty due to a systematic effect, enhancing the data spread.
Recently, a vast literature on primary morphological variables from MR images (Schnack et al., 2010; Pomponio et al., 2020; Bethlehem et al., 2022), and more complex ones (Fortin et al., 2018), explored the limitations of compiling images acquired with different scanners and protocols and its implication in statistical analysis results. The suggested solution for these limitations is a procedure called data harmonization, in which there is an “explicit removal of site-related effects in multi-site data” (Pomponio et al., 2020).
It is essential to add that one cannot, in general, distinguish the natural variance in healthy human anatomy (σNatural) and data acquisition random errors (σRandom). Also, it is not possible to distinguish between the uncertainty due to a systematic effect of the acquisition and the processing. In other words, any measurement setup can include an additional value at its final measure that is repeated for all measurements. Here, the acquisition effects are due to differences in MR structural images acquisition protocols and their parameters, while processing errors are due to differences in the software version, pipeline parameters, and computer used for processing (Gronenschild et al., 2012). In Fortin et al. (2018) acquisition and processing effects are described as scanner and site effects. However, in the case of the methodology of this manuscript, it is impossible to quantify the contribution of this systematic shift, from acquisition and processing, since we did not compare processing or acquisition parameters by controlling the other error sources. For this, one would need a diligent procedure to track all differences between acquisition and processing as well as multiple image of the same individual for each sample. This manuscript deals only with the first case, approximating the random errors due to acquisition and processing for all samples using three subsequent structural images of the subjects available only for one particular sample (refer to Section Supplementary Material Section 2). The total uncertainty associated with the variables (K, S, and I) and their components is then summarized in Equation (5), the error propagation formula with no correlation since we are treating errors from uncorrelated sources by definition.
At which σX, is the total uncertainty for K, S, and I, σNatural represents the so-called here Natural Fluctuation, σRandom, the random variance that can be evaluated from the repeated measures with the same subject at each Sample acquisition and processing steps, and σSample, for the uncertainty due to a systematic effect on each sample.
In this manuscript, we developed a simplistic harmonization method that allows a multi-site analysis, determining K, S, and I baseline values, their uncertainties, and the ratio of changing rates through the years for cross-sectional images of healthy controls. Although this work focused is these novel morphological variables, we also provided the baseline of usual morphological variables such as cortical thickness (T) and Gyrification Index (GI). A robust estimation framework of a baseline time function is important to provide future clinical studies with enough information to compare brain structure trajectories in case of a pathological investigation. We then, as a clinical application and extension of de Moraes et al. (2022) and Wang et al. (2019b), verify if Alzheimer’s Disease (AD) is similar to a premature accelerated aging brain in terms of the independent morphological variables, K, S, and I, by comparing their rates. Finally, we have determined the time evolution of the self-similarity dimension α proposed in the Mota and Herculano-Houzel (2015) model, interpreting all results within the framework of this theory.
2. Materials and Methods
Participants and data included in this analysis were either acquired by the Instituto D’Or de Pesquisa e Ensino (IDOR), Rio de Janeiro, Brazil; used at de Moraes et al. (2022) (approved by the Hospital Copa D’Or Research Ethics Committee under protocol number CAAE 47163715.0.0000.5249); shread datasets included at Wang et al. (2016, 2019b), as the Alzheimer’s Disease Neuroimaging Initiative (ADNI) database (adni.loni.usc.edu), the Human Connectome Project (HCP), Information eXtraction from Images (IXI), Nathan Kline Institute (NKI) and Open Access Series of Imaging Studies (OASIS) Wang et al. (2016, 2019b); or from open-access databases, as Ultra-high field adult lifespan (AHEAD) (Alkemade et al., 2020) and (Snoek et al., 2021). The total number of subjects is 3,095 healthy controls from 4 to 96 years old and 555 Alzheimer’s Disease subjects from 56 to 92 years old. Datasets’ demographics are summarized in Table 1. AHEAD and HCPr900 only include age range instead of the actual age due to local Ethics Committee rules. To overcome the lack of an exact number, we assumed the interval’s mean age as a reference.
Table 1. Summary for each dataset.
To investigate the uncertainty due to a systematic effect during image acquisition and variances in repeated measures performed in the same conditions, we included the first 50 subjects from the AOMIC ID1000 (Snoek et al., 2021), described in Table 1.
The structural images AHEAD, AOMIC ID1000, PIOP1, and PIOP2 were processed in FreeSurfer v6.0.0 (Fischl, 2012) without manual intervention at the surfaces (McCarthy et al., 2015). The original project from IDOR images is longitudinal. Those images in particular were processed with the longitudinal pipeline. The FreeSurfer localGI pipeline generates the external surface and calculates each vertex’s local Gyrification Index (localGI) (Schaer et al., 2008). Values of Cortical Thickness, Total Area, Exposed Area, and Local Gyrification Index were extracted with Cortical Folding Analysis Tool (Wang et al., 2019a). We defined as ROI the whole hemisphere, frontal, temporal, occipital, and lateral lobes (based on FreeSurfer definition of lobes). The lobes’ area measurements were corrected by their integrated Gaussian Curvature, removing the partition size effect and directly comparing lobes and hemisphere cortical folding, as described in Wang et al. (2019b).
Datasets’ demographics, acquisition, and processing information are summarized in Supplementary Material.
2.1. Statistical Data Analysis
Multiple comparisons of means were made with ANOVA, correlations between cortical folding independent variables, and age estimated with Pearson r and post-hoc evaluations with Tukey multiple comparisons of means, including comparisons within Diagnostic groups. The statistical significance threshold was set at α = 0.05, and multiple corrections (Bonferroni) were applied when needed. The Healthy Control subjects are used as a reference in the necessary normalization and harmonization procedures. A standard linear regression was used to obtain the self-similarity index α.
All analyses were done considering the whole cortical hemisphere and the lobes (refer to Supplementary Material), as most diseases imply local or non-uniform global structural damage. Linear Mixed Models (LMM), standard linear regression, and statistics derived from their results were analyzed with RStudio (v1.4.1717 and R v4.1.1). The Bayesian model comparison was developed in Python3 with the package PyMC3. Harmonization and deaging codes are available at https://zenodo.org/record/5879895 (de Moraes et al., 2021).
2.2. Harmonization
The data harmonization procedure is based on the time evolution of the basic morphological variables Y = T, AT, or AE, modeled as an exponential decay with the age tage. Using Linear Mixed Models (LMM) in the log-linear scale, we are able to perform a joint linear fit for each ROI, maintaining the same angular coefficient a among all samples but different linear coefficient bj, ROI according to Equation (6). We then subtracted the sample shift, given by the linear coefficients bj, ROI), directly obtained from the Sample Random Effects (log10(Variable) ~Age×ROI+(1|Sample:ROI)) in Healthy Controls. From this procedure, we directly calculate the variables K, S, and I and determine for all morphological variables of interest (also T, AT, AE, GI = At/Ae) their baseline and uncertainties. This harmonization procedure was validated as the most probable description of the data using a Bayesian approach. This simple harmonization procedure is equivalent to the approach presented in Bethlehem et al. (2022) for the special case of a purely linear function.
2.2.1. Bayesian Model Selection
In order to show that the harmonization procedure used was the most suitable, one can compare it with other models to describe the data. For example, the angular coefficient could also be a random effect for the LMM, varying across the samples. Using Bayesian model selection (Gregory, 2005), one can show that the LMM used in this study is the most probable description of the data.
Considering our dataset comprised by J samples where D(j) = {d1(j), d2(j), ⋯ , dNj(j), } denotes the data of a morphological variable in the j-th sample. The degree of belief in a model M that describes a normally distributed data for each age is given by the posterior probability:
The term μ(xi(j), A(j)β) is given by the model M proposed to describe the data, with the set {A(j)β} being the free parameters of the model defined within a range ΔA(j)β. A prior probability P(A(j)β|I) is assigned for all free parameters reflecting the prior knowledge about their values. The basic model of this study assumes that the log of morphological variables follows μ = ax+b(j), with the dispersion σ(j) for each sample. Uninformative priors were used for the free parameters. This integral can be easily done numerically using Python3 with PyMC3 package or similar. The advantage of this approach is that it naturally penalizes the insertion of unnecessary new parameters avoiding over-fitting.
With multiple models Mi, one defines the odds as . As , the probability of the model M1 is given by:
The alternative models considered by this study were: (i) μ = ax+b(j), with the same dispersion σ(j); (ii) μ = a(j)x+b(j), with the dispersion σ(j) for each sample; (iii) μ = a(j)x+b(j), with the same dispersion σ(j) for each sample; (iv) we also tested the simplest non-linear trend , with the dispersion σ(j) for each sample; we have found that the LMM used in the study is the most probable. The code implementing this approach is given by the authors.
It is important to note that all parameters calculated using the LMM approach can be obtained in this framework by a slight modification of Equation (8). Instead of marginalizing all the parameters, by leaving one out, we obtain P(Aj1|D) which gives directly the most probable value for the parameter Aj1 and its credible region.
3. Results
3.1. Multisite Harmonization
The analysis starts with data harmonization, a crucial step to combining information from multiple sites allowing a fair comparison among datasets. Figure 1 shows the result of the harmonization procedure and its reflection on the primary morphological variables. The LMMs are used to fit the data assuming a single angular coefficient for all samples and different linear coefficients due to the uncertainty caused by systematic effect. We considered the variables’ trend with age linear based on the correlation between T, AT, and AE with Age for the Healthy Control group calculated from our data. There was a statistically significant correlation for T (Pearson r = −0.6949; R2 = 0.48; p < 0.001), AT (Pearson r = −0.4711; R2 = 0.22; p < 0.001), and AE (Pearson r = −40.2005; R2 = 0.42; p < 0.001). The harmonization was validated as the most probable description of this data using the Bayesian approach.
Figure 1. Basic morphological variables through age for Healthy Controls (A–C) raw data (C–F) after harmonizing, removing the estimated residual from the Linear Mixed Model. This procedure results in the harmonized values of the novel morphological variables K, S, and I.
3.2. Baseline, Rate Estimates, and Alzheimer’s Disease Diagnostic Discrimination
An immediate consequence of the harmonization is defining a baseline value and its uncertainty for any morphological variable considered in this work. The LMM provides quantitative information about the trends of these variables through age. The estimated uncertainties values are summarized in Table 2, and summaries of linear mixed models are available in Supplementary Material. The so-called Natural Fluctuation (σNatural), representing human diversity, was considered the Residual Standard Deviation of the models. The approximation for the random component (σRandom) of the experimental error derived from the multiple-image acquisition and processing obtained using only 50 subjects from the AOMIC ID1000 dataset is fully described in Supplementary Material. We considered the Random Effect of the Sample category (intercept) Standard Deviation as the uncertainty caused by a systematic effect component of the experimental error (σSample). The σX is the resultant uncertainty, composed of the previously described errors, and is the square root sum of each one’s squared values. Thereby, σK = 0.026 (4.9%), σS = 0.14 (1.6%), and σI = 0.091 (0.9%) with the corresponding fraction of their healthy controls intercept.
Table 2. Summary of the estimated uncertainty components for each morphological variable of interest.
Considering Alzheimer’s Disease, we used the baseline (Figure 2) to verify if the age trends of K, S, and I are different between the diseased and Healthy Control groups. We extracted the rates of change for each variable and diagnostic, summarized in Table 3 and displayed in Figure 3. The slopes were compared with post-hoc pairwise mean comparisons.
Figure 2. Cortical folding variables (A) K, (B) S, and (C) I through age for Healthy Controls after harmonizing in the primary morphological variables T, AT , and AE (removing the estimated residual from the Linear Mixed Model).
Table 3. Changing rate per year for each variable (after harmonization) and diagnostic.
Figure 3. Fitted values extracted from the linear regression model after data harmonization. Bars represent the 95% CI. Complete summary of linear models in Supplementary Material. (A) Concerning K, Alzheimer’s Disease (AD) has a shallow slope, meaning small changes with Age. Compared to healthy controls, AD values of K are almost constant and similar to older subjects (AD slope p < 0.0001 and CTL slope p < 0.0001; pairwise comparison estimate −0.00090, p < 0.001). (B) For S, the AD and CTL patterns have similar intercepts, but statistically different slopes (AD slope p = 0.004 and CTL slope p < 0.0001; pairwise comparison estimate 0.0013, p = 0.004). (C) For I, that reflects brain volume, CTL has a decreasing volume with aging, while AD has a smaller slope (AD slope p < 0.0001 and CTL slope p < 0.0001; pairwise comparison estimate −0.0017, p < 0.001). The curves for the AD group were extended to early ages to ease the comparison.
There is a significant difference between Control and AD slopes in K (pairwise comparison estimate −0.00090, p < 0.001), S (pairwise comparison estimate 0.0013, p = 0.004), and I (pairwise comparison estimate −0.0017, p < 0.001) as expected from Figure 3. For K, the AD curve is virtually flat as if it had reached a plateau, while for I, we see a constant reduction in brain isometric volume, with a reduced intercept. Comparing the mean values within AD and CTL, there is a significant difference for K (pairwise comparison estimate 0.042, p < 0.001, representing 7.9% of the CTL mean) value, S (pairwise comparison estimate −0.14, p < 0.001, 1.5%), and I (pairwise comparison estimate 0.11, p < 0.001, 1.0%).
In order to complete the full description of the cortical morphology within the framework established by the Mota and Herculano-Houzel’s model, we separated the data from both Health controls and patients with AD per decade and fitted the self-similarity dimension α theorized to be 1.25. The result is shown in Figure 4.
Figure 4. Slope α with 95% CI derived from the Cortical Folding model from Mota and Herculano-Houzel through age for Healthy Controls and AD subjects. Data points from (0, 5] and [95, 100] years old were omitted due to the reduced data. The numbers on top of each point indicate the number of subjects at each interval.
As most diseases imply local or non-uniform global structural damage, we expanded the results of trajectories to the lobes. We included ROI as a fixed effect in the linear mixed model equation with harmonized variables: Variable~Age×Diagnostic×ROI+(1|Sample:Diagnostic:ROI). Rates are summarized in Supplementary Material. We discriminated AD and healthy aging (CTL) K trends with pairwise comparisons at the Frontal (p adj <0.0001), Occipital (p adj <0.0001), Parietal (p adj <0.0001), and Temporal lobes (p adj <0.0001), as expected from de Moraes et al. (2022). Despite the rates being smaller for AD than CTL, the results mimic the whole brain with Alzheimer’s Disease trajectories appearing as a plateau during aging. The Frontal and Parietal lobe shows the least folded pattern (K) in AD (within lobe comparison are described in Supplementary Material). Also, there are significant differences between S pairwise comparisons in CTL and AD at the Frontal (p adj = 0.00035) and Parietal lobes (p adj <0.0001); and among I, between CTL and AD at the Frontal (p adj <0.0001), Occipital (p adj = 0.0031), Parietal (p adj <0.0001), and Temporal lobes (p adj <0.0001)).
4. Discussion
We analyzed a combination of datasets aiming to establish the K, S, and I baseline through the human lifespan. The suggested model included complexity at two levels: accounting for heterogeneous acquisition and processing and heterogeneous brain structures across healthy and disease subjects. The hypothesis is that non-homogeneous methodology (inclusion criteria, acquisition, and processing) implies more significant uncertainty due to a systematic effect. Thereby, we suggested a harmonization procedure that allows the comparison of data from multisite reducing this kind of experimental error. This result is innovative in its simplicity and focuses on expanding our independent morphological variables knowledge. The main caveat of our procedure is the impossibility of accounting for all differences between the samples and using this to define a universal harmonization according to each acquisition protocol. Our method does not depend on one reference sample or a gold standard definition. The removal of the shift for each sample is carried out based on the available healthy control data, which can always be added to the analysis.
In other words, our method does not account for a possible uncertainty due to a global systematic effect that could shift all samples equally. By looking at the uncertainty in cortical thickness estimation, our results are compatible with previous work. Our findings of 0.02 mm (Table 2) of uncertainty in cortical thickness are in agreement with a variability less than 0.03 mm found in a test-retest analysis by Han et al. (2006). Moreover, they suggest that uncertainties of 0.15 mm across platforms of manufacturers and 0.17 mm across field strength, at this study, coupled with processing uncertainties and estimated in 0.14 mm. A most recent study by Frangou et al. (2020) suggests a Mean Inter-individual variation of 0.07 ± 0.06 mm for hemisphere’s cortical thickness at all age ranges in the study compared to 0.10 mm from our findings. This compatibility suggests that the possible global systematic effect is close to 0; thus, the harmonization could be regarded as being close to the universal one.
The typical values for the K, S, and I independent variables through the whole life span and their related uncertainties were obtained by this study. These curves not only can be used to discriminate pathological from healthy aging in cross-sectional harmonized data but also are a powerful tool to morphologically describe non-typical cortex. Estimating the uncertainties is crucial to constructing the baseline curves and was carefully developed by this study. We estimated three types of uncertainties with linear mixed models and repeated measures: the uncertainty due to a systematic effect (combination of multiple datasets), repeated measures, and the natural variation in the human species. The error from acquisition and processing were coupled at the individual and sample levels. We expect future studies to disassociate the uncertainties in detailed factors (test-retest analysis, across platforms, field strengths, and acquisition parameters with the same sample of subjects) to measure morphological variables of interest better.
We have demonstrated the usage of the baseline curves by analyzing AD data. In terms of K, a measure of axonal tension, the AD curve is virtually flat as if it had reached a plateau of reduced brain tension. Compared with Healthy Controls, this plateau seems to be the same experienced by the oldest subjects. Interpreting in terms of the Mota and Herculano-Houzel cortex model, AD can be regarded as a form of accelerated aging from the axonal tension perspective. On the other hand, the evolution of S, I, and even the self-similarity dimension α suggest that the AD cortex is not morphologically similar to the health one, with a very different shape, isometric volume, and time evolution. However, it is interesting to note that both S and I converge to a value compatible with the oldest health subjects. In the future, correlating these results with the elasto-mechanic properties of the AD cortex will help to understand better the physical limits of the Mota and Herculano-Houzel model and expand it theoretically.
A significant result that could only be achieved by combining multiple data is the value of the self-similarity dimension α through the human lifespan and its evolution through atypical conditions. We extended the analysis in de Moraes et al. (2021) and verified how the slope α behaves with healthy aging. We found that the slope is compatible with the theoretical prediction for subjects between 20 and 60 years old. Defining this range of applicability of the theory is essential to understanding the period of life in which the model’s basic assumptions are still valid. Thereafter, we confirmed the findings that the slope escapes the model after 60 years old, diverging from the theoretical value of 1.25. Before 15 years old, we have found a hint that the slope is also inferior to 1.25, in disagreement with the model. In addition, a non-linear pattern is suggested by the K baseline function before 20 years. One of the leading hypotheses to explain the deviation of the measured α from the expected value in both cases (before 15 and after 60 years old) is the breakdown of the homogeneity on the cortical surface. The suggestion of a homogeneity breakdown makes it particularly interesting to further study AD in the age range compatible with 1.25. Understanding the regimes outside the model is a fundamental step toward the theoretical understanding of cortical folding. Coupled with analysis of neuropathologies such as AD, these empirical results improve the cortical folding theory.
More than discriminating between disease and healthy subjects, understanding the pathology trajectory is vital to defining landmarks, indicating future injury, and collaborating to develop efficient treatment. Notably, in Alzheimer’s Disease and its prodromal form, Mild Cognitive Impairment, subsequent phases of damage may start at least 20 years before the diagnosis (Frozza et al., 2018). The structural damage is an initial state than cognitive alterations, such as memory and clinical functions (Jack et al., 2010b). As pointed out by Fjell et al. (2014), healthy aging also contributes to reductions in cortical thickness and understanding the overlap between healthy aging and pathological aging triggers is crucial to segregating both events. Nevertheless, the possible similarity of structural damage in Alzheimer’s Disease and aging is explored by multiple previous studies (Pacheco et al., 2015; Gutierrez Becker et al., 2018; Huizinga et al., 2018; Wang et al., 2021). We hypothesize that the structural changes due to AD are reflected in measures of the brain’s global structure and shape. Therefore, cortical folding derived variables K, S, and I are sensitive to this perturbation. We verified this hypothesis by analyzing the difference in aging trajectories between healthy and pathological brains.
At a more refined grain scale, the lobes, the results are compelling to confirm (Wang et al., 2021) suggestion that K, S, and I together would be a powerful tool to discriminate similarly but distinct events, such as AD and aging. In this case, the difference between events is in the cortical structure’s changing velocity or rate. Then with the rates for K, S, and I at the lobes, we could explore whether a lobe would unfold faster. For AD, the temporal lobe has the biggest unfolding (K) and shrinking rate (I), while healthy aging is more aggressive in the Parietal lobe, followed by the temporal and frontal lobes.
The natural expansion of this study is to include Mild Cognitive Impairment (MCI) subjects. As seen before, the MCI is an intermediary step between healthy aging and AD in the disease progression, cognitive impairment, and structural changes (Jack et al., 2010b; Frozza et al., 2018; de Moraes et al., 2022). Also, we intend to continue this study by applying this approach to a longitudinal dataset composed of MRI structural images with MCI subjects that eventually convert to AD and compare to the cortical folding components of non-converters. Including this data would allow having a more reasonable aging trajectory by having the transition period from healthy aging to pathological aging. Furthermore, we expect to explore K, S, and I function with age for the age range comprising newborns to 25 years old in future studies. This study will help better understand the value of the slope α at this age, revealing if there is an inflection point that can be connected to the development/aging transition in human brains (Davis, 2021), such as suggested by the data.
This study successfully achieved its goal of combining multiple samples with a simplistic harmonization procedure that adjusts K, S, and I values for the 3,650 subjects. Hereafter we estimated K, S, and I typical values and aging rates to discriminate against Alzheimer’s Disease subjects based on their changing rates. To validate the proposed methodology, we analyzed results for Cortical Thickness, which are in concordance with previously published findings. At a glance, the significant differences found in aging rates for the cross-sectional data are suggestive that a brain with Alzheimer’s Disease has premature and accelerated aging in terms of morphology.
Data Availability Statement
The original contributions presented in the study are publicly available. Data from previous publications are available from their repositories: Wang et al. (2016, 2019a). Data from AHEAD (Alkemade et al., 2020) and AOMIC (Snoek et al., 2021) datasets with the cortical folding variables extracted by these authors are available at https://doi.org/10.5281/zenodo.5750619 (de Moraes et al., 2021). Data acquired at IDOR does not have clearance for public sharing patients’ information that could lead to identification due to local Ethics Committee approval restrictions but will be shared upon request.
Ethics Statement
The studies involving human participants were reviewed and approved by Copa D’Or Research Ethics Committee under the protocol number CAAE 47163715.0.0000.5249. Written informed consent to participate in this study was provided by the participants’ legal guardian/next of kin.
Author Contributions
FM, VM, FT-M, and BM contributed to the conception and design of the study. FM organized the database. FM and VM performed the statistical analysis and wrote the first draft of the manuscript. FT-M and BM supervised and were responsible for the project administration and funding acquisition. All authors contributed to the article and approved the submitted version.
Funding
This study was funded by the Research Support Foundation of the State of Rio de Janeiro (FAPERJ), National Council for Scientific and Technological Development (CNPq), and intramural grants from Instituto D’Or de Pesquisa e Ensino (IDOR). BM was supported by Instituto Serrapilheira (Grant Serra-1709-16981) and CNPq (PQ 2017 312837/2017-8). ADNI data included in this work collection and sharing was funded by the Alzheimer’s Disease Neuroimaging Initiative (ADNI) (National Institutes of Health Grant U01 AG024904) and DOD ADNI (Department of Defense award number W81XWH-12-2-0012). ADNI is funded by the National Institute on Aging, the National Institute of Biomedical Imaging and Bioengineering, and through generous contributions from the following: AbbVie, Azheimer’s Association; Alzheimer’s Drug Discovery Foundation; Araclon Biotech; BioClinica, Inc.; Biogen; Bristol-Myers Squibb Company; CereSpir, Inc.; Cogstate; Eisai Inc.; Elan Pharmaceuticals, Inc.; Eli Lilly and Company; EuroImmun; F. Hoffmann-La Roche Ltd and its affiliated company Genentech, Inc.; Fujirebio; GE Healthcare; IXICO Ltd.; Janssen Alzheimer Immunotherapy Research & Development, LLC.; Johnson & Johnson Pharmaceutical Research & Development LLC.; Lumosity; Lundbeck; Merck & Co., Inc.; Meso Scale Diagnostics, LLC.; NeuroRx Research; Neurotrack Technologies; Novartis Pharmaceuticals Corporation; Pfizer Inc.; Piramal Imaging; Servier; Takeda Pharmaceutical Company; and Transition Therapeutics. The Canadian Institutes of Health Research is providing funds to support ADNI clinical sites in Canada. Private sector contributions are facilitated by the Foundation for the National Institutes of Health (www.fnih.org). The grantee organization is the Northern California Institute for Research and Education, and the study is coordinated by the Alzheimer’s Therapeutic Research Institute at the University of Southern California. ADNI data are disseminated by the Laboratory for Neuro Imaging at the University of Southern California.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Acknowledgments
We acknowledge the financial support from the Research Support Foundation of the State of Rio de Janeiro (FAPERJ), and National Council for Scientific and Technological Development (CNPq), intramural grants from Instituto D’Or de Pesquisa e Ensino (IDOR), Instituto Serrapilheira, and the included datasets: the Alzheimer’s Disease Neuroimaging Initiative (ADNI), and DOD ADNI, Amsterdam Ultra-high field adult lifespan (AHEAD), Amsterdam Open MRI Collection (AOMIC), the Human Connectome Project (HCP), Information eXtraction from Images (IXI), Nathan Kline Institute (NKI) and Open Access Series of Imaging Studies (OASIS). We are thankful to the research team and volunteers for participating in this research project.
Supplementary Material
1. One may understand this as similar to “Principal Component Analysis”, a change of basis on the data representation, but using scaling considerations rather than covariance eigenvectors to define the principal directions. The new set of variables are combinations of the usual geometrical variables used to describe the cortex. The principal component K is given by the power-law relation, which contains the correlation between the variables. The variables S and I are orthonormal to K.
2. One can also start requiring a plane perpendicular to K, impose the dimensionless requirement, resulting in the dimensionless vector S, and find I from K×S.
3. “Similar to batch effects in genomics (see Leek et al. (2010) for a review of batch effects), we use the term scanner effects in neuroimaging to refer to unwanted variation that is (1) non-biological in nature and (2) associated with differential scanning equipment or parameter configurations. Because different imaging sites use different physical scanners, site effects are one example of scanner effects”.
Alkemade, A., Mulder, M. J., Groot, J. M., Isaacs, B. R., van Berendonk, N., Lute, N., et al. (2020). The Amsterdam Ultra-high field adult lifespan database (AHEAD): a freely available multimodal 7 Tesla submillimeter magnetic resonance imaging database. Neuroimage 221, 117200. doi: 10.1016/j.neuroimage.2020.117200
Dickerson, B. C., Fenstermacher, E., Salat, D. H., Wolk, D. A., Maguire, R. P., Desikan, R., et al. (2008). Detection of cortical thickness correlates of cognitive performance: reliability across MRI scan sessions, scanners, and field strengths. Neuroimage 39, 10–18. doi: 10.1016/j.neuroimage.2007.08.042
Frangou, S., Modabbernia, A., Doucet, G. E., Papachristou, E., Williams, S. C., Agartz, I., et al. (2020). Cortical thickness trajectories across the lifespan: data from 17,075 healthy individuals aged 3-90 years. bioRxiv. 43, 431–451. doi: 10.1101/2020.05.05.077834
Gronenschild, E. H. B. M., Habets, P., Jacobs, H. I. L., Mengelers, R., Rozendaal, N., Os, J., et al. (2012). The effects of FreeSurfer version, workstation type, and macintosh operating system version on anatomical volume and cortical thickness measurements. PLoS ONE 7, e38234. doi: 10.1371/journal.pone.0038234
Han, X., Jovicich, J., Salat, D., van der Kouwe, A., Quinn, B., Czanner, S., et al. (2006). Reliability of MRI-derived measurements of human cerebral cortical thickness: the effects of field strength, scanner upgrade and manufacturer. Neuroimage 32, 180–194. doi: 10.1016/j.neuroimage.2006.02.051
Marcus, D. S., Fotenos, A. F., Csernansky, J. G., Morris, J. C., and Buckner, R. L. (2010). Open access series of imaging studies (OASIS): longitudinal MRI data in nondemented and demented older adults. J. Cogn. Neurosci. 22, 2677–2684. doi: 10.1162/jocn.2009.21407
Núñez, C., Callén, A., Lombardini, F., Compta, Y., Stephan-Otto, C., and Alzheimer’s Disease Neuroimaging Initiative (2020). Different cortical gyrification patterns in Alzheimer’s disease and impact on memory performance. Ann. Neurol. 88, 67–80. doi: 10.1002/ana.25741
Pacheco, J., Goh, J. O., Kraut, M. A., Ferrucci, L., and Resnick, S. M. (2015). Greater cortical thinning in normal older adults predicts later cognitive impairment. Neurobiol. Aging 36, 903–908. doi: 10.1016/j.neurobiolaging.2014.08.031
Pomponio, R., Erus, G., Habes, M., Doshi, J., Srinivasan, D., Mamourian, E., et al. (2020). Harmonization of large MRI datasets for the analysis of brain imaging patterns throughout the lifespan. Neuroimage 208, 116450. doi: 10.1016/j.neuroimage.2019.116450
Snoek, L., van der Miesen, M. M., Beemsterboer, T., van der Leij, A., Eigenhuis, A., and Steven Scholte, H. (2021). The Amsterdam Open MRI Collection, a set of multimodal MRI datasets for individual difference analyses. Sci. Data 8, 85. doi: 10.1038/s41597-021-00870-6
I know, I got off track.” 
