Developmental Trajectories of Anxiety and Depression Symptoms from Early to Middle Childhood: a Population-Based Cohort Study in the Netherlands

This study was embedded in the Generation R Study, an ongoing multi-ethnic population-based prospective cohort from fetal life onward in Rotterdam, the Netherlands. The Generation R Study is designed to identify early causes of normal and abnormal growth, development, and health. Its design has been previously described in detail (Kooijman et al. 2016). In short, all pregnant women living in Rotterdam, with an expected delivery date between April 2002 and January 2006, were invited to participate and contacted via obstetrician practices. The study was approved by the Medical Ethics Committee of the Erasmus Medical Center. Anonymity was guaranteed and written informed consent was obtained from children’s primary caregivers at each study phase.

In total, 9778 mothers enrolled in the study gave birth to 9749 live-born children. During the preschool period (0–4 years), parents of the 1166 children living outside the definite study area at birth were not approached as the logistics of postnatal follow-up were embedded in the municipal routine child care system. When all parents were contacted again when children were 5 years, consent for their children to participate was 85% of the original sample (Jaddoe et al. 2012). Data on children’s anxiety and depression symptoms were available at ages 1.5 years (n = 5223), 3 years (n = 4939), 5 to 7 years (n = 6210; hereafter referred to as ‘age 6’) and 10 years (n = 4938). In total, the study sample comprised of 7499 children with anxiety and depression data available for one or more research assessment(s) (76.2% response rate). For psychosocial and school-related outcomes, the study sample varied across outcomes (n = 4336 friendship quality, n = 4355 self-esteem, n = 3669 school performance, and n = 3857 school problems).

Categorical omega’s for internal consistency were computed with the MBESS R package version 4.4.3 (Kelley 2018). Descriptive analyses were conducted with SPSS version 21.0 (IBM Corp. 2016). All other analyses were performed with Mplus version 8 (Muthén and Muthén 2017). We used growth mixture modeling (GMM) to determine trajectories of anxiety and depression symptoms in early childhood (1.5–10 years). To be able to model these trajectories, average scores on the Anxious/Depressed scale of the CBCL/1 ½-5 and CBCL/6–8 from the four assessments were used when at least 75% of the items were completed. GMM is a statistical technique that can be used to find subgroups of change in individuals over time while allowing for within-class or trajectory variability (Fanti and Henrich 2010; Jung and Wickrama 2008).

To identify the best model fit, we compared models with different trajectory solutions with the use of both the Bayesian information criterion (BIC) and the Lo-Mendell-Rubin test (LMRT). A lower BIC score indicates a model with a better fit and the LMRT tests indicates whether a given model has a significantly better fit than a solution with one trajectory less. In addition, a plot of trajectories for each model was produced to see whether these solutions were distinctive in capturing the development of anxiety and depression symptoms. Each model was fitted with a linear slope and with a quadratic slope in addition. Classification accuracy was determined by both mean assignment probabilities for each trajectory and overall model entropy. Entropy with values approaching one indicate a clear delineation of trajectories (Celeux and Soromenho 1996). Because of the low internal consistency of the CBCL Anxious/Depressed subscale at age 1.5 years, we compared trajectories with four assessments from 1,5 to 10 years with trajectories based on the sample of children from the three assessments 3, 6, and 10 years (n = 7130).

We then examined gender invariance by determining whether the best fitting model in the total group was valid for both boys and girls. Invariance was examined by testing whether trajectories had the same growth parameters (mean, variance and covariance of the intercept, linear slope, and quadratic slope) across gender (Sterba et al. 2007). This was done by simultaneously estimating the models for both boys and girls and using Wald tests to examine differences in growth parameters between boys and girls.

After determining the final model, we examined the association between possible predictors of anxiety and depression symptom trajectories by entering the predictors as auxiliary variables in the growth mixture model. We simultaneously fitted the trajectories and the associations of predictors in a multinomial logistic regression model. As maternal age has been related to increased levels of psychological stress during pregnancy (Aasheim et al. 2012), maternal age was included as a predictor of anxiety and depression symptom trajectories as well.

Then, we examined the associations of anxiety and depression symptom trajectories with psychosocial and school-related outcomes using the three-step approach with an auxiliary model (Asparouhov and Muthén 2014). This approach uses classification probabilities together with the most likely trajectory to correct for classification errors when examining the association of the trajectories with outcomes. For continuous outcomes, differences in the mean and for dichotomous outcomes, the thresholds were compared across trajectories while adjusting for the predictors of these anxiety and depression symptom trajectories. A threshold model was used to test differences for dichotomous outcomes; thresholds are expected values that have to be exceeded for a person to be in a particular category (Muthén and Muthén 2017). Because of the high comorbidity between internalizing and externalizing problems (Thomas and Guskin 2001), all outcomes analyses were repeated while adjusting for externalizing problems at age 10 years. We used global Wald model tests for differences in outcomes across trajectories. When significant, subsequent follow-up Wald model tests were performed to determine which trajectories contributed to the overall difference in outcomes. This hierarchical approach reduces the number of independent tests.

Table 1 gives an overview of the descriptive characteristics of the children and their parents in the study sample and the original Generation R sample. Children from the study sample most often had a Dutch or non-Western ethnicity and the majority of the mothers had an intermediate or higher educational level. Most mothers and fathers reported low levels of parental psychopathology. Mean levels of anxiety and depression symptoms for children in the study sample were 0.15 (SD = 0.18) at 1.5 years, 0.13 (SD = 0.19) at 3 years, 0.18 (SD = 0.23) at 6 years, and 0.17 (SD = 0.20) at 10 years.

The sample characteristics of the study sample (n = 7499) were compared with the characteristics of children and their parents from the original sample who did not participate (n = 2250). Children in the study sample were more likely to be Dutch (χ² = 426.96, p < 0.001) than those not participating. Mothers of children from the study sample were more likely to have a high level of education (χ² = 689.50, p < 0.001) and be older (t = 25.20, p < 0.001). Psychopathology during pregnancy of participating mothers (t = −11.02, p < 0.001) and fathers (t = −5.99, p < 0.001) in the study, was lower compared to parents who were not participating.

Fit indices for models with 2 to 5 number of latent classes are presented in Table 2. The 4 class trajectory model solution with a quadratic slope had the best fit. Although the BIC for the model with 5 classes showed a slightly better fit, the 5 class solution did not lead to another distinct trajectory. Additionally, the LMRT value indicated that the more parsimonious model with 4 trajectory classes provided a better fit than the 5 class model. For the model with 4 latent classes, a model with a quadratic slope showed a better fit (BIC = −14,196.50) compared to a 4 class model with only a linear slope (BIC = −13,797.95). Classification accuracy of the 4 class model was good with an entropy value of 0.85 and average latent class probabilities ranging from 0.74 to 0.95. The same 4 class solution came forward for the sample that included measurements at 3, 6, and 10 years (n = 7130), supporting the inclusion of the measurement of anxiety and depression symptoms at 1.5 years despite lower internal consistency. As a 4 class model with a random linear slope and/or random quadratic slope did not converge for the gender-invariance analyses, a model with a random intercept only was chosen as a final model for subsequent analyses.

Figure 1 gives a graphical presentation of the four anxiety and depression symptom trajectories. We defined the steady low anxiety and depression slope over time (solid line) as the low trajectory. Most children were within this trajectory (82.4%). The decreasing anxiety and depression slope over time (dashed-dotted line) comprised 6.0% of the children (decreasing trajectory). The slope that shows increasing anxiety and depression symptoms up to 6 years and a decreasing slope from 6 to 10 years (dotted line), was defined as the preschool-limited trajectory and comprised 4.2% of the children. The increasing anxiety and depression slope over time (dashed line) comprised 7.4% of the children (increasing trajectory). Anxiety and depression symptoms for children in this trajectory increased to borderline clinical levels at age 10. To examine whether differences between the trajectories were driven by certain items of the CBCL Anxious/Depressed scale, we conducted in-depth analyses of discrepancies in anxiety and depression items between the trajectories. Eight overlapping items that were assessed in both versions of the CBCL showed the same pattern of trajectories as presented in Fig. 1. Thus, differences between trajectories were driven by age and the number of symptoms, instead of specific items.

With regard to differences in trajectories across gender, invariance held for the 4 class model solution. The 4 class model elicited the same type of trajectories for both boys (n = 3785) and girls (n = 3714). The means of the intercepts, linear slopes, quadratic slopes and variance of the intercept of the 4 trajectories were equal across gender (all p > 0.16). We, therefore, concluded that the anxiety and depression symptom trajectories were gender invariant. Small gender differences in the distribution of proportions across categories were observed within the trajectories. The percentage of boys was higher in the low trajectory (82.6%) and pre-school limited trajectory (3.9%) than the percentage of girls (77.3% and 2.9%, respectively). In contrast, the percentage of boys was lower in the increasing (7.3% boys vs. 8.4% girls) and decreasing (6.1% boys vs. 11.4% girls) trajectory compared to girls.

An overview of the associations between predictors and anxiety and depression symptom trajectories is presented in Table 3. Most strongly, children with a Western ethnicity or non-Western ethnicity had higher odds to be in the decreasing trajectory relative to the low trajectory, OR = 4.41, 95% CI [2.21, 8.82], OR = 8.79, 95% CI [5.10, 15.16], than Dutch children. Also, children with a non-Western ethnicity had higher odds to be in the pre-school limited trajectory, OR = 1.64, 95% CI [1.14, 2.34]. In contrast, children with a Western or non-Western ethnicity had lower odds to be in the increasing trajectory relative to the low trajectory, OR = 0.57, 95% CI [0.33, 1.00], OR = 0.68, 95% CI [0.48, 0.97]. Children of mothers with low educational level had higher odds to be in the decreasing, OR = 2.02, 95% CI [1.28, 3.16] as well as in the preschool-limited trajectory, OR = 1.57, 95% CI [1.02, 2.41], compared to children of mothers with high educational level. Maternal psychopathology was positively associated with the increasing trajectory, OR = 1.77, 95% CI [1.50, 2.10], the decreasing trajectory, OR = 1.62, 95% CI [1.36, 1.93], and the preschool-limited trajectory, OR = 1.80, 95% CI [1.55, 2.09]. A similar pattern was found for paternal psychopathology being positively associated with the increasing trajectory, OR = 1.24, 95% CI [1.05, 1.46], decreasing, OR = 1.24, 95% CI [1.07, 1.45], and preschool-limited trajectory, OR = 1.29, 95% CI [1.12, 1.49].

The associations of anxiety and depression symptom trajectories with friendship quality and self-esteem are presented in Table 4. Friendship quality varied between the four trajectories (Wald = 15.63, p = 0.001). Children in the increasing and preschool-limited trajectory had lower friendship quality compared to children in the low trajectory (d = −0.15). After adjusting for externalizing problems at 10 years, self-reported friendship quality did no longer differ across the four trajectories (Wald = 6.52, p = 0.089).

Further, different levels of self-esteem were observed between the trajectories (Wald = 181.51, p < 0.001). Children in the increasing trajectory had a lower self-esteem than children in the other three trajectories (d = −0.80 to −0.94). For self-esteem, the same results were found after adjusting for externalizing problems (Wald = 96.54, p < 0.001).

Table 4 also presents the results of the associations of anxiety and depression symptom trajectories with school performance and school problems. School performance of the children varied across trajectories (Wald = 14.05, p = 0.003). Children with a decreasing trajectory of anxiety and depression symptoms had a higher school performance (d = 0.17) compared to children in the low trajectory. Children within the increasing trajectory and the preschool-limited trajectory had worse school performance compared to children in the decreasing trajectory (d = −0.25, d = −0.26). Although the model test for school performance remained significant after adjusting for externalizing problems (Wald = 12.56, p = 0.006), school performance between children in the pre-school limited trajectory and decreasing trajectory no longer differed (Wald = 2.98, p = 0.084).

School problems of the children varied across trajectories (Wald = 57.15, p < 0.001). Children within the increasing trajectory had a lower threshold for school problems compared to children in the low and decreasing trajectory (d = −0.30, d = −0.36). Children within the preschool-limited trajectory had a lower threshold for school problems compared to children in the low class (d = −0.24). School problems varied across trajectories after adjusting for externalizing problems (Wald = 15.56, p < 0.001). However, the same thresholds for school problems were found for children in the increasing and decreasing trajectory (Wald = 3.76, p = 0.052).

One strength of the current study is that we identified predictors of anxiety and depression symptom trajectories that were measured when parents were expecting a child. Also, this is the first study that relates anxiety and depression symptom trajectories to broader psychosocial and school functioning. Moreover, we related these trajectories to the children’s own reports of their experienced friendship and feelings of self-worth. However, the use of shorter versions of these questionnaire resulted in a lower internal consistency of these outcomes. Therefore, replication of associations between the trajectories and outcomes is warranted in studies that use more extensive measures of children’s psychosocial and school functioning. A final strength is that most previous studies failed to take the uncertainty of trajectory membership into account when the most likely assigned trajectory for each individual is exported to other statistical software outside the growth mixture model. By conducting all analyses within the growth mixture model, our results are not hampered by possible misclassifications of class membership.

One limitation of the current study, and in line with previous studies, is that we only relied on reports of the primary caregiver for measuring child anxiety and depression symptoms. In addition, average scores of the Anxious/Depressed subscale 1.5–5 years and 6–18 years consisted of a different number of items which may hamper the comparability of the average scores across time. As already mentioned, internal consistency of the Anxious/Depressed subscale was low at the age of 1.5 years. Comparable reliability coefficients less than 0.70 have been found for most narrow-band scales of the CBCL in preschool children world-wide (Rescorla et al. 2011). Moreover, the factor structure of the CBCL provides culture-general taxonomic constructs (Ivanova et al. 2010) and the Anxious/Depressed subscale has criterion-related validity in preschool children, despite low internal consistency (Achenbach and Rescorla 2000). In our study, the same trajectory solution appeared when performing the analyses without this first measurement. Future studies may consider using a Factor Mixture Model with item-response information to identify item-level latent variables of the CBCL. Another limitation is the number of children for whom no outcomes measures were available. We should be cautious when generalizing the trajectories as children with a non-Western ethnicity, mothers with a low level of education, and parents who reported on levels of psychopathology during pregnancy were underrepresented in the outcome samples because of drop-out. In particular, missing data on paternal psychopathology was high which weakens our confidence in the results of this predictor. As previous research has shown that study drop-out does not result into different associations between predictors and outcomes in population-based longitudinal follow-up studies (Wolke et al. 2009), we are hesitant to speculate whether these findings would have been different for the whole cohort.

The current study has implications for prevention policies that aim to ward off problems with internalizing, psychosocial, and school functioning for young children in the general population. Our results showed that 7.4% of the children developed anxiety and depression symptoms within the borderline clinical range at age 10 and had lower self-esteem and poorer school-related outcomes. Our findings suggest that selective (i.e. asymptomatic population at higher risk) and indicated (i.e. high risk children with detectable symptoms foreshadowing clinical anxiety or depression) prevention programs are needed instead of universal prevention programs (Stockings et al. 2016; Vázquez-Bourgon et al. 2013; Werner-Seidler et al. 2017). In this study, parental psychopathology served as a risk factor for the development of children’s anxiety and depression symptoms. Although the risk to develop internalizing problems for children of parents with a mood or anxiety disorder is well known (Maciejewski et al. 2018; Weissman et al. 2016), general psychopathology should not be overlooked. Moreover, information about psychopathology experienced by the father should be taken into consideration as well. Prevention programs often target factors that are related to the maintenance of anxiety and depression symptoms instead of modification of risk factors (Lawrence et al. 2017). Selective and indicated prevention programs could benefit from addressing risk factors that predict the development of anxiety and depression symptoms, such as parental psychopathology or indicators of socioeconomic status.