Tracking Brain Development and Dimensional Psychiatric Symptoms in Children: A Longitudinal Population-Based Neuroimaging Study

Given that children who experience psychiatric problems at a young age are at an increased risk for impaired functioning and continued psychopathology later in life (1, 2), characterization of any associated neurodevelopmental features is crucial. Neuroimaging offers a unique window on in vivo brain development and the associated features of mental illness (3, 4). As the maturation of morphological (5) and white matter microstructural (6) features has been demonstrated with neuroimaging during childhood, the importance of examining emerging psychopathology in the context of typical brain development has been highlighted (7). However, information is limited on the exact interplay between the emergence of psychiatric problems and corresponding trajectories of macro- and microstructural neurodevelopment.

There is a vast literature on structural neuroimaging studies of psychopathology. Broadly speaking, externalizing disorders, such as attention deficit hyperactivity disorder (ADHD), have frequently been tied to anomalies in fronto-striatal/fronto-cerebellar circuitry (8, 9), and mood and anxiety disorders have been associated with anomalies in cortico-limbic circuitry (10). However, findings have been inconsistent, and multiple psychiatric disorders have been found to exhibit spatial overlap in alterations across a broad range of anatomical areas, including those with cortico-limbic and cortico-striatal components. While previous work has been highly informative, much of it is characterized by two notable limitations. First, most studies are cross-sectional, limiting the inferences that can be drawn about developmental processes. Specifically, it is unclear whether early neural anomalies are associated with later psychopathology, or if the reverse relationship also holds (i.e., early psychopathology is tied to later neural anomalies). Longitudinal research, which combines the collection of clinical and imaging data at baseline and follow-up, can help disentangle the temporality of this relationship (11). A second major limitation is that most studies have examined clinical samples, comparing cases and controls. However, many psychiatric symptoms exist on a continuum in the general population (12). Larger, population-based studies have demonstrated that symptoms, when considered dimensionally, vary with neurobiological features, providing further support for this framework (13, 14). However, to date, few studies have examined child psychiatric symptoms along a continuum in relation to longitudinal brain development (15, 16).

In the present study, we examined the bidirectional association of psychiatric problems with longitudinal gray and white matter microstructural development in a large sample of children from the general population. A dimensional approach was applied in the quantification of internalizing and externalizing problems, along with continuous measures of DSM symptom classes. We hypothesized that psychiatric problems along a continuum would be associated with altered anatomic and white matter microstructural development. In order to parse the precise direction of this relationship between brain and behavior over time, we used a cross-lagged panel model.

Method

Participants

This study is part of the Generation R Study, a population-based cohort study of maternal and child health from fetal life onward, in Rotterdam, the Netherlands (17). In addition to a behavioral assessment at ages 5–6 (18), a subsample of 1,070 children was recruited for MRI scanning (referred to as time 1) (19). As part of the age-10 assessment (referred to as time 2), 520 children who had a scan at time 1 also visited our research-dedicated imaging facility for MRI scanning at time 2. Figure 1 outlines the timeline of the various data collection efforts; a flow chart illustrating the exclusion of MRI data for both time points is depicted in Figure S1 in the data supplement that accompanies the online edition of this article. The final sample consisted of 845 usable T₁-weighted scans and 715 usable diffusion tensor imaging (DTI) data sets at time 1 and 480 T₁-weighted scans and 361 DTI data sets at time 2. The Medical Ethics Committee of the Erasmus Medical Center approved all study procedures, and all parents provided written informed consent.

Psychiatric Symptom Assessment

Child psychopathology was assessed using the Child Behavior Checklist (CBCL). The CBCL is a widely used 100-item inventory that provides parent-report information on a wide array of behavioral problems in young children (20). The instrument is reliable and valid, and it has been used internationally (20, 21). We utilized the broadband scales (internalizing and externalizing) and four DSM-oriented scales (the affective problem and anxiety problem scales, which correspond to DSM internalizing disorders, and the attention deficit hyperactivity problem and oppositional defiant problem scales, which correspond to DSM externalizing disorders). The DSM-oriented scales were developed to view the rated problems in the context of a formal diagnostic system (20), and they have been shown to correspond to actual clinical diagnoses (22). We administered the CBCL/1½–5 (the preschool version, for children ages 1½–5) even though some children were older than age 5 at time 1, as none of the children had yet received any formal schooling. Cronbach’s alpha was the same in 5-year-old children and children older than 5 years, indicating that problems were reliably measured in all children (23). At time 2, the CBCL/6–18 (the school-age version) was administered (24). For all analyses, raw scores (square-root transformed) were used to preserve the natural variation in the data from this nonclinical sample, as T-scores require truncation of values. The average age at the time 1 behavioral assessment was 5.9 years, and the average age at the time 2 assessment was 9.7 years. At baseline, the CBCL was administered prior to the MRI scan in all children, and at time 2, scanning preceded the behavioral assessment in a small number of participants (N=9). The percentages of participants meeting borderline and clinical cutoff thresholds for the broadband and DSM-oriented scales are presented in Table S1 in the online data supplement.

Image Acquisition

Prior to scanning, all children underwent a 30-minute mock scanning session for acclimation to the scanner environment (19). Data were acquired on 3-T scanners (at time 1 on a GE MR750, and at time 2 on a GE MR750w (General Electric, Milwaukee). Both systems utilized an eight-channel receive-only head coil. T₁-weighted structural images were acquired with an inversion recovery–prepared fast spoiled gradient recalled sequence. Diffusion MRI data were collected with three B=0 volumes and 35 noncollinear diffusion encoded volumes using an echo-planar imaging sequence (see the online data supplement for details).

Morphological Image Processing

Structural MRI data were processed through the FreeSurfer analysis suite, version 5.3 (25). Briefly, nonbrain tissue was removed, voxel intensities were normalized for B₁ inhomogeneity, whole-brain tissue segmentation was performed, and a surface-based model of the cortex was reconstructed. Global metrics of volume were extracted (e.g., total brain volume and subcortical volume), and a number of subcortical and cortical structures (amygdala, orbitofrontal cortex, etc.) were automatically labeled.

Diffusion Image Preprocessing

Diffusion image preprocessing was conducted using the FMRIB Software Library (FSL), version 5.0.9 (26), and the Camino diffusion MRI toolkit (27). Nonbrain tissue was removed, and diffusion images were corrected for eddy current–induced artifacts and translations/rotations resulting from head motion. In order to account for rotations applied to the diffusion data, the resulting transformation matrices were used to rotate the diffusion gradient direction table. The diffusion tensor was fitted at each voxel, and common scalar metrics (e.g., fractional anisotropy, mean diffusivity) were subsequently computed.

Probabilistic tractography was run on each subject’s diffusion data using the fully automated FSL plugin AutoPtx (28) (see the online data supplement). Connectivity distributions were estimated for 12 fiber bundles (Figure 2). Using the connectivity distributions, average fractional anisotropy and mean diffusivity values were then computed for each fiber bundle (29).

Image Quality Assurance

FreeSurfer reconstructions were visually inspected using a protocol similar to previously reported methods (30). Raw and processed diffusion image quality was assessed using a combination of automated and manual methods. Further details of the quality assurance procedure are available in the data supplement, and the flow chart in Figure S1 outlines the sample before and after data sets were excluded as part of quality assurance.

Statistical Analysis

Statistical analysis was conducted using the R statistical package, version 3.2.3 (31). The Lavaan package (32) was used to fit cross-lagged panel models to study the associations between longitudinal measures of brain and behavior. Cross-lagged panel models allow associations between two or more repeatedly measured variables to be investigated contemporaneously. Figure 3A depicts the general modeling strategy utilized. The first cross-lagged coefficient β_CL-1 represents the association between psychiatric problems measured at time 1 and brain metrics measured at time 2 that have been adjusted for baseline brain metrics measured at time 1. Similarly, the other cross-lagged coefficient, β_CL-2, represents the association between brain metrics measured at time 1 and psychiatric problems measured at time 2 that have been adjusted for baseline psychiatric problems measured at time 1. Cross-sectional associations between brain metrics and psychiatric problems are also modeled, although only the association at time 1 (coefficient β_CL-Baseline) is reported, given that the association at time 2 represents the correlation in the residual terms and is not straightforward to interpret. Lastly, autoregressive coefficients β_AR-CBCL and β_AR-MRI, representing the stability of psychiatric problems and brain metrics from time 1 to time 2, respectively, are modeled. While cross-lagged models can provide information on associations that relate to interindividual variability in two repeatedly measured variables, they do not provide an explicit metric of change. Therefore, significant cross-lagged associations were followed up with linear mixed-effects models in order to obtain an explicit interpretation of within-subject change. These models are described in detail in the online data supplement, and the correlations between behavioral and MRI metrics for time 1 and time 2 are presented in Tables S2 and S3.

A hierarchical approach was imposed in which broad/global metrics were first examined using cross-lagged panel models, in order to gain a comprehensive picture of the time 1 and time 2 associations. These models were subsequently followed by more refined metrics, homing in on the associations of interest using linear mixed-effects models (33). First, we tested whether broadband behavioral measures were associated with global MRI metrics (e.g., total brain volume and global fractional anisotropy). If a significant association (p<0.05) was observed with one of the broadband behavioral measures, follow-up associations between the corresponding DSM-oriented subscales and the MRI metrics were then tested. This approach was utilized in order to determine whether specific psychiatric traits account for any observed association between a broadband scale and brain metrics. Second, along similar lines, in order to determine whether the effects on the brain were global (i.e., widespread in the brain) or focal (i.e., limited to a particular set of regions or tracts), behavioral broadband scores showing a relationship with global MRI metrics were also tested for associations with changes in MRI metrics from more focal regions of interest. When examining individual regions rather than global metrics, each region of interest was z-transformed within the time point to control for confounding effects of MR scanner. Given the number of statistical tests examined with individual regions of interest, a false discovery rate correction was applied to control for type I error (34).

Global Macro- and Microstructural Brain Metrics

The anatomic metrics used were provided by FreeSurfer and included total brain, total cortical, subcortical, lateral ventricular, and white matter volume. For DTI data, rather than computing a simple average, latent factors were used to represent global DTI metrics (i.e., across multiple tracts) within the cross-lagged panel model (29, 32). Thus, all tract metrics were summarized by a single “global” factor. For each DTI scalar metric, the tracts depicted in Figure 2 were used to model the latent factors. Latent global DTI factors were modeled separately for time 1 and for time 2, given that the two assessments were acquired on different scanners, and were normally distributed. A visual depiction of the standardized factor loadings for the time 1 and time 2 global fractional anisotropy metric is available in Figure S2 in the online data supplement.

Covariates

All models were adjusted for age at assessment (behavioral and MRI, as the two assessments were conducted at different times), sex, and ethnicity (reference-coded with the Dutch population as the reference group).

Results

Table 1 summarizes the sample characteristics. Children were approximately 8 years of age at the time 1 MRI and 10 years of age at the time 2 MRI.

Neuroanatomical Macrostructure

Table 2 outlines associations between global cortical brain metrics and psychiatric problem scores. Coefficient labels in the model illustration (see Figure 3A) correspond to the headings used in the table. Cross-sectionally at baseline, higher broadband externalizing scores were associated with smaller total brain volume, cortical gray matter volume, white matter volume, and subcortical volume. Higher externalizing scores at time 1 were related to smaller subcortical volumes at time 2, after adjusting for time 1 volumes (Table 2, Figure 3B). Baseline externalizing scores predicted smaller increases in subcortical volume over time with linear mixed-effects models, and results remained consistent after adjusting for intracranial volume (see Table S4 in the data supplement). Interestingly, the path testing whether baseline neuroanatomical features predicted time 2 externalizing scores was nonsignificant.

For broadband internalizing scores, cross-sectionally at baseline, there were no significant associations with any of the macrostructural features. However, higher internalizing scores at time 1 were associated with smaller subcortical volumes at time 2 after adjusting for volumes at time 1 (Table 2). The subcortical association was consistent in linear mixed-effects models, even after further adjusting for intracranial volume (see Table S4). The path testing associations between global metrics of cortical morphology at time 1 and internalizing scores at time 2 was nonsignificant.

In order to better characterize the regional specificity of the significant subcortical volume analyses outlined above, linear mixed-effects models were run to predict change in the individual regions that comprise the total subcortical volume. After adjustment for multiple comparisons, psychiatric problems at time 1 were not related to changes in any of the subregions over time.

White Matter Microstructure

Externalizing scores at time 1 were not associated cross-sectionally with global DTI measures at time 1. However, externalizing scores at time 1 were negatively associated with global fractional anisotropy at time 2 after adjusting for global fractional anisotropy at time 1 (Table 2). Linear mixed-effects analyses were consistent, with higher baseline externalizing scores predicting smaller increases in global fractional anisotropy over time (Figure 4; see also Table S5 in the data supplement). Baseline global DTI measures were not associated with externalizing scores at time 2.

Internalizing scores at time 1 were also not cross-sectionally associated with global DTI measures at time 1. Broadband internalizing scores were negatively associated with global fractional anisotropy at time 2 after adjusting for global fractional anisotropy at time 1 (Table 2, Figure 3C). In linear mixed-effects models, higher baseline internalizing scores predicted smaller increases in global fractional anisotropy over time (see Figure 4 and Table S5). Similar to what was observed with the neuroanatomical features, time 1 global fractional anisotropy did not predict time 2 internalizing problems (Table 2).

Given that associations were observed between time 1 behavioral measures and time 2 DTI metrics, individual tracts were examined to determine whether there was any regional specificity in the white matter associations. Higher externalizing scores at time 1 were associated with smaller increases in fractional anisotropy in the superior longitudinal fasciculus (see Table S6 in the data supplement). Higher broadband internalizing scores at time 1 were associated with smaller increases in fractional anisotropy in the right cingulum and in the superior longitudinal fasciculus bilaterally (see Table S6).

The significant associations presented above were adjusted further for potential confounders (e.g., IQ and motion), and the results are presented in the data supplement.

DSM-Oriented Subscales

To further characterize the global associations, we examined the individual DSM-oriented subscales when broadband scales showed significant associations with changes in global MRI metrics. Higher baseline scores on both the DSM-oriented attention and oppositional defiant subscales predicted smaller increases in global fractional anisotropy over time (see Table S7 in the data supplement). Higher baseline affective subscale scores were associated with smaller increases in global fractional anisotropy (see Table S7).

Stability of Psychiatric Scores and Brain Metrics Over Time

Autoregressive coefficients in Table 2 demonstrate that psychiatric scores measured at time 1 are positively associated with those measured at time 2 (roughly 0.65 for externalizing scores and 0.56 for internalizing scores), suggesting stability in the measure and also some interindividual variability over time. Similarly, brain metrics at time 1 were positively associated with those measured at time 2, suggesting a relatively high stability with some interindividual variability over time, particularly in global fractional anisotropy.

Discussion

In this large population imaging study of children, we demonstrate a link between psychiatric problems along a continuum and a differential pattern of brain changes over time. Even in the general population, psychiatric problems were related to altered trajectories of both macro- and microstructural brain development. Interestingly, baseline brain metrics measured during childhood did not predict changes in psychiatric symptom ratings over time; instead, psychiatric problems at a young age predicted an altered pattern of brain changes over a 2½-year interval.

Consistent with the literature, this study shows that facets of psychopathology are related to smaller subcortical volumes and lower fractional anisotropy. The subcortex has been implicated to some degree in nearly all psychiatric disorders, and brain imaging studies have revealed evidence for its involvement in, for example, depression, ADHD, and obsessive-compulsive disorder (8, 35, 36). Similarly, numerous studies have demonstrated the potential role of white matter microstructure in psychopathology (10, 37). The present study expands on this literature by demonstrating that dimensionally assessed psychiatric problems are related to smaller changes in brain volume and fractional anisotropy over time. Most subcortical structures show increases in volume that peak during late childhood and into adolescence, followed by decreases in volume into adulthood (38). Within the relatively narrow age range of the present sample, we show associations that are consistent with this previously reported pattern, and additionally show attenuated increases in subcortical volume that are related to both externalizing and internalizing problems. Fractional anisotropy in white matter has shown largely linear increases over time within this restricted age range (6). Data from the present study fit with the literature, again additionally showing that these trajectories are potentially modified by the presence of psychiatric symptoms at an early age. Fractional anisotropy, but not mean diffusivity, showed a differential pattern of change. This could suggest differences in, for example, myelination or axonal packing.

Limited specificity was observed, as both externalizing and internalizing symptoms were associated with differential macro- and microstructural changes over time. Furthermore, three of the four DSM scales were associated with changes in white matter microstructure. As the two broadband domains are relatively highly correlated (20), it is possible that the overlapping variance in these scales best explains the differential changes in the brain. Alternatively, given the relatively young age of the sample, it is also possible that domain-specific patterns of development become more apparent at older ages. In terms of symptom assessment, this study used parent reports of children’s psychiatric symptoms. It is likely that combining data from multiple informants (parents, teachers, and the children themselves) would provide a more robust and accurate picture of their symptom profile, which might in turn help to elucidate more specific correlates with neuroimaging features (39). Along similar lines, although this study used high-resolution anatomical data and diffusion tensor imaging data, the two modalities were examined separately. It is possible that a multimodal approach in which both pieces of information were elegantly combined into a single analysis could help to disentangle unique neurobiological features in both categories of disorders (40).

In the context of psychopathology, the most widely applied model in the field of neuroimaging has focused on identifying associated underlying neurobiological substrates of potential etiological significance. However, we did not find that early brain metrics predicted changes in psychiatric problems over time. There are several potential explanations for this. As this is not a clinical sample, it is possible that such processes are not part of a continuum and are not present in subclinical presentations of psychiatric constructs in the general population. It is thus a priority to determine whether these findings hold in large cohorts that are enriched for psychiatric disorders. Alternatively, the underlying neurobiological predictors of the emergence and development of psychopathology may be spatially or temporally heterogeneous, particularly during brain development, and more sophisticated image analysis methods may be needed to better characterize them (e.g., machine learning) (41). It is also possible that the responsible underlying brain features are detectable at higher spatial resolution or perhaps using alternative imaging modalities (e.g., MR spectroscopy, MR perfusion, task-based or resting-state functional MRI).

The bulk of the findings in this study consisted of early measures of psychiatric problems predicting differential patterns of both macro- and microstructural brain development. This suggests that, in addition to the standard model of the brain shaping behavior discussed above, perhaps behaviors also shape the brain (42). While there are undoubtedly underlying neurobiological explanations for the emergence of psychopathology in children, the data in this study suggest that psychiatric symptoms in children may also contribute to some of the macro- and microstructural abnormalities reported in the literature—a potentially cascading interaction between psychopathology and the brain. Take the example of obsessive-compulsive disorder, where the symptomatology originates in a particular set of brain areas (e.g., fronto-striatal circuitry), but over time the symptoms themselves could modulate structural brain development (e.g., repetitive motor behaviors leading to alterations in the motor cortex) (42). A similar extension can be drawn, in anxiety disorders, for example, to the potential effect of increased hypothalamic-pituitary-adrenal axis activity on the brain through downstream hormonal exposure (43). A child with psychopathology may also experience his or her environment differently as a result of the disorder, which could consequently influence how the brain is shaped during development (e.g., reduced novelty seeking). The confounding effect of psychiatric treatment could also play a role (i.e., children with more psychiatric symptoms seek treatment, which subsequently influences brain development over time), although additional analyses in this study did not show that psychiatric medication status explained the association (see the data supplement). Future studies should not rule out the possibility that, in addition to a preceding causal factor, observed neuroimaging associations with a given disorder could also be a downstream effect of the disorder.

Despite the large population-based sample, longitudinal design, and dimensional assessment of psychiatric symptoms, this study has important limitations. First, although the data for this study are part of a larger population-based study, representativeness is important to discuss. While many features of this sample are representative of the catchment area (e.g., ethnicity), some potentially are less so. For example, the mean nonverbal IQ in this sample was a few points above the population average. Although models adjusted for IQ proved to be similar to unadjusted models (see the data supplement), IQ remains an important consideration in child psychiatric studies, as it is often intertwined with clinical diagnosis and is related to many neuroimaging features. Next, as a hierarchical approach was implemented, some associations at a more focal level could have been missed; future work will benefit from analyzing brain and behavior at a more focal level. Regarding data acquisition, imaging data were acquired on two different MRI scanners, possibly introducing problems with the longitudinal interpretation of results. However, acquisitions were made as similarly as possible (e.g., gradient table, head coil, etc.), and a number of steps were introduced to mitigate such problems, including within-scanner normalization of MRI metrics (e.g., latent factor modeling within scanner, and Z-score standardization), which should attenuate the effects of scanner difference. Additionally, the statistical models used in this study to test the association between psychiatric problems and change in macro- and microstructural brain features operate on a relative rather than an absolute scale, further mitigating concerns about the effect of scanner. Lastly, only linear models were tested, as data were acquired at only two points in time. Gray matter morphological studies have shown the importance of nonlinear trajectories, and future work with additional time points should address this in more detail.

In conclusion, we demonstrate that internalizing and externalizing problems are related to altered macro- and microstructural changes in a large sample of children from the general population. Tracking the emergence of psychopathology in children, both in terms of symptomatology and neurobiology, may not only help guide diagnosticians but also improve the selection and timing of treatments. It is important to appreciate that this study does not challenge the causal role of neural changes in the pathogenesis of psychiatric disorders. Nonetheless, our findings raise the intriguing possibility that some emergent neuroimaging features of psychopathology may be partly a consequence of the disorder.