Using the steps of the g-formula as described in the eAppendix 1; https://links.lww.com/EDE/B474, these models are then used to predict antidepressant purchasing () for a scenario in which no hypothetical intervention is performed on employment ( kept at empirically observed values ) termed and is also known as the “natural course scenario,” and to predict antidepressant purchasing in a scenario in which is changed so that no individual is unemployed (), denoted as and termed the “hypothetical intervention scenario.” These expectations represent the total number of person-years with at least one antidepressant purchase, divided by all person-years in the study. The effect of the hypothetical intervention is then estimated as
Each time-varying variable is predicted in sequence annually from t = 1995 to t = 2012 using covariate values from the year before. When using covariate sets 1 and 2, predictions are made for individuals by inserting covariate values into Eq. 1 and adding a random draw from . For binary and multinomial outcomes, predicted probabilities are then used to draw values (0 or 1) from binomial and multinomial distributions, respectively.
When using covariate set 3, predictions are made by inserting covariate values into Eq. 2, adding a random draw from , and drawing from binomial or multinomial distributions where applicable. However, it is important to note that when estimating with covariate set 3, the intercepts are estimated for all individuals, but only those individuals who vary on the time-varying variables contribute to the estimation of the corresponding and coefficients. Hence, when predicting using covariate set 3, the individuals in the sample (i = 1 … n) can be broken down into two groups: those who vary on any of the time-varying variables , referred to here as group Q, and those who do not vary on any of the time-varying variables , group O, with Q ∪ O = N. The sizes of both groups are termed and , respectively, and . If we use the empirical data to make predictions from our estimated models (i.e., the natural course scenario), then the average of a variable in the population at time t can be estimated as:
where hats represent that parameters are estimates from Eq. 2. In Eq. 3, the first term can be seen as containing the between-individual variation for the whole sample, and the second term as containing the within-individual variation for those who vary on any of .
Both between-individual and within-individual differences contribute to the population-averaged estimate of the level of an outcome in a scenario, such as . However, when including individual intercepts in the multivariable estimation, only within-individual differences contribute to the hypothetical intervention effect estimate , because ’s of both scenarios will cancel out. Furthermore, empirically, many individuals only change on a few time-varying covariates and thereby contribute to only a few of the estimated coefficients and in any of the multivariable models. The models with individual-level fixed-effect intercepts therefore have different effective sample sizes for each coefficient (Table 2). In addition, effects can only occur for those individuals whose covariates change as a consequence of the hypothetical intervention, i.e., those for whom . However, not all who belong to the group experience a change on each covariate empirically and therefore do not contribute to the corresponding and coefficients, but our approach transports the estimated and coefficients to all individuals who have .
The models with individual-level fixed-effect intercepts (Eq. 3) were estimated using the “plm” package in R.19 Annotated R code performing the g-formula for sets 2 and 3 is available in eAppendices 2 and 3; https://links.lww.com/EDE/B475 and https://links.lww.com/EDE/B476, respectively.
Subgroup and Sensitivity Analyses
We perform additional analyses to gain further insight into our results. First, we estimate effects by education (secondary and tertiary levels) and sex, because previous research has shown that effects can differ between these groups.5
Second, introducing individual-level intercepts into the modeling procedure reduces the sample that contributes to effect estimation. To gain some insight into the degree to which the change in sample changes our results, we additionally perform the g-formula with covariate set 2 on the population of individuals who changes their employment status, i.e., all individuals who remain constantly employed or unemployed throughout follow-up are removed from this sample.
Third, the estimation of all three covariate sets is also performed with hospitalization because of injuries and accidents (International Classification of Diseases, 10th Revision codes: S00-T79, T90-T98) instead of antidepressant purchasing. Compared with antidepressant purchases, this outcome is less likely to be sensitive to differences in treatment-seeking behavior between social strata.
A sensitivity analysis investigating the strength of potential (unobserved) time-constant confounding can be found in the eAppendix 1; https://links.lww.com/EDE/B474.
Of the study cohort, 18% have lower secondary as their highest education, 48% higher secondary, 24% lower tertiary, and 9% higher tertiary. The registered language of 94% of the cohort is Finnish, 4% Swedish, and 2% other. The percentage of men is 50%.
Antidepressants are purchased by 3% of cohort members in 1995, increasing over time (as individuals age) to 11% in 2012, with 8% on average. In 1995, 16% of all individuals are registered as unemployed, decreasing to 6% at the end of follow-up, with 10% on average. Around 11% of unemployment spells last 5 consecutive years or more. From 1995 to 2012, the percentage employed increases from 71% to 84%, with the discrepancy between the decline in unemployment and the rise in employment caused by changes in other categories (retired or other). In 1995, the average yearly taxable income is €22,090, rising to €35,600 in 2012.
g-Formula Results: No Unemployment Versus Natural Course
Our hypothetical intervention sets all unemployed person-years (10.2% of all person-years) to employed. We use three covariate sets to estimate the effect of this scenario on antidepressant purchases (Table 1). The g-formula’s natural course scenario from the richest measured covariate set (set 2) appears to approximate the empirical data adequately (see eAppendix 1; https://links.lww.com/EDE/B474). Covariate sets that do not include individual-level fixed-effect intercepts (sets 1 and 2) all estimate some population-level reduction in annual person-years with antidepressant purchases (Figure 2). For the set with the most measured covariates (set 2), both time constant (sex, language, and education) and time varying (income, household status, other drug purchases, and previous antidepressant purchases), the hypothetical intervention reduces antidepressant purchasing by 5%. Comparing covariate sets 1 and 2 shows that including time-varying variables attenuates effect estimates.
Covariate set 3 adjusts for unobserved time-constant confounders by including individual intercepts. This covariate set did not find an effect, suggesting that reducing unemployment among those who have experienced unemployment may not reduce antidepressant purchases.
Models additionally including 2- and 3-year lags for all covariates and outcomes were also fitted but did not substantially alter our conclusions. The results of all multivariable models used in our g-formula procedures can be found in eAppendix 4; https://links.lww.com/EDE/B477.
Subgroup and Sensitivity Analyses
The differences found between the covariate sets largely persist in each stratum of our subgroup analysis (Table 3). In the sets without individual intercepts (1 and 2), effect estimates of the hypothetical intervention are stronger for the low educated compared with the high educated and are closer to a null association for women than for men. Reflecting the overall analysis, null associations are found for all subgroups when including individual intercepts (set 3).
Roughly half (25,355) of all individuals in the sample experienced at least one change in their employment status during follow-up. Performing the g-formula with covariate set 2 on this population showed a population-averaged reduction in antidepressant purchasing of 3.8% (95% confidence interval [CI] = 0.6%, 7.0%).
Finally, we find that the effect of our hypothetical (no unemployment) intervention is similar for hospitalization because of injury or accident, as compared with antidepressant purchasing. The effect is slightly stronger than for antidepressant purchasing, but the analysis including individual-level fixed-effect intercepts (set 3) also found no effect despite high precision (Figure 3).
The estimated effect of unemployment on antidepressant purchasing may be biased by time-varying, intermediate, and time-constant confounding. The parametric g-formula is one of the few methods that can account for these sources of bias but has previously required that all relevant confounders be measured. For the first time (to our knowledge), we combine the g-formula with methods to adjust for unobserved time-constant confounding. We estimate the effect of a hypothetical intervention which sets all unemployed person-years (10.2% of all person-years) to employed on antidepressant purchasing. The covariate sets that do not include individual-level fixed-effect intercepts estimate a reduction in the number of person-years with antidepressant purchasing when compared with the natural course. The reduction is estimated to be 5% points (95% CI = 2.5%, 7.4%), population-averaged, in the covariate set that includes time-varying covariates and which most closely follows our assumed directed acyclic graph. However, when including individual-level fixed-effect intercepts to remove unobserved time-constant confounding, the estimate becomes −0.1% points (95% CI = −1.8%, 1.5%).
Strengths and Limitations
Our data constitute an 11% random sample of register data, extracted and anonymized by Statistics Finland following data protection regulations. To limit the influence of potential cohort effects, we choose to follow a closed cohort with a small age range (30–35 years followed until 47–52 years). Furthermore, the age range was chosen because the vast majority of participants have finished education by the age of 30, and because of the increased probability of leaving the labor market for reasons other than unemployment at ages of 50 and above.
Since missingness on covariates is very small (<1%) we do not use missing data imputation methods. The exception is employment status, which has 4.6% missingness. Using a multiple imputation procedure for this variable, including all time-constant and time-varying covariates, does not substantially alter the findings of this study. Furthermore, our natural course scenario of covariate set 2 closely approximates empirical antidepressant purchases, and all other time-varying covariates, which indicates that our models are not grossly misspecified.5,16
Antidepressant purchases do not perfectly indicate the presence of depression. In Finland, only about a third of those with major depression or anxiety disorder received antidepressant treatment.20 However, the specificity of antidepressant purchases is about 90%, and the majority of individuals purchasing psychotropic medication are likely to have a mental health condition.20–23 Furthermore, prior research shows that psychotropic drug purchases are associated with depression and suicide, and react to adverse life events such as divorce and workplace downsizing.20–28 In addition, psychotropic medication purchases require a prescription and are thus based on clinical assessments.
Individuals who become employed may obtain better access to occupational health services, which may increase the probability of a depression diagnosis. It is also possible that higher socioeconomic status individuals and individuals with partners are more likely to receive care.5 If this were the case, the results of the hypothetical intervention of this study (Figure 2) may be negatively biased. However, it is contested whether unmet need for treatment is higher among individuals in low socioeconomic positions, and Finnish population-based studies indicate that income is unrelated to whether or not an individual with a mental health disorder receives treatment.20,29–33
The causal claims of this study rely on three fundamental assumptions: consistency, positivity, and no unobserved confounding.34–36 The consistency assumption requires that there are no different versions of the exposure that have different effects on the outcome. Strictly speaking, different types of employment may have different effects on mental health, which affects the external validity (generalizability) of our study. We reflect on this in later paragraphs. The positivity assumption requires that observed treatment levels vary within confounder strata. In our empirical data, employment status varies within the strata of all measured covariates. Importantly, our study adjusts for unobserved time-constant confounding using individual-level fixed-effect intercepts and observed time-varying confounding. To the best of our knowledge, we are the first to do so using the parametric g-formula. It is possible that by restricting our analysis (through the use of individual intercepts) to the group that has experienced changes in employment status, the influence of time-varying confounders becomes more prominent. Thereby, adding individual-level intercepts may compromise internal validity, potentially biasing effect estimates even if we restrict our inference to the population of changers.37 We include a large set of potentially important time-varying confounders; however, it remains possible that some time-varying confounding is unobserved. For example, unobserved health shocks may influence both employment status and antidepressant purchasing. Although we are not able to adjust for this, we believe this form of bias in our study is small, since the hypothetical intervention effect is also small.
Unobserved Time-constant Confounding
When including the individual intercepts into the modeling procedure, the hypothetical intervention effect becomes null. As indicated by our subgroup analysis, the null effect is also found across educational and sex strata. This could indicate the presence of strong time-constant factors that affect both unemployment and antidepressant purchasing. The sensitivity analysis that performs the g-formula without individual intercepts on the subsample of those who change their employment status showed only a 1.2% point lower reduction in antidepressant purchasing compared with the ordinary set 2 estimate. This may indicate that only a small part of the change in the estimate between the g-formulas with and without individual intercepts is due to the effect being estimated for a qualitatively different population. However, this is not certain; due to having multiple time-varying covariates, most individuals who do not vary on employment still vary on other variables, and thereby contribute to the coefficients of those variables. These contributions are now also lost from this sensitivity analysis’ estimate.
In the sensitivity analysis where we replace antidepressant purchasing with injuries and accidents as the outcome variable, conclusions reflect those of the main analysis. Not adjusting for individual intercepts results in a hypothetical intervention effect that reduces injuries and accidents, but adjusting for individual-level intercepts results in a null estimate. This finding is important because unemployment and hospitalization due to injuries and accidents are likely to be partly affected by time-constant factors similar to those that also affect unemployment and depression (or antidepressant purchasing), such as early problem behavior, childhood adversity, and personality traits such as neuroticism and lack of conscientiousness.38–40 At the same time hospitalization due to injuries and accidents as an outcome is less sensitive to selection through care-seeking behavior, compared with antidepressant purchasing. Therefore, these findings imply that differences in care-seeking do not explain our results.
The Unemployment Effect
The association between employment status and antidepressant purchasing may be explained by time-constant factors that make individuals both more prone to be unemployed and to purchase antidepressants, though it may also be explained by a change in meaning and generalizability of the employment effect in the model with individual-level intercepts. Time-constant factors may include poor childhood physical and mental health, adverse childhood socioeconomic conditions or personality traits such as neuroticism and lack of conscientiousness.38–40 The findings of our study corroborate those of other studies that found no effects of unemployment on mental health when including individual intercepts to adjust for unobserved time-constant confounding.41,42 A third study using individual intercepts found only weak estimated effects of unemployment on health, and large selection into unemployment and poor health.43 However, it is important to note that adding individual intercepts to the modeling procedure also changes the meaning of the coefficients in the multivariable models, and thereby the results of a g-formula procedure that uses these coefficients. When including individual-level fixed-effect intercepts, the only observations that contribute to the estimation of relevant coefficients are those where the exposure varies over time.17 In our study, on average 10.2% of person-years were unemployed, and roughly half of all individuals were continuously employed throughout the study. Individuals in the ages of 30–55 years who spent a year unemployed are likely to be qualitatively different from individuals who have not spent any time unemployed (as measured annually). The type of employment that they can obtain may also differ qualitatively. Individuals who switch regularly between unemployment and employment are more likely to experience precarious employment. Various studies have found that those in precarious employment, such as temporary employment or employment with weak social protections, and those transferring to such jobs have higher mental health risks than those in nonprecarious employment.44–47 This is not to imply that the estimates of the covariate set with individual intercepts should be disregarded. Some employment interventions may result in less fulfilling and more precarious forms of employment and may therefore be better approximated by the g-formula with individual intercepts.
To our knowledge, this is the first study–of any exposure or outcome–to combine individual-level fixed-effect estimation with the parametric g-formula in order to adjust for unobserved time-constant confounding. We use this method to study the impact of unemployment on antidepressant purchasing. At the population level, without fixed-effect intercepts, we find a substantial reduction in antidepressant purchasing in a scenario where everyone is employed. However, when we include individual-level fixed-effect intercepts, we estimate a null effect. We argue that this null effect arises not only due to the adjustment for unobserved time-constant confounding but also due to a change in the meaning and generalizability of the effect of employment in models with individual-level fixed-effect intercepts.
1. McLaughlin KA. The public health impact of major depression: a call for interdisciplinary prevention efforts. Prev Sci. 2011;12:361–371.
2. Lépine JP, Briley M. The increasing burden of depression. Neuropsychiatr Dis Treat. 2011;7(suppl 1):3–7.
3. World Health Organization. The Global Burden of Disease: 2004 Update. 2004. Available at: http://www.who.int/entity/healthinfo/global_burden_disease/GBD_report_2004update_full.pdf
. Accessed 1 August 2018.
4. Mathers CD, Loncar D. Projections of global mortality and burden of disease from 2002 to 2030. PLoS Med. 2006;3:e442.
5. Bijlsma MJ, Tarkiainen L, Myrskylä M, Martikainen P. Unemployment
and subsequent depression: a mediation analysis using the parametric G-formula. Soc Sci Med. 2017;194:142–150.
6. Jefferis BJ, Nazareth I, Marston L, et al. Associations between unemployment
and major depressive disorder: evidence from an international, prospective study (the predict cohort). Soc Sci Med. 2011;73:1627–1634.
7. Paul KI, Moser K. Unemployment
impairs mental health: meta-analyses. J Vocat Behav. 2009;74:264–282.
8. Buckley JP, Keil AP, McGrath LJ, Edwards JK. Evolving methods for inference in the presence of healthy worker survivor bias. Epidemiology. 2015;26:204–212.
9. Wagenaar AF, Kompier MA, Houtman IL, van den Bossche SN, Taris TW. Employment contracts and health selection: unhealthy employees out and healthy employees in? J Occup Environ Med. 2012;54:1192–1200.
10. Steele F, French R, Bartley M. Adjusting for selection bias in longitudinal analyses using simultaneous equations modeling: the relationship between employment transitions and mental health. Epidemiology. 2013;24:703–711.
11. De Stavola BL, Daniel RM, Ploubidis GB, Micali N. Mediation analysis with intermediate confounding: structural equation modeling viewed through the causal inference lens. Am J Epidemiol. 2015;181:64–80.
12. Vanderweele TJ, Vansteelandt S, Robins JM. Effect decomposition in the presence of an exposure-induced mediator-outcome confounder. Epidemiology. 2014;25:300–306.
13. Jain P, Danaei G, Robins JM, Manson JE, Hernán MA. Smoking cessation and long-term weight gain in the Framingham Heart Study: an application of the parametric g-formula for a continuous outcome. Eur J Epidemiol. 2016;31:1223–1229.
14. Keil AP, Edwards JK, Richardson DB, Naimi AI, Cole SR. The parametric g-formula for time-to-event data: intuition and a worked example. Epidemiology. 2014;25:889–897.
15. Robins JM. A new approach to causal inference in mortality studies with a sustained exposure period—application to control of the healthy worker survivor effect. Math Mod. 1986;7:1393–1512.
16. Vangen-Lønne AM, Ueda P, Gulayin P, Wilsgaard T, Mathiesen EB, Danaei G. Hypothetical interventions to prevent stroke: an application of the parametric g-formula to a healthy middle-aged population. Eur J Epidemiol. 2018;33:557–566.
17. Gunasekara FI, Richardson K, Carter K, Blakely T. Fixed effects analysis of repeated measures data. Int J Epidemiol. 2014;43:264–269.
18. Lin SH, Young J, Logan R, Tchetgen EJ, VanderWeele TJ. Parametric mediational g-formula approach to mediation analysis with time-varying exposures, mediators, and confounders. Epidemiology. 2017;28:266–274.
19. Croissant Y, Millo G. “Panel Data Econometrics in R: The plm Package.” Journal of Statistical Software 27;1–43. 2008). Available at: http://doi.org/10.18637/jss.v027.i02
20. Sihvo S, Isometsä E, Kiviruusu O, et al. Antidepressant utilisation patterns and determinants of short-term and non-psychiatric use in the Finnish general adult population. J Affect Disord. 2008;110:94–105.
21. Hämäläinen J, Isometsä E, Sihvo S, Kiviruusu O, Pirkola S, Lönnqvist J. Treatment of major depressive disorder in the Finnish general population. Depress Anxiety. 2009;26:1049–1059.
22. Hämäläinen J, Isometsä E, Sihvo S, Pirkola S, Kiviruusu O. Use of health services for major depressive and anxiety disorders in Finland. Depress Anxiety. 2008;25:27–37.
23. Piek E, van der Meer K, Hoogendijk WJ, Penninx BW, Nolen WA. Most antidepressant use in primary care is justified; results of the Netherlands Study of Depression and Anxiety. PLoS One. 2011;6:e14784.
24. Kivimäki M, Tabák AG, Lawlor DA, et al. Antidepressant use before and after the diagnosis of type 2 diabetes: a longitudinal modeling study. Diabetes Care. 2010;33:1471–1476.
25. Metsä-Simola N, Martikainen P. Divorce and changes in the prevalence of psychotropic medication use: a register-based longitudinal study among middle-aged Finns. Soc Sci Med. 2013;94:71–80.
26. Magnusson Hanson LL, Westerlund H, Chungkham HS, Vahtera J, Sverke M, Alexanderson K. Purchases of prescription antidepressants in the Swedish population in relation to major workplace downsizing. Epidemiology. 2016;27:257–264.
27. Thielen K, Nygaard E, Andersen I, et al. Misclassification and the use of register-based indicators for depression. Acta Psychiatr Scand. 2009;119:312–319.
28. Tiihonen J, Lönnqvist J, Wahlbeck K, Klaukka T, Tanskanen A, Haukka J. Antidepressants and the risk of suicide, attempted suicide, and overall mortality in a nationwide cohort. Arch Gen Psychiatry. 2006;63:1358–1367.
29. Hämäläinen J, Isometsä E, Laukkala T, et al. Use of health services for major depressive episode in Finland. J Affect Disord. 2004;79:105–112.
30. Kivimäki M, Honkonen T, Wahlbeck K, et al. Organisational downsizing and increased use of psychotropic drugs among employees who remain in employment. J Epidemiol Community Health. 2007;61:154–158.
31. Starkes JM, Poulin CC, Kisely SR. Unmet need for the treatment of depression in Atlantic Canada. Can J Psychiatry. 2005;50:580–590.
32. von Soest T, Bramness JG, Pedersen W, Wichstrøm L. The relationship between socio-economic status and antidepressant prescription: a longitudinal survey and register study of young adults. Epidemiol Psychiatr Sci. 2012;21:87–95.
33. Hansen DG, Søndergaard J, Vach W, et al. Socio-economic inequalities in first-time use of antidepressants: a population-based study. Eur J Clin Pharmacol. 2004;60:51–55.
34. Rehkopf DH, Glymour MM, Osypuk TL. The consistency assumption for causal inference in social epidemiology: when a rose is not a rose. Curr Epidemiol Rep. 2016;3:63–71.
35. Petersen ML, Porter KE, Gruber S, Wang Y, van der Laan MJ. Diagnosing and responding to violations in the positivity assumption. Stat Methods Med Res. 2012;21:31–54.
36. Greenland S, Robins JM. Identifiability, exchangeability and confounding revisited. Epidemiol Perspect Innov. 2009;6:4.
37. Frisell T, Öberg S, Kuja-Halkola R, Sjölander A. Sibling comparison designs: bias from non-shared confounders and measurement error. Epidemiology. 2012;23:713–720.
38. Almlund M, Duckworth AL, Heckman JJ, Kautz TD. Personality Psychology and Economics [Internet]. National Bureau of Economic Research. Report No.: 16822. 2011. Available at: http://www.nber.org/papers/w16822
. Accessed 1 August 2018.
39. Klein DN, Kotov R, Bufferd SJ. Personality and depression: explanatory models and review of the evidence. Annu Rev Clin Psychol. 2011;7:269–295.
40. Goodman A, Joyce R, Smith JP. The long shadow cast by childhood physical and mental problems on adult life. Proc Natl Acad Sci U S A. 2011;108:6032–6037.
41. Schmitz H. Why are the unemployed in worse health? The causal effect of unemployment
on health. Labour Econ. 2011;18:71–78.
42. Böckerman P, Ilmakunnas P. Unemployment
and self-assessed health: evidence from panel data. Health Econ. 2009;18:161–179.
43. Tøge AG, Blekesaune M. Unemployment
transitions and self-rated health in Europe: a longitudinal analysis of EU-SILC from 2008 to 2011. Soc Sci Med. 2015;143:171–178.
44. Virtanen M, Kivimäki M, Joensuu M, Virtanen P, Elovainio M, Vahtera J. Temporary employment and health: a review. Int J Epidemiol. 2005;34:610–622.
45. Quesnel-Vallée A, DeHaney S, Ciampi A. Temporary work and depressive symptoms: a propensity score analysis. Soc Sci Med. 2010;70:1982–1987.
46. Rodriguez E. Marginal employment and health in Britain and Germany: does unstable employment predict health? Soc Sci Med. 2002;55:963–979.
47. Jang SY, Jang SI, Bae HC, Shin J, Park EC. Precarious employment and new-onset severe depressive symptoms: a population-based prospective study in South Korea. Scand J Work Environ Health. 2015;41:329–337.