Charting brain growth in tandem with brain templates for schoolchildren

Brain growth charts and age-normed brain templates are essential resources for researchers to eventually contribute to the care of individuals with atypical developmental trajectories. The present work generates age-normed brain templates for children and adolescents at one-year intervals and the corresponding growth charts to investigate the influences of age and ethnicity using a common pediatric neuroimaging protocol. Two accelerated longitudinal cohorts with the identical experimental design were implemented in four centers in the United States and China. Anatomical magnetic resonance imaging (MRI) of typically developing school age children (TDC) was obtained up to three times at nominal intervals of 1.25 years. The protocol generated and compared populationand age-specific brain templates and growth charts, respectively. A total of 674 Chinese pediatric MRI scans were obtained from 457 Chinese TDC and 190 American pediatric MRI scans were obtained from 133 American TDC. Populationand age-specific brain templates were used to quantify warp cost, the differences between individual brains and brain templates. Volumetric growth charts for labeled brain network areas were generated. Shape analyses of cost functions supported the necessity of age-specific and ethnicity-matched brain templates, which was confirmed by growth chart analyses. Analyses revealed volumetric growth differences between the two ethnicities primarily in lateral frontal and parietal areas, regions which are most variable across individuals in regard to structure and function. Ageand ethnicity-specific brain templates facilitate establishing unbiased pediatric brain growth charts, indicating the necessity of the brain charts and brain templates generated in tandem. These templates and growth charts as well as related codes have been made freely available to the public for open neuroscience (https://github.com/zuoxinian/CCS/tree/master/H3/GrowthCharts).


Introduction
Growth charts are an invaluable resource for pediatric clinicians. They are essential for screening the developmental status of individuals and monitoring for early detection of abnormal growth [1]. Deviations from normative age-expected values trigger evaluations searching for underlying abnormal factors, facilitating early detection and intervention when required. Extending this approach to the evaluation of neurodevelopmental status has been impeded by the lack of reliable brain growth charts. Magnetic resonance imaging (MRI) is increasingly being employed to map human brain development. Anatomical MRI (aMRI) can capture developmental changes in brain morphology [2,3], comprising full-brain geometrical transformations (e.g., cortical thinning and surface expansion) [4,5]. For example, changes in cortical thinning trajectories have been linked with inter-individual differences in IQ in children and adolescents [6]. Such developmental effects in brain structure have also been shown to be detectable across adulthood [7] and are supported by brain network studies using diffusion-weighted (dMRI) and restingstate functional (rfMRI) imaging methods [8,9], providing a framework for quantifying multimodal brain development at the population level [10,11]. Although sparse, efforts to translate developmental trajectories into GCs have begun to be initiated for neuropsychiatric conditions [12][13][14], which are hypothesized to have abnormal neurodevelopmental origins [15,16].
Despite the promise of developmental population neuroscience [17,18], a number of key issues must be addressed prior to establishing brain GCs for clinical use. First, reliability of MRI-based measurements must meet clinical standards for capturing intra-and inter-individual differences [19,20]. Core anatomic MRI measures (e.g., volume, cortical thickness, surface area) currently meet this standard [21], but most dMRI and rfMRI measures do not, due to multiple confounds and substantial measurement error [22,23]. This suggests aMRI-derived measures could provide the bases for developing reliable and valid imaging GCs [20]. Second, MRI samples of brain development cohorts for building GCs are currently limited. Largescale brain development cohorts are fundamental for charting growth [24][25][26], but obtaining longitudinal assessments across multiple centers with the same protocols is rare [27,28]. Third, the necessity of deriving GCs for height, weight and head circumference for specific populations or countries is established [29,30], and is likely even more critical for building brain GCs, given the neurodevelopmental diversity potentially associated with differences in ethnicity and culture. Fourth, detection of differences could be biased by using inappropriate templates, as has been documented for developmental ages [31][32][33][34][35], indicating the need for agespecific brain templates. Finally, despite the public health importance of creating normative GCs [12,[36][37][38][39], whether and how brain templates affect the estimation of brain morphological GCs remains unknown.
This study was designed to begin to fill this gap by developing an integrative pipeline for generating brain templates and GCs in schoolchildren. Volumetric measurements were quantified with aMRI of 603 school-age scans from two accelerated longitudinal cohorts with the same experimental design obtained in the United States (Enhanced Nathan Kline Institute Rockland Sample -eNKI sample, 190 scans) [40] and China (Chinese Color Nest Project -CCNP, 413 scans) [11,41], respectively. Standard brain templates were constructed for each year of age that can serve as a resource for generating volumetric GCs of tissue compartments (white or gray matter), lobes and large-scale brain networks. These brain templates and GCs were further validated across two cultures to offer an initial normative reference for studies of brain development during school age.  [11], which is a five-year accelerated longitudinal study designed to delineate normative trajectories of brain development of Chinese children [41]. Each participant in CCNP is invited to undergo scans three times at nominal 15-month intervals. The age and sex of overall MRI scans are listed in Table 1.

MRI
For the eNKI Rockland Sample (from 2013 to 2017), a total of 190 scans from 133 TDC were included after the same QCP applied for CCNP samples in final analyses. CCNP and eNKI datasets are both accelerated longitudinal designs, and were initially designed with matched age span and imaging resolution. Participants in the CCNP and eNKI sample who had a history of neurological or mental disorder, family history of such disorders, organic brain diseases, physical contraindication to MRI scanning, a total Child Behavior Checklist (CBCL) T-score higher than 70, or a   Wechsler Intelligence Scale for Children IQ standard score lower than 80 were   excluded. CCNP and eNKI projects obtained Institutional Review Board approval from IPCAS and NKI, respectively. Written informed assent and consent were obtained from both participants and their parents/guardians. Sex and ethnicity information was collected from self-and parent/guardian report. Details of the PKU131 sample has been reported previously [42,43] and also can be found in ADHD200 consortium, The CCNP-SWU413 and eNKI samples, which were matched on age and imaging resolution and used the identical experimental design (Table 2 and Figure 1) were employed for growth chart modeling. As few children were in the two youngest ages (n=7 at age 6, n=22 at age 7 in CCNP sample, n=4 at age 6, n=3 at age 7 in eNKI sample), the two age groups were combined.
Two additional pediatric neuroimaging datasets from Weifang Medical University [45] and Zhejiang University [46,47], including 84 structural MRI scans, were employed to validate (validation data) the necessity of constructing age-and ethnicity-matched MRI templates using brain deformation cost function (see age, sex in Table 3.

Analysis
To chart brain growth, we developed a standard pipeline consisting of customized brain template construction [48,49], robust imaging registration and growth chart estimation (Fig. S1 online). To incorporate atlas information, we also performed a two-step protocol of image registration (steps 5-7 in Fig. S1 online). Detailed analysis procedures can be found in Supplementary Materials. Ethnicity and age are the two major variables addressed in this work. We examined a 2 (ethnic levels, CCNP vs. eNKI) × 11 (age levels) within-subject design to test template effects in registration. The 11 age levels ranged from 6 to 17 years old; ages 6 and 7 were combined to increase sample size. This generated 22 ethnicity-and age-specific templates. Individual brains from the validation dataset were then registered to these 22 templates, with 22 corresponding registration deformations calculated for each subject. To assess the influences of the brain template shapes on the registration across different age groups, we performed affine registration on the IPCAS pediatric templates to match the size of the Montreal Neurological Institute (MNI152) adult brain template. We also include the Chinese pediatric brain templates constructed by Zhao et al. [35], which are transformed to match the shape and size of the MNI152 template. We denoted the former as IPCAS-Affine templates and the latter as BNU templates. A total of 4 groups of warp deformation were calculated to determine the optimal template for Chinese pediatric population.
For modeling brain morphological growth charts, ethnicity was also considered as a potential confounder. Due to the absence of ethnic-specific templates, the MNI152 template has been usually utilized as the default in most previous studies. To test how an ethnicity-unspecified template affects morphological estimation, we performed a 2 × 2 mixed design with ethnicity as the between-group factor (CCNP vs. eNKI) and appropriateness as the withingroup factor (individual brains registered to ethnic-matched and ethnicmismatched template brains). Two registrations were performed for each participant. Specifically, for a child from the CCNP sample, ethnicity-appropriate registration refers to using an IPCAS age-appropriate template while ethnicity-inappropriate refers to using an NKI age-appropriate template, and similarly for youth in the eNKI sample. To achieve more detailed characteristics on brain growth, we employed the 400-unit areal parcellation in MNI space, namely localglobal parcellation (LGP) [50]. It combined both local and global properties of intrinsic functional connectivity and parcellated the human cerebral cortex into 400 parcellation elements (parcels) with almost equal surface areas. This parcellation has two advantages: (1) clear brain network attribution [51], and (2) usability validated for studies across different stages of the human lifespan (6 -85 years) [50]. Individual-level LGPs were extracted based on the above two registrations for each subject and their GCs were modeled, generating four GCs for each brain area: (1) CCNP2IPCAS (CCNP brains registered to the IPCAS templates), (2) eNKI2NKI (eNKI brains registered to the NKI template), (3) CCNP2NKI (CCNP brains registered to the NKI template),( 4) eNKI2IPCAS (eNKI brains registered to the IPCAS template). We hypothesized that the former two charts would be more appropriate than the latter two charts.
Quantile regression was employed to build brain GCs [52]. Analysis was performed using GAMLSS implemented in R (version 3.4.3) [53]. Two models were applied to explore developmental trajectories. In one model, volume data from all participants was utilized for modeling GCs, while in the other, growth curves were modeled separately for boys and girls. This analytic strategy has been employed by the World Health Organization (WHO) and the Centers for Disease Control and Prevention (CDC) to delineate height and weight growth charts [54][55][56].
To quantitatively estimate the diversity in GCs attributable to ethnicity, normalized variance (NV) was calculated across 400 brain parcels in MNI space [50] contrasting CCNP and eNKI samples, with the NV values calculated as follows: where V is a vector referring to the parcel unit volume at every age point estimated in the last step; the standard deviation of the differences between the two samples was calculated to characterize the degree of chart shape dispersion across different ages. It was further normalized by the mean parcel volume across the two cohorts. A large NV quantifies the shape difference of the two GCs while a small NV indicates that the two GCs share similar shapes.

Results
The constructions of the IPCAS and NKI brain templates were based on 674 and 190 MRI scans, respectively; detailed sex and age decompositions are shown in Table 1 were utilized for GC modeling. Standard brain and tissue probability templates, and GCs for brain tissues (gray matter (GM) and white matter (WM) and cerebrospinal fluid (CSF)) and 400 areas (all areas assigned to the seven brain networks identified by Yeo et al. [51]) were generated in the current study. The axial slices of agespecific templates are illustrated in Figure 2a. The upper four rows were built from the CCNP cohort while the lower four rows were constructed from the eNKI cohort.
Clear differences in tissue spatial profiles across childhood and adolescence are observable, and ethnic differences in brain shape can be appreciated. To further examine whether these template morphological differences induced by ethnicity changed with ages, we segmented the age-specific brain templates into GM, WM and CSF and demonstrated their volume changes across ages between ethnicities ( Figure 2b). The ethnic differences in age-related changes of the three brain tissues are observable for specific age ranges while their overall patterns are similar: decreasing GM volumes, increasing CSF volumes, and increasing then stabilizing WM volumes.
Multi-comparisons reveal significant differences between each pair of groups; details are shown in the violin plots in Figure 3, right panel and Table 4. As illustrated in Figure 3, small deformations occurred when sources and targets of brain registration were approximately age matched, but the smallest deformation occurred when individual brain images were registered to templates corresponding to 1-year younger. Negative matching ages refers to templates based on younger samples, whether CCNP or eNKI.
The two templates with the same brain size of the MNI152 adult template, IPCAS-Affine and BNU template, exhibited comparable costs of registration deformation with each other. However, they showed much higher deformation costs than the eNKI and IPCAS templates (Fig. S2 online). The deformation cost curves of IPCAS and IPCAS-Affine exhibited identical shapes. The BNU templates were constructed with samples younger than 13 years old, and thus only a slightly decreasing trend were observed at this age span, which is almost matched to the IPCAS-Affine templates at the same age span. The only difference between IPCAS and IPCAS-Affine templates is the template size, thus it suggests that the improper size of the template can introduce larger registration deformation. The fact that warp deformation induced by IPCAS-Affine templates is even larger than that induced by eNKI templates indicates that the developmental factor (i.e., size) may have more impacts on registration.

Differences in brain growth charts between cohorts
Growth charts including global metrics like intracranial volume (ICV), GM volume, WM volume and CSF are depicted in Figure 4 (first row) for CCNP and eNKI samples with all subjects combined and separately for males and females in Fig. S3 online. Overall, brain volume was larger in the eNKI than in the CCNP sample.
Intracranial volume and WM volume grow rapidly during childhood and slowly in adolescence in the eNKI dataset. In the CCNP dataset, few age-related differences were observed for ICV and WM. For GM, both datasets exhibited linear decreases.
CSF volume exhibited an inverted pattern from other tissues, with larger volume in CCNP participants, rapid enlargement during childhood and stable in adolescence.
By contrast, in the eNKI dataset, CSF volume increased continuously throughout childhood and adolescence. Developmental trajectories at the level of brain lobes are shown in Figure 4 (second row), which exhibits a similar linear decreasing pattern across samples; sex differences in growth charts are shown in Fig. S4 online.
At a more resolved scale, GCs of regional brain volumes were compared.
The similarity of trajectories between CCNP and eNKI cohorts was estimated for each area and depicted in Figure 5. Large differences were mostly observed in association cortex while primary cortex exhibited similar developmental trajectories.
To better summarize the distribution of NVs among common brain networks, the pie graphs of regional NVs that are larger than their mean value are also shown in

Growth charts are driven by developmental underpinnings of brain templates
Potential effects of using age-specific brain templates on GCs were evident. Such effects are reflected in results dramatically driven by age-specific brain templates.
For instance, we consider an area located in the right superior parietal gyrus (RH_Cont_Par_5 of the frontal-parietal network in Schaefer et al. [50]) ( Figure 5).
Its growth curve exhibited relatively distinct patterns between CCNP and eNKI cohorts (NV = 0.29) when individual brains were registered to ethnic-and ageappropriate templates (right panel). Growth patterns were almost inverted when individual brains were deformed to mismatched ethnic-templates (left and middle panels). This implies that the GC patterns were largely driven by the developmental underpinnings of the brain templates used. The observation that extraction of the areal metric largely depended on the target templates used for registration held generally across the whole brain.

Discussion
Growth charts built with existing big datasets provide clear evidence that it is necessary to estimate brain morphological properties as well as their corresponding brain templates within specific age groups. This is especially true for pediatric neuroimaging studies of school-age participants, who undergo considerable changes in brain morphology. Both cross-sectional and longitudinal applications of growth charts are facilitated by choosing proper templates to define regional areas in individual brains.
Previous studies [33,34,57] have demonstrated that different ethnicities and ages increase deformation costs associated with morphing anatomical regions of individual brains, which if done poorly, can result in mismatches in brain segmentation tissue profiles [33,34]. This is supported by our findings that even with an identical brain atlas, morphological metrics can differ substantially when registered to different brain templates. Ideally, these metrics should be identical. In practice, the method with lower registration errors or costs is preferred. However, the desirability of population-and age-specific brain templates for modeling growth charts has not been prioritized. In the single exception, group differences in registration errors relating to ethnic and developmental factors were tested [34]. In that study, the authors applied group-level comparisons in which registrations were divided into appropriate and inappropriate groups for estimating template effects, with paired t-tests or variance analysis performed to assess ethnicity differences in template registration. The utility and generalizability of their templates was limited by the samples from a single site using broader age intervals (2 years). Moreover, the relationship between registration errors and age in pediatric samples has yet to be examined and quantified. Therefore, we strongly recommend creating custom brain templates from a homogenous population across the appropriate age ranges to improve registration performance [58] and the accuracy of brain growth charts in pediatric MRI studies.
More samples are included; more representative brain templates are constructed.
In adult population, Yang et al. [59] first revealed that the stability of brain template is improved by increasing sample size, suggesting a large sample size is preferred in constructing brain template. This could be even more complicated in pediatric population as brain morphology keeps changing during childhood and adolescence for higher inter-individual differences. It still remains an open question whether an optimal age span or sample size existed for the pediatric template construction and how developmental factors affect its stability. In the current study, we have used the available samples as many as possible to increase the sample size. So far as we known, CCNP is the largest longitudinal Chinese pediatric MRI dataset. To address the potential limitation of sample size, we also included two other public available Chinese pediatric datasets (SMU130 and PKU131). The sample size used to construct Chinese pediatric brain templates in the current study is the largest among all the similar studies. The number of scans (n = 674) is nearly 5 times and 13 times larger than the previous studies by Xie et al. [34] (n = 138) and Luo et al. [32] (n = 53), respectively. The templates constructed by Zhao et al. [35] are based upon an elegant algorithm but also with limited sample sizes (n < 55 per group) for the six age groups (6-12 years). With the accumulation of publicly shared Chinese pediatric MRI datasets, we will update our brain templates and apply more advanced algorithm to improve the template quality in future.
In our analyses, age and ethnicity accounted for the most variance across individuals. Ethnicity plays a critical role in shaping brain morphology [57]. For instance, significant volumetric differences were observed between Chinese and Caucasian adult brain templates [33,57], reflecting rounder global shape and shorter axial distance in Chinese adults. Dynamic neurodevelopmental factors, namely age factors, affect brain maturation, suggesting such brain morphological differences should also be observable during childhood or adolescence. Both of these effects were revealed in our results. As demonstrated in Figures 2 and 3, we observed significant shape distinctions between IPCAS and NKI brain templates ( Figure 2).
Moreover, the corresponding registration costs were also related to national origin and stage of development ( Figure 3). This was particularly well illustrated by using age-specific brain templates for the longitudinal CCNP and eNKI samples to model growth charts of brain volume (Figure 4), which were similar (peaking at 12-13 years of age but differing in specific details) to the inverted-U shaped curves observed in previous studies [60,61]. The eNKI sample exhibited larger brain volumes and more accelerated increases during childhood than the CCNP sample ( Figure 2). Differences in such a fundamental morphological characteristic may lead to increases of registration errors related to age and ethnic differences ( Figure   3). Several factors may account for the observation that the minimum deformation value occurred at when individual brain images were registered to the templates corresponding to one-year younger. First, the validation dataset only included participants below 12 years old and were most collected about 10 years ago.
Social-economic development has an inevitable impact on human health, especially on brain maturation [62]. The validation dataset may have a systematic bias comparing with the dataset collected in recent years, this is also the reason why the validation dataset was not included in the construction of brain templates.
Second, the limited sample size may also result in large inter-subject variability, which may be another factor that leads to a negative matching age in registration deformation.
At a more resolved scale, the areas showing large developmental divergences in growth charts between CCNP and eNKI samples are almost located within highlevel hierarchical cortical regions such as default mode network and frontal parietal network. Other areas situate at the boundaries between networks also show moderate divergences. This phenomenon may reflect several facts. First, from a methodological aspect, the parametric (linear or quadratic) models were applied to estimate growth charts in most previous studies. We employed a semi-parametric model without setting a strong prior to construct growth charts, which is more sensitive to the sample data. This may be less stable with a limited sample size; however, the linear trajectories were still captured as shown in Figure 5 for the visual network's 19 parcel situated as well as other areas with low NVs. Large-scale measurements like brain lobe and network level often show monotonic trend in trajectory [63], which is also observed and replicated in our study. However, few studies have focused on growth charts at such a detailed parcellation level, and smaller brain areas may be highly valuable for revealing more complex patterns of higher growth variability in specific functional areas. Cognitive abilities and selfconsciousness are rapidly developed during the childhood and adolescence period.
The dynamic developmental trajectories of brain morphology may reflect the potential integration and differentiation of cortical function. This may also be the reason why areas located within unimodal regions like visual cortex and somatosensory cortex exhibited almost the identical trajectories. Second, we could speculate that potential socio-cultural and economic status may also account for this divergence. It is well known that socio-economic status is closely associated with human health [62,64], i.e. height of children and adolescents improved in tandem with socioeconomic development in China [65]. Although it is not fully understood that how socio-economic status affects brain morphology during development, it has a huge impact on human brain maturation. In summary, the divergence in growth charts between ethnicities might be an indicator of the differences in socio-culturaleconomic status, which needs to be further explored in future work.
Many factors can affect the construction of pediatric growth charts, including the data preprocessing pipeline [66], modeling methods [11], and site effects [63].
Dynamic developmental trajectories might be confounded by image registration errors if inappropriate brain templates are employed. However, ethnic-and agespecific brain templates have not been used in previous developmental studies.
This is partially because small samples are insufficient to construct such templates and few developmental studies have focused on such areal scales (small parcels) [50]. In the present study, we showed that for growth chart modeling, use of improper brain templates can substantially distort the estimations of underlying morphological development, making conclusions questionable. Hence, we believe that the construction of age-specific brain templates and developmental trajectories or growth charts should be performed in tandem.

Limitation and future work
Several limitations must be considered regarding the application of pediatric templates and the interpretation of brain growth charts. The proportion of males and females was balanced in most age groups except for 15-, 16-and 17-year-olds.
Given sex differences in brain development [61,64,67], sex-specific templates should also be constructed. Also, the age intervals used to define templates in the present studies were defined provisionally due to the lack of more detailed evidence on brain development. Image variance introduced by different scanners/sites is another potential confounding factor reducing the accuracy of growth charts estimation. It is still unclear that whether and how systemic bias in scanner and parameters would affect the estimation of growth charts, however, incorporation of multi-center datasets in the future would greatly facilitate the generalization of growth charts application. Nonetheless, the age-specific brain templates generated in the current study can facilitate the estimation of more precise changes in dynamical growth pattern of human brain morphology during development. Regional volume was employed in this protocol to demonstrate age and ethnicity effects on brain templates and growth charts. How such effects can be generalized to other metrics of human brain morphometry (see reference [68] for a review) and function is worth considering.

Conclusions
Age-and ethnicity-specific brain templates are required to establish in tandem with unbiased brain growth charts for pediatric population. Specifically, using a large neuroimaging dataset of Chinese and United States pediatric brain images, we demonstrated for the first time that greater age mismatching of templates introduces larger registration deformations. Age-specific templates can improve the accuracy of image registration between individual pediatric brain images, thereby facilitating more reliable and accurate human brain mapping studies in healthy and clinical pediatric populations [36,69]. By modeling growth charts, we found that differences across western and eastern samples were decreased when examined at large-scale levels, including tissue classes of brain lobe volumes. At more fine-grained levels of spatial resolution, ethnic differences in cortical surface area indices became markers of developmental divergence, particularly in association cortex, which exhibits greater flexibility, morphological variability and hemispheric asymmetry [70]. To be noted, it is the growth pattern rather than the absolute estimation of brain morphological metric that is emphasized in the current report. The identification of deviation from typical brain growth pattern is the main target of growth charts application.

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. White helped to polish the paper. All authors discussed the results and commented on the manuscript.         Similarities of brain growth charts between CCNP and eNKI samples. The first and third row shows the 400 brain parcellation units, with parcels colored according to the NV values and the Yeo2011 seven-network organization. The second row shows similarities of the brain growth charts between CCNP and eNKI samples as measured by NV. The second row depicts growth charts with largest and smallest NVs. The pie graph demonstrates the NV ratio of networks. The bottom row shows the template effect in growth charts modeling, take parcel Cont_Par_5 for instance, CCNP2NKI indicates the subjects from CCNP sample registered to NKI template, eNKI2IPCAS indicates the subjects from NKI sample registered to IPCAS template. Registration to same template drives the pattern of growth charts similar.