This is a follow-up report to The differential impact of major life events on cognitive and affective wellbeing (Kettlewell et al., 2020).


1 Aims

This document describes the heterogeneity of mental health (MHi-5) trajectories around life events.

The objective is to identify distinct trajectories or response-types in the average response to life-events, and challenge assumptions that an “average” response exists. If distinct trajectories exist then the response to major life events is heterogenous in the population, and inferring individuals are resilient to major life events on the basis of the marginal response in the population may be committing the ecological fallacy and not be justified. Some preliminary investigation of factors predicting the subgroups associated with distinct trajectories was also explored.


1.1 Background

Kettlewell et al. (2020) reported the average response to the worst life events such as bankruptcy or getting fired was \(U\)-shaped, and concluded that most people (on average) are resilient to negative life events. Here I test that assumption by estimating whether latent classes of trajectories exist within the “average” trajectory using growth mixture modeling (GMM). It is possible that most people are not resilient, and instead the average trajectory contains multiple trajectories within it. These sub trajectories may or may not exhibit the same trend as the average but nevertheless combine to produce an average U-shaped trajectory.

Since then two studies of the HILDA survey, using different methods, have provided evidence for heterogeneous responses to major life events. Quantile regression of life satisfaction scores showed differences between trajectories at different levels (quantiles) of life satisfaction (O’Leary, Li, Gupta, & Blackaby (2020)). And finite mixture modelling of MHI-5 scores revealed the three latent trajectories (Etilé, Frijters, Johnston, & Shields (2021)), shown below.

Etile et al 2021 Figure 3

Model estimated trajectories in MHI-5 scores in response to a life event (at Year 0)
Model estimated trajectories in MHI-5 scores in response to a life event (at Year 0)

 

In both studies, the differences between trajectories can best be described as a difference in the degree of the same \(U\)-shaped function, rather than representing distinctly different shaped trajectories. Such results support the original inference that individuals are resilient to major life events (e.g., Etilé et al., 2021).


1.2 Growth Mixture Models (GMM)

Standard growth models (Latent Growth Curve Analysis or LGCA) estimate a single trajectory that averages the individual trajectories of all participants in the sample. Individual differences are captured by estimating a random coefficient representing the variability around the slope and intercept. Categorical or continuous features can be used to predict individual differences in slope or intercept (e.g., risk factors). LGCA is only useful when all individuals are expected to change in the same direction over time (the degree of change can vary).

GMMs extend LGCA by allowing for multiple latent classes, each with its own growth trajectory. It estimates class-specific means and variances for intercepts and slopes. Non-linear slopes can be included using polynomial terms (quadrative and cubic trends), or by relaxing the functional form entirely using splines. Individuals are probabilistically assigned to classes based on their data. GMM captures heterogeneity in the marginal trajectory, and allows for within-class variability so individuals within a class can still differ somewhat.


1.2.1 Notes of caution:

Multiple authors have pointed out that most (all?) widely-used methods to identify latent trajectories in psychology to date are typically not assessed against a ground truth, or provide reasonable controls for overfitting. In other words, such methods will tend to “discover” latent trajectories where none exist. We should be very careful when interpreting results of methods which rely solely on model fit statistics (e.g., \(R^2\), AIC) or impose arbitrary cut-offs (e.g., \(p < .05\)), especially when using large datasets.

  • Erosheva, E. A., Matsueda, R. L., & Telesca, D. (2014). Breaking bad: two decades of life-course data analysis in criminology, developmental psychology, and beyond. Annual Review of Statistics and Its Application, 1, 301-332.
  • Sonya K. Sterba , Ruth E. Baldasaro & Daniel J. Bauer (2012) Factors affecting the adequacy and preferability of semiparametric groups-based approximations of continuous growth trajectories, Multivariate Behavioral Research, 47:4, 590-634.
  • Bauer DJ. 2007. Observations on the use of growth mixture models in psychological research. Multivar. Behav.Res. 42(4):757–86.
  • Bauer DJ, Curran PJ. 2004. The integration of continuous and discrete latent variable models: Fourth level potential problems and promising opportunities. Psychol. Methods 9(1):3–29.
  • Hoeksma, J. B., & Kelderman, H. (2006). On growth curves and mixture models. Infant and Child Development: An International Journal of Research and Practice, 15(6), 627-634.
  • Bauer DJ, Curran PJ. 2003a. Distributional assumptions of growth mixture models: Implications for overextraction of latent trajectory classes. Psychol. Methods 8(3):338–63.


2 Methods


2.1 Subjective wellbeing

Two measures of subjective wellbeing were used in the present report: MHI-5 scores representing mental health (or affective wellbeing), and life-satisfaction scores representing happiness (or cognitive wellbeing).


2.1.1 MHI-5 scores

To be added…


2.1.2 Life-satisfaction

To be added…


2.2 Growth Mixture Model

To identify distinct patterns of change over time, we employed growth mixture modelling using the lcmm package in R (Proust-Lima, Philipps, & Liquet, 2017). This approach estimates growth trajectories while accounting for individual variability and class-specific parameters. Models with two and three latent classes were fit to the data, allowing for the identification of subgroups characterized by differing temporal profiles. Class membership was probabilistically assigned based on posterior probabilities, and model comparisons were guided by fit indices and classification diagnostics.

Each model included polynomial terms (i.e., linear, quadratic and cubic trends) to distinguish the trajectory of each latent class. Polynomial terms are a common choice due to their flexibility and ease of implementation, however they impose some restrictions on the functional form (e.g., they do not allow for extrapolation or prediction beyond the observed data range). In our case the restricted functional form is not a disadvantage since our aim is descriptive, and the restriction will facilitate convergence during model estimation. We selected a data range of ±3 years around the event, which was deemed sufficient to capture the variation due to anticipation and reaction to the event. Years beyond this range are likely to reflect other events or factors (rendering any extrapolation unreliable and of little value, see Kettlewell et al 2020 and also Fritjers et al). We also restricted the analysis to 3rd order terms because higher order polynomial terms result in parameters that are difficult to interpret, and are highly contingent on data handling features such as centering choices (Dalal & Zickar, 2012).

We compared models with the inclusion of each additional trajectory \(j\) by examining three key fit indices: AIC computed as -2L+2P (lower is better); BIC computed as -2L+P log(N) (lower is better); and Entropy computed as \(1+\sum[\pi_{ig} \times log(\pi_{ig})]/(N \times log(G))\) where \(\pi_{ig}\) is the posterior probability that subject \(i\) belongs to class \(g\) (the closer to one the better).

To evaluate the success of classification in the winning model, we examined the distribution of posterior probabilities for each latent class. Specifically, we quantified the mean posterior probabilities (MPP) of each class to indicate the certainty of class membership. A mean posterior probability greater than 70% was considered an acceptable level of certainty. We also examined the percentage of individuals whose probability of membership exceeded a similar arbitrary threshold (e.g., > 0.7) for each class, which indicates similar patterns of change within trajectories. This percentage provided an index of separation or discrimination between classes, and both the MPP and individual percentage statistics help quantify the discriminative ability and internal reliability of the classification.


3 Results


3.1 Bankruptcy

A major financial loss or bankruptcy (lefnw) was the most severe life event reported in HILDA according to the depth of the U-shaped function (Kettlewell et al., 2020).

Initially a linear multilevel model (LMM) of MHI-5 scores was fit using the lcmm package, with fixed and random effects for linear, quadratic and cubic polynomial terms and no mixtures, so this model assumed a homogeneous response. This model was used as a baseline comparison for subsequent growth mixture models (GMM) with 2, 3 and 4 classes.


Figure 1. LMM results: Predictions (points) and observed means (lines ±95%CI shaded)

*Figure legend*: Predicted (points) and mean observed MHI-5 score (lines and ±95%CI shaded) before and after bankruptcy at time 0

Figure legend: Predicted (points) and mean observed MHI-5 score (lines and ±95%CI shaded) before and after bankruptcy at time 0

The MHI-5 trend around bankruptcy revealed the familiar U-shaped form previously reported by Kettlewell et al. (2020), which was interpreted as resilience (see also Etilé et al., 2021). Note the model predictions (points) did not entirely fall inside the expected 95%CI (shaded region) of the observed data, likely due to the restrictive functional form of the polynomial terms.


3.1.1 Model comparison

Each GMM successfully converged according to the default thresholds in lcmm. The three key fit statistics across 1- to 4-class models were compared below.

Figure 2. GMM fit comparisons

*Figure legend*: Fit statistics over models spanning 4 latent classes. Lower scores better for AIC and BIC. Scores closer to unity (1) better for entropy

Figure legend: Fit statistics over models spanning 4 latent classes. Lower scores better for AIC and BIC. Scores closer to unity (1) better for entropy

The results of the model comparison indicated the single class LLM was preferred by AIC and entropy, while a 3-class GMM was preferred by BIC.


3.1.2 Class discrimination

The discriminative performance of the 3-class GMM was assessed by examining the distribution of the posterior probability of class memberships, as shown below.

Table 1. Posterior probabilities of classification (3-class GMM)

##  
## Posterior classification: 
##   class1  class2 class3
## N 355.00 2282.00  81.00
## %  13.06   83.96   2.98
##  
## Posterior classification table: 
##      --> mean of posterior probabilities in each class 
##         prob1  prob2  prob3
## class1 0.7862 0.1657 0.0481
## class2 0.0531 0.8963 0.0506
## class3 0.0527 0.2142 0.7331
##  
## Posterior probabilities above a threshold (%): 
##          class1 class2 class3
## prob>0.7  66.76  89.40  53.09
## prob>0.8  54.65  81.81  35.80
## prob>0.9  35.49  67.22  24.69
## 

The three classes distinguished by the GMM were dominated by a single class (class2 was 83.96% of cases). The mean posterior probablities of each class were greater than 0.7 and acceptable (lowest MPP was 0.7331). However the proportion of individuals assigned to a class with prob>0.7 was low (e.g., class3 53.09%). This indicates a substantial proportion of individuals were not strongly classified into class1 or class3, perhaps due to overlapping trajectories or too many classes in the model.

The predicted vs observed trajectories of the two and three class models are shown below to further judge class discrimination.


Figure 3. GMM results: Predictions (points) and observed means (lines ±95%CI shaded) of each class

*Figure legend*: Weighted subject-specific predicted and mean MHI-5 scores of class 1 (Black), class 2 (Red), and class 3 (Green)*Figure legend*: Weighted subject-specific predicted and mean MHI-5 scores of class 1 (Black), class 2 (Red), and class 3 (Green)

Figure legend: Weighted subject-specific predicted and mean MHI-5 scores of class 1 (Black), class 2 (Red), and class 3 (Green)

In both the 2-class and 3-class GMMs, distinct trajectories of mental health were revealed over the time period. In the 3-class GMM, the 3rd class (green) did not display a \(U\)-shaped recovery after the event at time 0, and indicating that resilience may not be a universally shared response to bankruptcy in the population.


3.1.3 Predicting classes

The class discrimination performance indicated the response to bankruptcy may be heterogeneous, however there was only weak evidence the individuals in the smaller two classes were well distinguished. Adding class predictors to the model may help distinguish the smaller classes and improve the evidence for true heterogeneity of the response.

We fit a 3-class growth mixture model, which included 3 latent classes as well as predictors for class membership. As the GMM assumes class membership does not vary over time we selected stationary variables as class predictors, including sex, chronic illness and highest education level. We also calculated the average socioecononimc index (SEIFA), average age and average calendar years over the period observed for each person as stationary variables. Finally we included dummy covariates for other events as fixed effects to explain variability due to concomitant life events.

Table 2. Posterior probabilities of classification with class predictors

##  
## Posterior classification: 
##    class1 class2 class3
## N 1881.00 320.00 517.00
## %   69.21  11.77  19.02
##  
## Posterior classification table: 
##      --> mean of posterior probabilities in each class 
##         prob1  prob2  prob3
## class1 0.9030 0.0431 0.0539
## class2 0.1402 0.7498 0.1100
## class3 0.1248 0.0900 0.7852
##  
## Posterior probabilities above a threshold (%): 
##          class1 class2 class3
## prob>0.7  88.46  58.75  66.15
## prob>0.8  80.81  44.06  51.84
## prob>0.9  69.86  27.81  34.82
## 

The 3-class model with class predictors produced similar class discrimination as above (3.1.2), however the size of the smallest classes increased so the largest class represented less than 7 in 10 people (69.21%). The individual posterior probabilities in the smallest class were relatively improved (53.09% increased to 58.75%).


Figure 4. Predicted and observed trajectories for the 3-class GMM with predictors

*Legend*: Three latent class trajectories of the mean predicted (points) and mean observed response (lines and ±95%CI shaded) to bankruptcy. Black: class 1; Red: class 2; Green: class 3

Legend: Three latent class trajectories of the mean predicted (points) and mean observed response (lines and ±95%CI shaded) to bankruptcy. Black: class 1; Red: class 2; Green: class 3

The trajectories of all three classes were well-distinguished and the model predictions (points) fell within the expected 95%CI (shaded region) of the observed data. While there were differences in the overall level of mental health between classes, as well as differences in the extent of recovery, all three classes showed a common \(U\)-shaped functional form.


Table 3. Fixed estimates of class-membership

(the class of reference is the last class) 

                     coef      Se    Wald p-value
intercept class1  1.03173 2.99766   0.344 0.73071
intercept class2 -0.32457 0.23466  -1.383 0.16661
chronic class1   -1.14240 0.09552 -11.960 0.00000***
chronic class2   -0.10783 0.11943  -0.903 0.36656
edu.L class1      0.42316 2.66311   0.159 0.87375
edu.L class2     -0.22727 0.26913  -0.844 0.39842
edu.Q class1      0.75848 5.98532   0.127 0.89916
edu.Q class2     -0.20923 0.33188  -0.630 0.52841
edu.C class1     -1.04017 8.03291  -0.129 0.89697
edu.C class2      0.01808 0.39547   0.046 0.96353
SEIFA class1      0.23881 0.08127   2.938 0.00330**
SEIFA class2      0.04090 0.10625   0.385 0.70029
age class1        0.89698 0.09972   8.995 0.00000***
age class2       -0.12844 0.12265  -1.047 0.29497
period class1    -0.39324 0.08226  -4.780 0.00000***
period class2    -0.22761 0.11139  -2.043 0.04101*

The standardized class-membership coefficients revealed some significant differences between classes (expressed as the log odds of membership relative to class 3). The predictor with the largest standardized coefficient was chronic illness which significantly reduced the log odds of membership in class 1 (log odds -1.14 [0.09], p < .001). Conversely the log odds of membership in class 1 significantly increased with higher socioeconomic index (p < .01) and age (p < .001), and decreased with calendar year (p < .001). The log odds of membership in class 2 also decreased with calendar year (p < .05).

Thus relative to class 3 (green), class 1 (black) were less likely to be chronically ill, older and enjoy a higher socioeconomic status. Class 2 was slightly more likely to be observed at an earlier time point.


3.1.4 Life-satisfaction

Figure 5. Model fit comparisons

*Figure legend*: Fit statistics over models spanning 4 latent classes. Lower scores better for AIC and BIC. Scores closer to unity (1) better for entropy

Figure legend: Fit statistics over models spanning 4 latent classes. Lower scores better for AIC and BIC. Scores closer to unity (1) better for entropy


Table 4. Posterior probabilities of classification

##  
## Posterior classification: 
##   class1  class2
## N 267.00 2748.00
## %   8.86   91.14
##  
## Posterior classification table: 
##      --> mean of posterior probabilities in each class 
##         prob1  prob2
## class1 0.8520 0.1480
## class2 0.0225 0.9775
##  
## Posterior probabilities above a threshold (%): 
##          class1 class2
## prob>0.7  77.53  98.07
## prob>0.8  68.16  96.29
## prob>0.9  55.43  93.16
## 

The estimated classes are dominated by class2 (91.14%), while the other class comprises less than 10 percent of the respondents. The mean posterior probabilities in each class are greater than 0.8 (supragood), and a substantial majority of individual’s posterior probabilities in each class are above 0.7 (lowest class1 77.53%).


Figure 6. Predicted and observed trajectories for 1- and 2-class models

*Figure legend*: Latent class trajectories of the mean predicted (points) and mean observed response (lines & ±95%CI shaded) to bankruptcy. Black: class 1; Red: class 2*Figure legend*: Latent class trajectories of the mean predicted (points) and mean observed response (lines & ±95%CI shaded) to bankruptcy. Black: class 1; Red: class 2

Figure legend: Latent class trajectories of the mean predicted (points) and mean observed response (lines & ±95%CI shaded) to bankruptcy. Black: class 1; Red: class 2

The impact of bankruptcy on life satisfaction among the majority class (red 91.14%) was relatively diminished, however both classes showed a \(U\)-shaped functional form. Many of the predicted estimates fall outside the expected 95%CI range of the observed mean in each class suggesting the model was not precise, which may be improved by adding class predictors (see above).


3.2 Divorce

Separation or divorce (lesep) was one of the top three most severe life events reported in HILDA (Kettlewell et al., 2020).

A baseline linear multilevel model (LMM) of MHI-5 scores was fit, with fixed and random effects for linear, quadratic and cubic polynomial terms and no mixtures. This was compared against subsequent growth mixture models (GMM) with 2, 3 and 4 classes.

Sex was included as a class predictor in each GMM.


The default convergence threshold was reached for 2- and 3-class GMMs only, so the 4-class model was excluded from the results.

Figure 7. Model fit comparisons

The 2-class GMM had the lowest BIC while the entropy of the 1-class model (LMM) was nearest unity.


Table 5. Posterior probabilities of classification

##  
## Posterior classification: 
##   class1 class2
## N 2884.0  529.0
## %   84.5   15.5
##  
## Posterior classification table: 
##      --> mean of posterior probabilities in each class 
##         prob1  prob2
## class1 0.9424 0.0576
## class2 0.1819 0.8181
##  
## Posterior probabilities above a threshold (%): 
##          class1 class2
## prob>0.7  94.00  72.40
## prob>0.8  89.18  58.41
## prob>0.9  82.04  42.53
## 


Figure 8. Predicted and observed trajectories for 1-, and 2-class models

*Figure legend*: Latent class trajectories of the mean predicted (points) and mean observed response (lines & ±95%CI shaded) to divorce. Black: class 1; Red: class 2

Figure legend: Latent class trajectories of the mean predicted (points) and mean observed response (lines & ±95%CI shaded) to divorce. Black: class 1; Red: class 2

The trajectory of both classes improved after divorce but only class 1 (black) was \(U\)-shaped. The trajectory of class 2 suggests the mental health of individuals in that class improved over the period, perhaps in anticipation of the divorce. The fixed effect of sex (log odds ratio of 0.717±0.12 SE, \(p < .001\)) indicated males were twice as likely (OR 2.05) to be members of class 1 than class 2, which is partially consistent with other research using HILDA showing that unhappy women tend to anticipate divorce and are happier after divorce.


4 Conclusions

This report establishes some preliminary evidence of heterogeneity in response to life events, however there were very few examples in the small set of life events studied (so far) of deviations from a \(U\)-shaped response. That is, there is little to suggest that some level of resilience is not commonly shared in the population. Nevertheless studying more events may reveal additional latent trajectories which can test this assumption.



5 Appendix


5.0.1 Latent Class Growth Analysis (LCGA) (Nagin, 1999)

https://www.tqmp.org/RegularArticles/vol05-1/p011/p011.pdf

LCGA is a semi-parametric technique used to identify distinct subgroups of individuals following a similar pattern of change over time. Non-linear trajectories can be modelled by including parametric terms (e.g., linear, quadratic, cubic). Heterogeneity is summarised by a finite set of unique polynomial functions each corresponding to a distinct trajectory. Instead of modelling individual differences with random coefficients, LCGA fixes the slope and intercept to equality across individuals within a trajectory. This means the variance in growth parameters are fixed to zero, and there is no additional variation within class - just random error. See Nagin 2005 for details.

See R package lcmm and specify no random effect (i.e., random=~-1 in lcmm or hlme)

Other forms of growth curve modelling also use fixed functional forms (linear, quadratic and cubic trends) to identify trajectories, including the standard growth curve models or group-based models (Nagin 1999, 2005, 2010).


References

Dalal, D. K., & Zickar, M. J. (2012). Some common myths about centering predictor variables in moderated multiple regression and polynomial regression. Organizational Research Methods, 15(3), 339–362. https://doi.org/10.1177/1094428111430540
Etilé, F., Frijters, P., Johnston, D. W., & Shields, M. A. (2021). Measuring resilience to major life events. Journal of Economic Behavior & Organization, 191, 598–619. https://doi.org/10.1016/j.jebo.2021.09.004
Kettlewell, N., Morris, R. W., Ho, N., Cobb-Clark, D. A., Cripps, S., & Glozier, N. (2020). The differential impact of major life events on cognitive and affective wellbeing. SSM-Population Health, 10, 100533. https://doi.org/doi.org/10.1016/j.ssmph.2019.100533
O’Leary, N., Li, I. W., Gupta, P., & Blackaby, D. (2020). Wellbeing trajectories around life events in australia. Economic Modelling, 93, 499–509. https://doi.org/10.1016/j.econmod.2020.08.021
Proust-Lima, C., Philipps, V., & Liquet, B. (2017). Estimation of extended mixed models using latent classes and latent processes: The r package lcmm. Journal of Statistical Software, 78, 1–56. https://doi.org/10.18637/jss.v078.i02