Evaluating a Method to Estimate Mediation Effects With Discrete-Time Survival Outcomes

The utility of evaluating mediation effects spans across research domains. The model facilitates investigation of underlying mechanisms of event timing and, as such, has the potential to help strengthen etiological research and inform intervention work that incorporates the evaluation of mediating variables. In order for the analyses to be maximally useful however, it is critical to employ methodology appropriate for the data under investigation. The purpose of this paper is to evaluate a regression-based approach to estimating mediation effects with discrete-time survival outcomes. We empirically evaluate the performance of the discrete-time survival mediation model in a statistical simulation study, and demonstrate that results are functionally equivalent to estimates garnered from a potential-outcomes framework. Simulation results indicate that parameter estimates of mediation in the model were statistically accurate and precise across the range of examined conditions. Type 1 error rates were also tolerable in the conditions studied. Adequate power to detect effects in the model, with binary X and continuous M variables, required effect sizes of the mediation paths to be medium or large. Possible extensions of the model are also considered.

Introduction

Mediation hypotheses have become increasingly pervasive across a wide range of research areas in recent years, as investigators use statistical mediation analysis to elucidate developmental pathways of behavior (e.g., Selig and Preacher, 2009; George et al., 2014; Gonzalez and MacKinnon, 2018), inform etiological underpinnings of disease symptoms and other health behaviors (e.g., Ranby et al., 2009, 2011), and enhance the evaluation of prevention/intervention programs (e.g., Petrosino, 2000; Donaldson, 2001; Fairchild and MacKinnon, 2014; Fairchild and McDaniel, 2017). In reaction to this sustained interest in use, methodologists have developed and analyzed a wide array of techniques to examine statistical mediation across a variety of data types, as well as have increased investigations into improving causal interpretation of model parameters. The use of discrete-time survival analysis to model time-to-event occurrence has also increased with applications across a wide range of developmental and behavioral outcomes such as youth alcohol and other drug involvement onset, first relapse following drug or alcohol treatment, early initiation of sexual activity, school expulsion, and timing of death (e.g., Capaldi et al., 1996; Nichols et al., 1999; Xie et al., 2005; Petras et al., 2011;Kerr et al., 2015; Atherton et al., 2016; Beydoun et al., 2016; Mason et al., 2017; Kilburn et al., 2018). To date, however, little methodological work has been conducted to study methods appropriate for examining mediation effects with discrete-time survival outcomes. This is a salient outstanding need.

Previous work has considered the examination of mediation effects with continuous time survival data (e.g., Tein and MacKinnon, 2003; Lange and Hansen, 2011; Tchetgen Tchetgen, 2011; Gelfand et al., 2016). The estimation and interpretation of model parameters across continuous vs. discrete-time survival analyses differ in critical ways, however, engendering the need to understand and develop discrete-time survival methods in their own right. These differences manifest in large part due to the “tied survival data” (e.g., Chalita et al., 2002) that materializes when one measures time-to-event occurrence in discrete epochs such as days, months, semesters, or years. Such periodic data collection, which is routine in academic settings and clinical trials, allows multiple individuals to experience an event within the same period, thus yielding grouped survival times. The ensuing focus of analysis in discrete-time models centers on the conditional probability of an individual experiencing the event of interest in a given time period, given that they have not experienced the event prior to that time (i.e., the hazard probability; Masyn, 2014). This focus is in contrast to an emphasis on the usual hazard rate in continuous-time survival models (e.g., Cox, 1972).

Fairchild et al. (2015) discussed key strengths of investigating mediation effects in discrete-time survival data with respect to ameliorating the nationwide public health problem of substance use in youth. Indeed, a brief survey of recent applied work suggests that there is a clear and present interest in examining mediation effects with discrete-time survival outcomes¹, particularly with research questions relating to different facets of substance use onset in youth. For example, Atherton et al. (2016) examined mediation effects in discrete-time survival data in a sample of Mexican-Origin youth to understand how familism impacted timing to onset of substance use through three proposed mediators: access to substances, intent to use substances, and association with deviant peers. In another study related to timing to substance use initiation, Mason et al. (2017) explored how cumulative contextual risk at the time of birth predicted timing to onset of substance use through childhood peer marginalization, aggression, and behavior problems. Kerr et al. (2015) also examined a substance use model that considered how parental use of marijuana in adolescence impacted time to onset of child marijuana use indirectly via a variety of risk and protective factors at both the child and parent-level. There have been applications beyond the substance use literature as well (Reis et al., 2011; Garner and Hunter, 2013; Hill et al., 2013; Lambert et al., 2013; Slater and Henry, 2013; Beydoun et al., 2016).

Varied tests of mediation were used throughout these applications. Although a handful of examples used contemporary, recommended approaches to test for mediation (e.g., bias-corrected bootstrapped confidence intervals, and potential outcomes approach), the majority of examples employed methods that have been shown to be limited in a variety of circumstances (see MacKinnon et al., 2002; Fritz and MacKinnon, 2007). That is, some of the applications used a point estimate and normal theory standard error estimator to test mediation but did not have sufficient sample size to invoke asymptotic efficiency (Sobel, 1982). Other studies assessed mediation via the causal steps (Baron and Kenny, 1986) and did not report a point estimate of the mediated effect, precluding tests of its statistical significance. An abundant amount of methodological work has demonstrated serious shortcomings with the latter approach to testing mediation effects (e.g., Shrout and Bolger, 2002; MacKinnon and Fairchild, 2009). Namely, the causal steps approach as originally described lacks sufficient power to detect mediation effects in a variety of circumstances (e.g., MacKinnon et al., 2002). This problem is particularly serious in complete mediation models, where (in the population) there is no direct effect of an independent variable on the outcome while controlling for the mediator. Indeed, Fritz and MacKinnon (2007) illustrated that the causal steps required approximately n = 21,000 subjects for adequate power to detect mediation when component paths of the mediated effect were small in size and the direct effect was zero. Others have also written on issues of power associated with this requirement of the test (e.g., Kenny and Judd, 2014; O’Rourke and MacKinnon, 2015)².

Taken as a set, the above examples serve to highlight the need for well-tested, uniform methods for estimating, and testing mediation effects on discrete time survival outcomes in the applied literature. Fairchild et al. (2015) introduced a framework to test mediation hypotheses with univariate discrete time survival data (i.e., non-recurring events) that incorporates current methodological recommendations from the broader mediation analysis literature. The authors demonstrated both a regression-based approach to analysis, as well as a structural equation modeling-based approach that estimates model parameters within a mixture-modeling framework. Neither variant, however, has yet to be evaluated statistically. Given the utility of examining mediation hypotheses with discrete-time survival outcomes, and the current lack of guidelines to do so, the formal evaluation of such a method is merited.

Current Study

Our goal in this paper is to evaluate statistical properties of the discrete-time survival mediation (DTSM) model. Specifically, we examine the models in a statistical simulation study to examine power, Type 1 error, and accuracy of model parameter estimates. Additionally, in line with an increased focus on causal frameworks to conduct statistical mediation analysis, we demonstrate analytical equivalency (under linearity assumptions) of the Fairchild et al. (2015) approach to that of Imai et al. (2010) potential outcomes framework for estimating the average causal mediated effect (ACME) and illustrate that simulation results of the two approaches are functionally indistinguishable. The intention of this study is to document the statistical properties of the DTSM model and to develop recommendations for its use in applied research. To best achieve study goals, we first orient readers to the DTSM model (for a more comprehensive overview of the model and component parts see Fairchild et al. (2015). We then briefly explain the potential outcomes framework for examining mediation effects, before describing the design, and findings of the statistical simulation study. We end by providing general recommendations and considerations for applied researchers, as well as by discussing limitations of our study.

Discrete-Time Survival Mediation Analysis

Mediation analysis allows a researcher to examine how one or more intermediate variables conveys the effect of a predictor to an outcome of interest, enhancing one’s understanding of theoretical relations in a model (Judd and Kenny, 1981; Baron and Kenny, 1986; MacKinnon et al., 2007a; MacKinnon, 2008). Discrete-time survival analysis supports research questions regarding both the “if” and “when” of an outcome by modeling the probability of event occurrence over discrete intervals of time (e.g., Allison, 1982; Singer and Willett, 1991; Singer and Willett, 1993; Willett and Singer, 2004; Masyn, 2014). This method lends insight into temporal windows of critical risk by providing information on risk for event occurrence in each interval. Drawing on mediation analysis and discrete-time survival analysis, the DTSM model allows testing of hypotheses about the mechanisms of event timing.

In line with conventional mediation models, the DTSM model parses the overall effect of an independent variable on the outcome (here the hazard probability of event occurrence) into a direct effect and an indirect effect. The direct effect in the DTSM model captures the influence of a predictor on the hazard probability of event occurrence controlling for the mediator, and the indirect effect captures the influence of the predictor on the hazard probability of event occurrence through the mediator.

With a single predictor, mediator, and event time, two equations define the DTSM model under a proportional hazard odds assumption that specifies invariant effects of X and M (see Figure 1).

FIGURE 1

Figure 1. Path diagram for the discrete-time survival mediation model with proportional odds constraint imposed for both the effects of M and X on the event history indicators, where X = the independent variable, M = the mediator variable, η = the latent propensity for event occurrence, and e₁-e_j = binary indicators of event occurrence at each time period.

$$M_i={\mathrm{\beta }}_{01}+a{X}_i+{\mathrm{\epsilon }}_{1i}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(1)$$

$$\mathrm{\eta =logit}\left(P_h\left(t_j|X_i,M_i\right)\right)={\mathrm{\beta }}_{02}+c^{\prime }{X}_i+b{M}_i+{\mathrm{\epsilon }}_{2i}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(2)$$

where P_h represents the hazard probability of the event outcome, and η is the logit of the latent propensity for event occurrence³. The hazard probability describes the probability of an individual experiencing the event in a particular time interval provided that individual is still event-free at the beginning of the interval. The logit hazard is defined as a linear weighted combination of the predictor X and mediator M. Indirect effects in the model can be defined in line with tracing rules from the path analysis literature (e.g., Wright, 1921; Alwin and Hauser, 1975), such that a parameter estimate of the mediation effect can be computed by multiplying the a coefficient from equation (1) and the b coefficient from equation (2). We consider the product of coefficients estimator with bootstrapped confidence limits in our study, as prior methodological work has shown that this is a desirable approach to conducting statistical mediation analysis with a variety of different data types (e.g., Fritz and MacKinnon, 2007; MacKinnon et al., 2007b; Preacher and Hayes, 2008; Taylor et al., 2008)⁴.

The model can be estimated using full-information maximum likelihood procedures as used in conventional logistic regression models (Singer and Willett, 1993; Muthén and Masyn, 2005). This estimation allows for missing data in the form of non-informative right censoring. Individuals with right-censored event times are individuals who are still “event free” at the time they cease to be observed in the study. FIML estimation can incorporate the partial information available about the event times for right-censored individuals into the likelihood function for the sample. Null hypothesis significance testing of mediated effects can be conducted by estimating asymmetric confidence limits for ab either via bootstrapping (Efron and Tibshirani, 1993) or the distribution of the product of two random variables (Aroian, 1947; Meeker et al., 1981).

We consider the simplest case of both X and M as time invariant variables with a time invariant effect for exposition purposes, but more complicated parameterizations of the model are possible. It is straightforward to incorporate time-varying predictors into the DTSM model, as well as to test for evidence of time-varying effects in model estimation. Though the DTSM model can be estimated in an SEM-based mixture-modeling framework to accommodate model extensions (see Fairchild et al., 2015), we limit our scope to the regression-based framework presented in Fairchild et al. (2015) in the interest of simplicity. This facilitates comparison of simulation results to Imai et al. (2010) causal mediation analysis framework. Note however, that in the simplest case as demonstrated here, maximum likelihood estimates of the hazard probabilities in the latent variable parameterization of the model are statistically equivalent to maximum likelihood estimates of the binary event indicators in an observed variable regression approach (see Muthén and Masyn, 2005).

Potential Outcomes Approach to DTSM

In theory, DTSM implies causal relations among variables in the model. However, as with other mediation models, real-world study design elements will often preclude investigators from drawing causal inferences even in the presence of longitudinal data. Though random assignment of the X variable allows the a path from the model to be treated as a causal estimate, the b and c’ coefficients cannot be interpreted as causal effects unless there is also random assignment of M. Given these limitations, methodologists have worked on developing analytical corrections to improve causal inference in mediation models (Imai et al., 2010; Coffman, 2011; Jo et al., 2011; VanderWeele, 2011). These methods permit one to move beyond associational estimators that consider how variables relate, toward causal estimators that instead evaluate how a given variable of interest causes another (in this case through a third, mediating variable). The frameworks also provide machinery to conduct sensitivity analyses to assess robustness of model effects against omitted confounders (e.g., VanderWeele, 2015). Though there are alternative frameworks available, such as Didelez’s (2018) recent work stemming from the decision theoretic perspective, we focus on the causal mediation approach presented by Imai et al. (2010), which utilizes counterfactuals and is an extension of the potential outcomes framework originally presented by Rubin (1974).

Suppose X is a binary variable, with X = 1 for treatment and X = 0 for control. Imai et al. (2010) define the causal mediation effect based on the potential outcomes framework as:

$${\mathrm{\delta }}_i(t)\equiv Y_i(t,M_i(1))-Y_i(t,M_i(0)),\mathrm{\hspace{0.17em}}for\mathrm{\hspace{0.17em}}t=0,\mathrm{\hspace{0.17em}}1\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(3)$$

where δ(x) represents the difference between the two potential outcomes, changing the mediator from M_i(1) to M_i(0) while holding the treatment constant at t. This quantity has also been termed the natural indirect effect in the literature. Note that for M_i(1) and M_i(0) only one of the quantities will actually be observed, and thus unit-level causal effects are not identifiable. If the unit participates in the control group _(Xi=0), then we will observe M_i(0) but not M_i(1). In contrast, if the unit participates in the treatment group _(Xi=1), then we will observe M_i(1) but not M_i(0). Accordingly, δ(0) reflects how the outcome would change if the mediator value was adjusted from the value that was observed under the control group M_i(0) to the value that would be observed if the participant was assigned to the treatment group M_i(1), and δ(1) reflects how the outcome would change if the mediator value was adjusted from the value that was observed under the treatment group M_i(1) to the value that would be observed if the participant was assigned to the control group M_i(0).

The average causal mediation effect (ACME) across participants in a sample is identifiable and given by:

$${\overline{\mathrm{\delta }}}_i(t)\equiv {\rm E}\left({\mathrm{\delta }}_i(t)\right)\equiv {\rm E}\left(Y_i(t,M_i(1))-Y_i(t,M_i(0))\right),\mathrm{\hspace{0.17em}}for\mathrm{\hspace{0.17em}}t=0,\mathrm{\hspace{0.17em}}1\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(4)$$

Under the sequential ignorability assumption, we can rewrite the equations (1) and (2) as:

$${\rm E}(M_i|X_i)={\mathrm{\beta }}_{01}+a{X}_i\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(5)$$

and

$${\rm E}(\log it(P_h(t|X_i,M_i)))={\mathrm{\beta }}_{02}+c^{\prime }{X}_i+b{M}_i\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(6)$$

Assuming that there is no interaction between the treatment and the ACME, the ACME is identified and given by $\overline{\mathrm{\delta }}$(1) = $\overline{\mathrm{\delta }}$(0) = ab. In logit form, the product of coefficients method is equivalent to the causal mediated effect under the potential outcomes framework.

Methods

To assess the statistical performance of the DTSM model, we conducted a Monte Carlo simulation study. Specifically, we conducted an experiment in which we had knowledge of population parameters by random sampling from probability distributions specified a priori. With access to these true parameter values, we defined and evaluated different statistical properties of model estimates. We conducted two variants of the study to demonstrate equivalency of the mediated effect estimates across the product of the coefficients and potential outcomes frameworks. The first study variant was conducted in Mplus Version 6.11 (Muthén and Muthén, 1998–2012) to estimate the product of the coefficients approach in a conventional, observed variable, linear structural model as defined by equations (1) and (2). The second study variant was conducted in R (R Core Team, 2013) using the “mediation” R package for Causal Mediation Analysis Tingley et al. (2014) to estimate the ACME as defined above. The latter approach proceeded in two steps. First, we specified the model for the conditional distribution of the mediator M given predictor variable X and fit a least squares regression through the linear model function. We also fit a model for the conditional distribution of the outcome Y given predictor variable X and mediator M by specifying a generalized linear mixed effect model, with a logit link to model the discrete-time survival data via a random effect component. After fitting the mediator and outcome models, we specified the fitted objects as well as the names of the treatment and mediating variables to the ‘mediate’ function to compute the estimated ACME.

Monte Carlo Study Design

Data Generation

We generated values of the exogenous X variable for the simulation studies with the Mersenne Twister pseudorandom number generator (Matsumoto and Nishimura, 1998). The default random number generator in R 2.13.1, the Mersenne Twister algorithm constructs a generally uniform distribution of data with a period of 2^(19937-1). Formal tests of statistical randomness have shown that the number sequences the Mersenne Twister generates are sufficiently random. We considered the DTSM model with a dichotomous X variable⁵, a specified that the data be randomly sampled from a binomial Bernoulli distribution, B(1,0.5), with probability mass function given by:

$$p\left(x\right)=\{\begin{array}{cc}{p}^x{(1-p)}^{1-x}& x\in \left\{0,1\right\}\\ 0& x\notin \left\{0,1\right\}\end{array}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(7)$$

We generated data for the M variable using equation (2), such that M was defined by a linear regression on X with β₀₁ = 0 and ε₁ randomly drawn from a standard normal distribution, N(0,1). We generated event histories under a proportional hazard odds model with a constant baseline hazard probability as given by:

$$\mathrm{Pr}\left(e_{ji}|X_i,M_i\right)=\frac{1}{1+\exp \left(-\left({\mathrm{\beta }}_{02}+c^{\prime }{X}_i+b{M}_i\right)\right)}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(8)$$

where β₀₂ was fixed to values that corresponded to specific baseline hazard probability values, e.g., β₀₂ = -2.94 corresponds to a baseline hazard probability of 0.05. Using a different random seed for each time interval, we randomly sampled values y_ij from a standard logistic distribution, L(0,1), with probability density function given by:

$$f(y)=\frac{\exp \left(-y\right)}{{\left(1+\exp \left(-y\right)\right)}^2}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(9)$$

If _{yij>−(β02+c′X+bM)}, then e_ji = 1, otherwise e_ji = 0.

After all e_ji values had been generated, values were recoded such that if, for a given case i, the first e_ji = 1 occurred at j=m, all remaining e_ji (i.e.,e_(m+1)i,…,e_Ji) were coded as missing, to yield non-informative fixed right-censoring at the end of the final time period.

Simulation Parameters

We varied 6 factors in the Monte Carlo study (see Table 1). We considered two variants of the number of measurement occasions associated with the event outcome to investigate how waves of measurement may affect DTSM parameter estimates. We also compared three different sample size conditions to ascertain how performance of the DTSM model fared across a range of sample sizes observed in the social sciences, and evaluated a variety of effect sizes for the a, b, and c’ parameters in the DTSM model to observe how the magnitude of different relations influenced properties of the ab parameter estimate in both full mediation and partial mediation contexts. Note that values of the b and c’ parameters reflect partially standardized odds ratios. Finally, we investigated two different baseline hazard conditions to assess in what ways baseline risk impacted performance of the DTSM model. Crossing all levels of each of these factors in a full factorial design yielded 384 unique experimental conditions. Replicating the model estimation process R = 500 times (e.g., MacKinnon et al., 1995; Nylund et al., 2007; Beyersmann et al., 2009) for each unique parameter combination created a sample of 192,000 datasets for analysis. A 6-way 2 × 3 × 4 × 4 × 2 × 2 ANOVA was conducted on the data to analyze the influence of each factor, and all possible interactions among the factors, on outcome variables under study. Due to the extremely large sample size of the experiment and corresponding inflation of statistical significance, we examined the practical significance of effects via η² to determine whether a given factor or interaction was meaningful. Specifically we interpreted the effect if η² ≥ 0.01, indicating that at least 1% of the variance in a given outcome was attributable to the effect of interest (Cohen, 1988).

TABLE 1

Table 1. Simulation factors and corresponding levels of each factor.

Simulation Outcome Measures

We evaluated four statistical properties to assess performance of the DTSM model in the simulation study. First, we looked at the accuracy of parameter estimates by assessing the mean relative bias of the ab estimator. We defined mean relative bias by ascertaining the average deviation of the estimator’s expected value from the true population parameter as a proportion of the true value in each unique parameter combination.

$$\widehat{RB}\left(\mathrm{\theta }\right)=\frac{\widehat{E}(\widehat{\mathrm{\theta }})-\mathrm{\theta }}{\mathrm{\theta }}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(10)$$

where

$$\widehat{E}\left(\widehat{\mathrm{\theta }}\right)=\frac{1}{R}{\displaystyle \sum }_{r=1}^R\left({\widehat{\mathrm{\theta }}}^{(r)}\right)\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(11)$$

${\widehat{\mathrm{\theta }}}^{\left(r\right)}$ refers to the sample estimate of $\widehat{a}\widehat{b}$ in a given dataset r, θ refers to the population value of ab in the experimental condition, r refers to the replication number, and R is the total number of replications. Values of relative bias ≤0.10 were considered desirable (Muthén et al., 1987). For those parameter combinations where the true value of ab was zero and relative bias was therefore undefined, we evaluated accuracy of parameter estimates by simply assessing the magnitude of unstandardized bias for the condition:

$$\widehat{B}\left(\mathrm{\theta }\right)=\widehat{E}\left(\widehat{\mathrm{\theta }}\right)-\mathrm{\theta }\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{(12)}$$

We then examined the mean squared error (MSE) as a measure of overall accuracy and precision for the ab estimator in the DTSM model:

$$\widehat{MSE}\left(\mathrm{\theta }\right)=\frac{1}{R}{\displaystyle \sum _{r=1}^R\left[{({\widehat{\mathrm{\theta }}}^{(r)}-\mathrm{\theta })}^2\right]=\widehat{var}\left(\widehat{\mathrm{\theta }}\right)+{\left(\widehat{B}\left(\mathrm{\theta }\right)\right)}^2}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(13)$$

where

$$\widehat{var}\left(\widehat{\mathrm{\theta }}\right)=\frac{1}{R}{\displaystyle \sum _{r=1}^R\left[{\left({\widehat{\mathrm{\theta }}}^{(r)}-\widehat{E}(\widehat{\mathrm{\theta }})\right)}^2\right]}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}(14)$$

Though often used as a comparative metric across two or more estimators (where the estimator with the smallest MSE is supported), we analyzed the metric in an absolute sense here to accompany our evaluation of bias in investigating how the simulation factors impact variability in the parameter estimates. As with both measures of bias, smaller values of MSE are desirable.

Finally, we evaluated statistical power and Type 1 error rates to evaluate accuracy of hypothesis testing of the ab estimate in the DTSM model. Specifically, we conducted significance testing by using the non-parametric percentile bootstrap to construct asymmetric confidence limits for the mediated effect. A resampling technique, the percentile bootstrap forms the 95% asymmetric confidence limits by empirically constructing the ab sampling distribution from the data and identifying the 2.5th and 97.5th quantiles of the distribution. The percentile bootstrap constructs an empirical sampling distribution of ab by taking repeated draws with replacement of size=n from the data and estimating $\widehat{a}\widehat{b}$ in each drawn sample. By using the data as a pseudo-population and creating a large number of bootstrap samples in which the parameter is estimated, the resampling technique effectively approximates the sampling distribution of the test statistic. In the studies conducted here, we utilized n = 1000 bootstrap draws (e.g., Fisher et al., 2016; Chowdhury et al., 2017; Fang et al., 2019).

For experimental conditions where the true value of ab was non-zero, we defined power as the proportion of simulation replications where the 95% asymmetric confidence interval did not include zero. We evaluated the results against the nominal 1-β = 0.80 criterion (Cohen, 1988), such that the statistical test was characterized as underpowered in conditions where empirical power was <0.80. For experimental conditions where the true value of ab was zero, we defined Type 1 error as the proportion of simulation replications where the 95% asymmetric confidence interval did not include zero. We evaluated results against the nominal α = 0.05, such that inflated Type 1 error rates were identified for those conditions where empirical α > 0.05.

Results

Accuracy and Precision of Parameter Estimates

Because results across the ab and ACME estimators were either identical or only differed in the third decimal place for all effects, we just plot results associated with the ab estimator. We report both estimates in the text when they differ. When X is binary and M is a continuous variable, mean relative bias of the ab estimator and ACME in the DTSM model (when ab≠0) were both 0.013, indicating that parameter estimates were within 1.3% of the true population values across simulation conditions. There were no main effects or interactions of practical significance on the relative bias outcome. In those parameter combinations where the true value of ab = 0 or ACME = 0, mean unstandardized bias of the estimators was 0.00007. There were no main effects or interactions of practical significance on the unstandardized bias outcome. These results indicate that the ab estimator and ACME in the DTSM model are unbiased when X is binary and M is continuous.

With regard to MSE across both the ab and ACME estimators, there were several practically significant main effects and interactions. The largest effect, that accounted for 20.9% of the variance in MSE, was a main effect of the size of the a path, η² = 0.209. We also found significant main effect of sample size, η² = 0.100 as well as the baseline hazard, η² = 0.025. However, there were significant interactions associated with these variables that warrant further discussion. First, there was a practically significant two-way interaction of sample size and the magnitude of a,η² = 0.065 (ab) and η² = 0.066 (ACME). The decrease in MSE associated with larger sample size was of greater magnitude for those conditions with a larger a path (see Figure 2). We also found a significant interaction of the number of waves of data collected and the size of the a path, η² = 0.011 (ab) and η² = 0.010 (ACME). Here the increase in MSE associated with an increase in the magnitude of the a path was mitigated by increasing waves of data collected (see Figure 3). That is, the increase in MSE with increasing size of the a coefficient was more negligible with 8 waves of data than with 4 waves of data. In addition, we detected a significant interaction of the baseline hazard and the size of the a path, η² = 0.027 (ab) and, η² = 0.028 (ACME). The increase in MSE associated with an increase in the magnitude of the a path was moderated by the size of the baseline hazard (see Figure 4). When the baseline hazard was greater, the increase in MSE associated with the a path was smaller.

FIGURE 2

Figure 2. Two-way interaction of sample size and the magnitude of a on MSE.

FIGURE 3

Figure 3. Two-way interaction of number of waves of data collected and the magnitude of a on MSE.

FIGURE 4

Figure 4. Two-way interaction of baseline hazard and the magnitude of a on MSE.

Finally, there was a practically significant main effect of the effect size of the b path, η² = 0.025 (ab) and η² = 0.024 (ACME). MSE of the mediated effect increase with increasing effect size of the b path (see Figure 5). Across both studies, average MSE was 0.002 for conditions that had a b path equal to 1.00 but increased to an average MSE of 0.004 when the b path was equal to 1.8.

FIGURE 5

Figure 5. Main effect of the effect size of the b path on MSE.

Power and Type I Error

The mean Type I error rate for the ab estimator in the DTSM model with binary X and continuous M was α = 0.043; similarly, the mean Type I error rate for the ACME was α = 0.039. There were no practically significant main effects or interactions among simulation factors of practical significance on the outcome. Evaluated against the nominal α = 0.05 criterion, the ab estimator and ACME in the DTSM model had acceptable Type I error rates across parameter combinations in these simulation studies.

There were several practically significant main effects and two-way interactions on the power to detect effects in the DTSM model. Specifically, there were significant main effects of sample size, η² = 0.077 (ab) and η² = 0.075 (ACME), the baseline hazard, η² = 0.015 (ab) and η² = 0.014 (ACME), the size of the b path, η² = 0.250 (ab) and η² = 0.257 (ACME), and the size of the a path, η² = 0.024 (ab) and η² = 0.026 (ACME). However, these effects should be interpreted in the context of their respective, practically significant interactions.

We found a practically significant two-way interaction between sample size and the size of the a path, η² = 0.021 (ab) and η² = 0.020 (ACME), as well as between sample size and the b path, η² = 0.027 (see Figures 6, 7). Power to detect ab increased with increasing sample size but this was moderated by the impact of the effect size of the a or b coefficient. The increase in power associated with increasing sample size was most marked when the size of a or b was small.

FIGURE 6

Figure 6. Two-way interaction between sample size and the size of the a path on power.

FIGURE 7

Figure 7. Two-way interaction between sample size and the size of the b path on power.

We also noted a practically significant two-way interaction between the baseline hazard and the size of the b path, η² = 0.017 (see Figure 8). In general, greater power was associated with a baseline hazard, however, this effect was moderated by the size of the b path. In conditions in which the b path was small, power rates were more disparate in conditions with a baseline hazard of 0.05 vs. 0.20. As the effect size of b increased, power rates in conditions with baseline hazard of 0.05 vs. 0.20 were more similar.

FIGURE 8

Figure 8. Two-way interaction between baseline hazard and the size of the b path on power.

Discussion

Discrete-time survival mediation analyses can strengthen behavioral research by improving one’s understanding of predictors of discretized periods of risk over time as well as further illuminating through what mechanisms different predictors impact the timing of an event. Further, DTSM analyses can be used to assess the impact of prevention or intervention programs on time-to-event outcomes. By incorporating the DTSM model into program evaluation efforts, researchers can enhance their understanding of how prevention or intervention programs achieve (or fail to achieve) effects on time-to-event outcomes.

This article evaluated the statistical performance of a DTSM model originally presented by Fairchild et al. (2015) in a Monte Carlo simulation study. It was shown that the conventional product of coefficients estimator and the ACME yielded comparable model results in the simulation conditions under study. Results showed favorable statistical properties of the model across a variety of different conditions. Specifically, MSE and Type 1 error for the model were tolerable across a range of different sample sizes, effect sizes, baseline hazard conditions, and number of measurement waves. Adequate power to detect mediation effects, however, required that effect sizes of the model parameters be medium or large. Note that considering the performance of a statistical test for mediation inherently involves null hypothesis significance testing, which is not without limitations. Given these limitations, and in line with other authors (e.g., Sullivan and Feinn, 2012), we encourage readers to report effect size estimates in any statistical model to assess the magnitude of an effect under study. Alternatively, researchers may consider different strategies entirely to investigate mediation effects (e.g., design approaches; Pirlott and MacKinnon, 2016; or Bayesian analyses; Enders et al., 2013; Miočević et al., 2018).

Limitations and Model Extensions

As this was the first examination of the DTSM model, we considered the simplest case of a single mediator model with time-invariant predictors that satisfied a proportional hazard odds assumption. These constraints are necessarily a limitation of the study that could be explored in subsequent research. Accommodating time-varying predictors and non-proportionality in the DSTMed model involves multiple a, b, and c’ paths, which necessarily introduces additional issues in estimating the mediated effect.

Future work might consider evaluating different latent variable model extensions of the DTSM model that permit examination of more complicated longitudinal processes. One such example is a finite mixture model that examines joint event history and growth processes (e.g., Malone et al., 2012). Such models permit the examination of reciprocal effects across the two processes, and may provide a good foundation for exploring multifaceted developmental theories. The estimation of non-parametric frailties to model unobserved heterogeneity in the survival process is also possible. Such random factors have been shown to be important in continuous survival models (e.g., Henderson and Oman, 1999), and may be particularly relevant to multivariate survival models where recurring events and/or competing risks are of interest.

Recommendations and Concluding Remarks

The utility of evaluating mediation effects spans across research domains, as such analyses strengthen prevention science and etiological research more generally. In order for the analyses to be maximally useful however, it is critical to employ methodology appropriate for the data under investigation (e.g., logistic or probit models for categorical outcomes; multilevel models for clustered data). The DTSM model presented in this article intends to facilitate the investigation of mediation effects with discrete-time survival data. By evaluating discrete-time mediation effects, we can improve our understanding of discretized periods of risk for myriad behavioral outcomes, as well as target prevention and intervention work more effectively. Interested researchers may consult Fairchild et al. (2015) for applied example syntax and a detailed interpretation of model parameters in the substantive context, for the model specification evaluated here⁶.

Simulation results indicate that parameter estimates of mediation in the DTSM model were statistically accurate and precise across the range of examined conditions. Type 1 error rates of the ab estimator and ACME were also tolerable in the conditions studied. Thus, substantive researchers can feel confident in using the DTSM model to make statistical inferences in circumstances where the model is adequately powered. With binary X variables, we cannot advocate that the model be used when the relation between X and M is posited to be small in effect, however. Though it may be possible that the model achieves adequate power with a small effect size of a in sample sizes larger than n = 1000, no such conditions were considered in this study precluding any comment on the matter.

Finally, as noted by Pratschke et al. (2016), a linear specification of the DTSM model using an applicable link function (such as the logit) allows for a causal interpretation of ab in this context. That is, the product of coefficients will yield an estimate of the natural indirect effect on a logit scale, lending support for ab as a viable estimator of mediation in linearized discrete-time survival models when there are no interaction effects. When incorporating non-linear relationships into the model, however, we encourage researchers to estimate ACME effects. The product of coefficients diverges from the causal estimator in this case, and ab no longer retains a causal interpretation (a limitation of the estimator). Fortunately, both methods to estimate mediation effects are readily available in a wide variety of software packages.

Author Contributions

AF designed, wrote, and approved all contributions to the study. CC ran all analysis for the work. HM and DS helped to revise the study. AG ran earlier analyses that contributed to the work. KM helped to design the study.

Funding

This study was supported in part by the National Institute on Drug Abuse grant R01DA030349, on which AF served as PI.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpsyg.2019.00740/full#supplementary-material

DATA SHEET S1 | Partial R Code for Simulation Study.

Footnotes

1. ^The focus of this paper is on mediation models, where a discrete-time onset variable serves as the outcome. Note that there are other types of mediation models where putative mediators can take the form of a discrete-time onset variable. Kilburn et al. (2018) examined one such model in the context of a randomized controlled trial to understand the impact of a conditional cash transfer intervention on reducing intimate partner violence for young women in South Africa.
2. ^Beyond issues associated with statistical power, the causal steps notion that there needs to be a significant overall relation of the predictor and outcome variables for mediation to exist is fallible from a variety of additional perspectives. Though testing the overall relationship between a predictor and outcome is important in any research study, it is a separate test from mediation and indirect effects can exist in the absence of a significant overall effect. This might occur in multiple mediator models where different mediated effect pathways have varying signs or when the sign of the mediated effect is opposite in value from the direct effect, yielding inconsistent mediation models.
3. ^Note: we use a latent underlying response variable parameterization, with latent variable η, as estimated in Mplus.
4. ^Previous work has noted the fallibility of estimating mediation with a difference in coefficients approach for both categorical and time-to-event outcomes (e.g., Lapointe-Shaw et al., 2018).
5. ^We specified a dichotomous X variable to align with a randomized controlled trial scenario where there would be treatment and control groups.
6. ^A simulation code template file for the current study is provided for interested readers in the Supplementary Materials to this manuscript.

References

Allison, P. D. (1982). “Discrete-time methods for the analysis of event histories,” in Sociological Methodology, ed. S. Leinhardt (San Francisco, CA: Jossey-Bass), 61–98. doi: 10.2307/270718

CrossRef Full Text | Google Scholar

Alwin, D. F., and Hauser, R. M. (1975). The decomposition of effects in path analysis. Am. Sociol. Rev. 40, 37–47.

Google Scholar

Aroian, L. A. (1947). The probability function of the product of two normally distributed variables. Ann. Math. Stat. 18, 265–271. doi: 10.1214/aoms/1177730442

PubMed Abstract | CrossRef Full Text | Google Scholar

Atherton, O. E., Conger, R. D., Ferrer, E., and Robins, R. W. (2016). Risk and protective factors for early substance use initiation: a longitudinal study of mexican-origin youth. J. Res. Adolesc. 26, 864–879. doi: 10.1111/jora.12235

PubMed Abstract | CrossRef Full Text | Google Scholar

Baron, R. M., and Kenny, D. A. (1986). The moderator-mediator variable distinction in social psychological research: conceptual, strategic, and statistical considerations. J. Personal. Soc. Psychol. 51, 1173–1182. doi: 10.1037/0022-3514.51.6.1173

CrossRef Full Text | Google Scholar

Beydoun, M. A., Beydoun, H. A., Mode, N., Dore, G. A., Canas, J. A., Eid, S. M., et al. (2016). Racial disparities in adult all-cause and cause-specific mortality among us adults: mediating and moderating factors. BMC Public Health 16:1113. doi: 10.1186/s12889-016-3744-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Beyersmann, J., Latouche, A., Buchholz, A., and Schumacher, M. (2009). Simulating competing risks data in survival analysis. Stat. Med. 28, 956–971. doi: 10.1002/sim.3516

PubMed Abstract | CrossRef Full Text | Google Scholar

Capaldi, D. M., Crosby, L., and Stoolmiller, M. (1996). Predicting the timing of first sexual intercourse for at-risk adolescent males. Child Dev. 67, 344–359. doi: 10.2307/1131818

PubMed Abstract | CrossRef Full Text | Google Scholar

Chalita, L. V. A. S., Colosimo, E. A., and Demétrio, C. G. B. (2002). Likelihood approximations and discrete models for tied survival data. Commun. Stat. Theory Methods 31, 1215–1229. doi: 10.1081/STA-120004920

CrossRef Full Text | Google Scholar

Chowdhury, M., Wu, C., and Modarres, R. (2017). Local Box–Cox transformation on time-varying parametric models for smoothing estimation of conditional CDF with longitudinal data. J. Stat. Comput. Simul. 87, 2900–2914. doi: 10.1080/00949655.2017.1347656

CrossRef Full Text | Google Scholar

Coffman, D. L. (2011). Estimating causal effects in mediation analysis using propensity scores. Struct. Equ. Modeling 18, 357–369. doi: 10.1080/10705511.2011.582001

PubMed Abstract | CrossRef Full Text | Google Scholar

Cohen, J. (1988). Statistical Power Analysis for the Behavioral Sciences, 2nd Edn. Hillsdale, NJ: Lawrence Earlbaum Associates.

Google Scholar

Cox, D. R. (1972). Regression models and life-tables (with discussion). J. R. Stat. Soc. Ser. B 34, 187–220. doi: 10.2307/2985181

CrossRef Full Text | Google Scholar

Donaldson, S. I. (2001). “Mediator and moderator analysis in program development,” in Handbook of Program Development for Health Behavior Research and Practice, ed. S. Sussman (Thousand Oaks, CA: Sage Publications, Inc), 470–500.

Google Scholar

Efron, B., and Tibshirani, R. J. (1993). An Introduction to the Bootstrap, Volume 57 of Monographs on Statistics and Applied Probability. New York, NY: Chapman & Hall.

Google Scholar

Enders, C. K., Fairchild, A. J., and MacKinnon, D. P. (2013). Using Bayesian parameter simulation to evaluate mediation effects with missing data. Multivariate Behav. Res. 48, 340–369. doi: 10.1080/00273171.2013.784862

PubMed Abstract | CrossRef Full Text

Fairchild, A. J., Abara, W. E., Gottschall, A. C., Tein, J. Y., and Prinz, R. J. (2015). Improving our ability to evaluate underlying mechanisms of behavioral onset and other event occurrence outcomes: a discrete-time survival mediation model. Eval. Health Prof. 38, 315–342. doi: 10.1177/0163278713512124

PubMed Abstract | CrossRef Full Text | Google Scholar

Fairchild, A. J., and MacKinnon, D. P. (2014). “Using mediation and moderation analyses to enhance prevention research,” in Advances in Prevention Science, Vol 1. Defining Prevention Science, eds Z. Sloboda and H. Petras (New York, NY: Springer), 537–555.

Google Scholar

Fairchild, A. J., and McDaniel, H. L. (2017). Best (but oft-forgotten) practices: mediation analysis. Am. J. Clin. Nutr. 105, 1259–1271. doi: 10.3945/ajcn.117.152546

PubMed Abstract | CrossRef Full Text | Google Scholar

Fang, J., Wen, Z., and Hau, K. T. (2019). Mediation effects in 2-1-1 multilevel model: evaluation of alternative estimation methods. Struct. Equ. Modeling 1–16. doi: 10.1080/10705511.2018.1547967

CrossRef Full Text | Google Scholar

Fisher, A., Caffo, B., Schwartz, B., and Zipunnikov, V. (2016). Fast, exact bootstrap principal component analysis for p > 1 million. J. Am. Stat. Assoc. 111, 846–860. doi: 10.1080/01621459.2015.1062383

PubMed Abstract | CrossRef Full Text | Google Scholar

Fritz, M. S., and MacKinnon, D. P. (2007). Required sample size to detect the mediated effect. Psychol. Sci. 18, 233–239. doi: 10.1111/j.1467-9280.2007.01882.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Garner, B. R., and Hunter, B. D. (2013). Examining the temporal relationship between psychological climate, work attitude, and staff turnover. J. Subst. Abuse Treat. 44, 193–200. doi: 10.1016/j.jsat.2012.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Gelfand, L. A., MacKinnon, D. P., DeRubeis, R. J., and Baraldi, A. N. (2016). Mediation analysis with survival outcomes: accelerated failure time vs. proportional hazards models. Front. Psychol. 7:423. doi: 10.3389/fpsyg.2016.00423

PubMed Abstract | CrossRef Full Text | Google Scholar

George, M. W., Fairchild, A. J., Cummings, E. M., and Davies, P. T. (2014). Marital conflict in early childhood and adolescent disordered eating: emotional insecurity as an explanatory mechanism. Eat. Behav. 15, 532–539. doi: 10.1016/j.eatbeh.2014.06.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Gonzalez, O., and MacKinnon, D. (2018). Statistical Mediation in Lifespan Developmental Analyses. Available at: http://oxfordre.com/psychology/view/10.1093/acrefore/9780190236557.001.0001/acrefore-9780190236557-e-358 (accessed March 13, 2018).

Henderson, R., and Oman, P. (1999). Effect of frailty on marginal regression estimates in survival analysis. J. R. Stat. Soc. B 61, 367–379. doi: 10.1111/1467-9868.00182

CrossRef Full Text | Google Scholar

Hill, K. G., Bailey, J. A., Hawkins, J. D., Catalano, R. F., Kosterman, R., Oesterle, S., et al. (2013). The onset of STI diagnosis through age 30: results from the seattle social development project intervention. Prev. Sci. 15 Suppl. 1:S19–S32. doi: 10.1007/s11121-013-0382-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Imai, K., Keele, L., and Yamamoto, T. (2010). Identification, inference and sensitivity analysis for causal mediation effects. Stat. Sci. 25, 51–71. doi: 10.1214/10-STS321

CrossRef Full Text | Google Scholar

Jo, B., Stuart, E. A., MacKinnon, D., and Vinokur, A. D. (2011). The use of propensity scores in mediation analysis. Multivariate Behav. Res. 46, 425–452. doi: 10.1080/00273171.2011.576624

PubMed Abstract | CrossRef Full Text | Google Scholar

Judd, C. M., and Kenny, D. A. (1981). Process analysis: estimating mediation in treatment evaluations. Eval. Rev. 5, 602–619. doi: 10.1177/0193841X8100500502

CrossRef Full Text | Google Scholar

Kenny, D. A., and Judd, C. M. (2014). Power anomalies in testing mediation. Psychol. Sci. 25, 334–339. doi: 10.1177/0956797613502676

PubMed Abstract | CrossRef Full Text | Google Scholar

Kerr, D. C., Tiberio, S. S., and Capaldi, D. M. (2015). Contextual risks linking parents’ adolescent marijuana use to offspring onset. Drug Alcohol Depend. 154, 222–228. doi: 10.1016/j.drugalcdep.2015.06.041

PubMed Abstract | CrossRef Full Text | Google Scholar

Kilburn, K. N., Pettifor, A., Edwards, J. K., Selin, A., Twine, R., MacPhail, C., et al. (2018). Conditional cash transfers and the reduction in partner violence for young women: an investigation of causal pathways using evidence from a randomized experiment in South Africa. J. Int. AIDS Soc. 21, 47–54. doi: 10.1002/jia2.25043

PubMed Abstract | CrossRef Full Text | Google Scholar

Lambert, S. F., Bettencourt, A. F., Bradshaw, C. P., and Ialongo, N. S. (2013). Early predictors of urban adolescents’ community violence exposure. J. Aggress. Maltreat. Trauma 22, 26–44. doi: 10.1080/10926771.2013.743944

CrossRef Full Text | Google Scholar

Lange, T., and Hansen, J. V. (2011). Direct and indirect effects in a survival context. Epidemiology 22, 575–581. doi: 10.1097/EDE.0b013e31821c680c

PubMed Abstract | CrossRef Full Text | Google Scholar

Lapointe-Shaw, L., Bouck, Z., Howell, N. A., Lange, T., Orchanian-Cheff, A., Austin, P. C., et al. (2018). Mediation analysis with a time-to-event outcome: a review of use and reporting in healthcare research. BMC Med. Res. Methodol. 18:118. doi: 10.1186/s12874-018-0578-7

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKinnon, D. P. (2008). Introduction to Statistical Mediation Analysis. Mahwah, NJ: Erlbaum.

Google Scholar

MacKinnon, D. P., and Fairchild, A. J. (2009). Current directions in mediation analysis. Curr. Direct. Psychol. 18, 16–20. doi: 10.1111/j.1467-8721.2009.01598.x

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKinnon, D. P., Fairchild, A. J., and Fritz, M. S. (2007a). Mediation analysis. Annu. Rev. Psychol. 58, 593–614. doi: 10.1146/annurev.psych.58.110405.085542

CrossRef Full Text | Google Scholar

MacKinnon, D. P., Lockwood, C. M., Brown, C. H., Wang, W., and Hoffman, J. M. (2007b). The intermediate endpoint effect in logistic and probit regression. Clin. Trials 4, 499–513. doi: 10.1177/1740774507083434

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKinnon, D. P., Lockwood, C. M., Hoffman, J. M., West, S. G., and Sheets, V. (2002). A comparison of methods to test mediation and other intervening variable effects. Psychol. Methods 7, 83–104. doi: 10.1037/1082-989X.7.1.83

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKinnon, D. P., Warsi, G., and Dwyer, J. H. (1995). A simulation study of mediated effect measures. Multivariate Behav. Res. 30, 41–62. doi: 10.1207/s15327906mbr3001_3

PubMed Abstract | CrossRef Full Text | Google Scholar

Malone, P. S., Northrup, T. F., Masyn, K. E., Lamis, D. A., and Lamont, A. E. (2012). Initiation and persistence of alcohol use in United States Black, Hispanic, and White male and female youth. Addict. Behav. 37, 299–305. doi: 10.1016/j.addbeh.2011.11.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Mason, W. A., Patwardhan, I., Smith, G. L., Chmelka, M. B., Savolainen, J., January, S. A. A., et al. (2017). Cumulative contextual risk at birth and adolescent substance initiation: peer mediation tests. Drug Alcohol Depend. 177, 291–298. doi: 10.1016/j.drugalcdep.2017.03.045

PubMed Abstract | CrossRef Full Text | Google Scholar

Masyn, K. (2014). “Advances in event history and survival models for prevention science,” in Advances in Prevention Science, Vol 1. Defining Prevention Science, eds Z. Sloboda and H. Petras (New York, NY: Springer), 537–555. doi: 10.1007/978-1-4899-7424-2_23

CrossRef Full Text

Matsumoto, M., and Nishimura, T. (1998). Mersenne twister: a 623-dimensionally equidistributed uniform pseudo-random number generator. ACM Trans. Modeling Comput. Simul. 8, 3–30. doi: 10.1145/272991.272995

CrossRef Full Text | Google Scholar

Meeker, W. G. Jr., Cornwell, L. W., and Aroian, L. A. (1981). Selected Tables in Mathematical Statistics, Volume VII: The Product of Two Normally Distributed Random Variables. Providence, RI: American Mathematical Society.

Google Scholar

Miočević, M., Gonzalez, O., Valente, M. J., and MacKinnon, D. P. (2018). A tutorial in bayesian potential outcomes mediation analysis. Struct. Equ. Modeling 25, 121–136. doi: 10.1080/10705511.2017.1342541

PubMed Abstract | CrossRef Full Text | Google Scholar

Muthén, B., Kaplan, D., and Hollis, M. (1987). On structural equation modeling with data that are not missing completely at random. Psychometrika 52, 431–462. doi: 10.1007/BF02294365

CrossRef Full Text | Google Scholar

Muthén, B., and Masyn, K. (2005). Discrete-time survival mixture analysis. J. Educ. Behav. Stat. 30, 27–58. doi: 10.3102/10769986030001027

CrossRef Full Text | Google Scholar

Muthén, L. K., and Muthén, B. O. (1998–2012). Mplus User’s Guide, 7th Edn. Los Angeles, CA: Muthén & Muthén.

Google Scholar

Nichols, J. D., Ludwin, W. G., and Iadicola, P. (1999). A darker shade of gray: a year-end analysis of discipline and suspension data. Equity Excell. Educ. 32, 43–55. doi: 10.1080/1066568990320105

CrossRef Full Text | Google Scholar

Nylund, K. L., Asparouhov, T., and Muthén, B. O. (2007). Deciding on the number of classes in latent class analysis and growth mixture modeling: a monte carlo simulation study. Struct. Equ. Modeling 14, 535–569. doi: 10.1080/10705510701575396

CrossRef Full Text | Google Scholar

O’Rourke, H. P., and MacKinnon, D. P. (2015). When the test of mediation is more powerful than the test of the total effect. Behav. Res. Methods 47, 424–442. doi: 10.3758/s13428-014-0481-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Petras, H., Masyn, K. E., Buckley, J. A., Ialongo, N. S., and Kellam, S. (2011). Who is most at risk for school removal? A multilevel discrete-time survival analysis of individual- and context-level influences. J. Educ. Psychol. 103, 223–237. doi: 10.1037/a0021545

CrossRef Full Text | Google Scholar

Petrosino, A. (2000). Mediators and moderator in the evaluation of programs for children: current practice and agenda for improvement. Eval. Rev. 24, 47–72. doi: 10.1177/0193841X0002400102

PubMed Abstract | CrossRef Full Text | Google Scholar

Pirlott, A. G., and MacKinnon, D. P. (2016). Design approaches to experimental mediation. J. Exp. Soc. Psychol. 66, 29–38. doi: 10.1016/j.jesp.2015.09.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Pratschke, J., Haase, T., Comber, H., Sharp, L., de Camargo Cancela, M., and Johnson, H. (2016). Mechanisms and mediation in survival analysis: towards an integrated analytical framework. BMC Med. Res. Methodol. 16:27. doi: 10.1186/s12874-016-0130-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Preacher, K. J., and Hayes, A. F. (2008). Asymptotic and resampling strategies for assessing and comparing indirect effects in multiple mediator models. Behav. Res. Methods 40, 879–891. doi: 10.3758/BRM.40.3.879

CrossRef Full Text | Google Scholar

R Core Team (2013). R: A Language and Environment for Statistical Computing. Vienna: R Foundation for Statistical Computing.

Google Scholar

Ranby, K. W., Aiken, L. S., MacKinnon, D. P., Elliot, D., Moe, E. L., McGinnis, W., et al. (2009). A mediation analysis of the ATHENA intervention for female athletes: prevention of athletic-enhancing substance use and unhealthy weight loss behaviors. J. Pediatr. Psychol. 34, 1069–1083. doi: 10.1093/jpepsy/jsp025

PubMed Abstract | CrossRef Full Text | Google Scholar

Ranby, K. W., MacKinnon, D. P., Fairchild, A. J., Elliot, D. L., Kuehl, K. S., and Goldberg, L. (2011). The PHLAME (promoting healthy lifestyles: alternative models’ effects) firefighter study: testing mediating mechanisms. J. Occup. Health Psychol. 16, 501–513. doi: 10.1037/a0023002

PubMed Abstract | CrossRef Full Text | Google Scholar

Reis, O., Dornte, M., and von der Lippe, H. (2011). Neuroticism, social support, and the timing of first parenthood: a prospective study. Personal. Individ. Differ. 50, 381–386. doi: 10.1016/j.paid.2010.10.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Rubin, D. B. (1974). Estimating causal effects of treatments in randomized and nonrandomized studies. J. Educ. Psychol. 66, 688–701. doi: 10.1037/h0037350

CrossRef Full Text | Google Scholar

Selig, J. P., and Preacher, K. J. (2009). Mediation models for longitudinal data in developmental research. Res. Hum. Dev. 6, 144–164. doi: 10.1080/15427600902911247

CrossRef Full Text | Google Scholar

Shrout, P. E., and Bolger, N. (2002). Mediation in experimental and nonexperimental studies: new procedures and recommendations. Psychol. Methods 7, 422–445. doi: 10.1037/1082-989X.7.4.422

PubMed Abstract | CrossRef Full Text | Google Scholar

Singer, J. D., and Willett, J. B. (1991). Modeling the days of our lives: using survival analysis when designing and analyzing longitudinal studies of duration and timing of events. Psychol. Bull. 110, 268–290. doi: 10.1037/0033-2909.110.2.268

CrossRef Full Text | Google Scholar

Singer, J. D., and Willett, J. B. (1993). It’s about time: using discrete-time survival analysis to study duration and the timing of events. J. Educ. Stat. 18, 155–195. doi: 10.2307/1165085

CrossRef Full Text | Google Scholar

Slater, M. D., and Henry, K. L. (2013). Prospective influence of music-related media exposure on adolescent substance-use initiation: a peer group mediation model. J. Health Commun. 18, 291–305. doi: 10.1080/10810730.2012.727959

PubMed Abstract | CrossRef Full Text | Google Scholar

Sobel, M. E. (1982). Asymptotic confidence intervals for indirect effects in structural equation models. Sociol. Methodol. 13, 290–312. doi: 10.2307/270723

CrossRef Full Text | Google Scholar

Sullivan, G. M., and Feinn, R. (2012). Using effect size—or why the P value is not enough. J. Grad. Med. Educ. 4, 279–282. doi: 10.4300/JGME-D-12-00156.1

PubMed Abstract | CrossRef Full Text | Google Scholar

Taylor, A. B., MacKinnon, D. P., and Tein, J. Y. (2008). Tests of the three-path mediated effect. Organ. Res. Methods 11, 241–269. doi: 10.1177/1094428107300344

PubMed Abstract | CrossRef Full Text | Google Scholar

Tchetgen Tchetgen, E. J. (2011). On causal mediation analysis with a survival outcome. Int. J. Biostat. 7:33. doi: 10.2202/1557-4679.1351

PubMed Abstract | CrossRef Full Text | Google Scholar

Tein, J.-Y., and MacKinnon, D. P. (2003). “Estimating mediation effects with survival data,” in New Developments on Psychometrics, eds H. Yanai, A. O. Rikkyo, K. Shigemasu, Y. Kano, and J. J. Meulman (Tokyo: Springer-Verlag Tokyo Inc), 405–412.

Tingley, D., Yamamoto, T., Hirose, K., Keele, L., and Imai, K. (2014). Mediation: R package for causal mediation analysis. J. Stat. Softw. 59, 1–38. doi: 10.1037/a0036434

PubMed Abstract | CrossRef Full Text | Google Scholar

VanderWeele, T. J. (2011). Causal mediation analysis with survival data. Epidemiology 22, 582–585. doi: 10.1097/EDE.0b013e31821db37e

PubMed Abstract | CrossRef Full Text | Google Scholar

VanderWeele, T. J. (2015). Explanation in Causal Inference: Methods for Mediation and Interaction. Oxford: Oxford University Press.

Google Scholar

Willett, J. B., and Singer, J. D. (2004). “Discrete-time survival analysis,” in The SAGE Handbook of Quantitative Methodology for the Social Sciences, ed. D. Kaplan (Thousand Oaks, CA: Sage Publications), 199–211.

Google Scholar

Wright, S. (1921). Correlation and causation. J. Agric. Res. 20, 557–585.

Google Scholar

Xie, H., McHugo, G. J., Fox, M. B., and Drake, R. E. (2005). Substance abuse relapse in a ten-year prospective follow-up of clients with mental and substance use disorders. Psychiatr. Serv. 56, 1282–1287. doi: 10.1176/appi.ps.56.10.1282

PubMed Abstract | CrossRef Full Text | Google Scholar