As reported on Twitter, my paper on causal mediation analysis with multiple mediators has been published online:
I have already discussed the theoretical background of multiple mediators in other posts (find them here, here, here, and here). Hence, and as promised, I will instead focus on the code that you need to replicate the results of the paper.
I will start with the post-treatment confounding case. First, and to provide a very general example, the Cs are the pre-treatment covariates (C1=gender; C2=age, C3=ethnicity), the X is the treatment, the M1 is normative alignment with the police, the M2 is duty to obey, and the Y is willingness to cooperate with the police.
As discussed in detail in the paper, you need to derive the squared terms of the mediators and the possible interaction between the mediator and the treatment.
rename C1 c1 rename C2 c2 rename C3 c3 rename X x rename Y y rename M1 m gen m2=m*m gen xm=x*m rename M2 ma gen ma2=ma*ma gen xma=x*ma
After creating the new variables, the only thing left is to run the analysis. Again, the model specification is outlined in the paper. The Monte Carlo simulations are intentionally set to a very high number to assure reliable effect size estimates. By contrast, the number of bootstraps for the standard errors could be conceivably increased, but unless you have access to a cluster-computing system, it is not advisable.
*Full number of bootstraps and MINSIM #delimit ; gformula y x m m2 ma ma2 c1 c2 c3 xma, mediation outcome(y) exposure(x) mediator(m) post_confs(ma) base_confs(c1 c2 c3) obe control(m:0) commands(y:regress, m:regress, ma:regress) equations(y:x m m2 ma ma2 c1 c2 c3 xma, m:x ma c1 c2 c3 xma, ma:x c1 c2 c3) derived(m2 ma2 xma) derrules(m2:m*m, ma2:ma*ma, xma:x*ma) minsim samples(300) moreMC simulations(100000) replace seed(215); #delimit cr
Finally, in case you want to replicate the sensitivity analysis, you can find the code for that below:
*Sensitivity analysis - Defining sens_rho.ado. capture program drop sens_rho program define sens_rho, rclass version 12 preserve cap matrix drop Psi sem (y <- x ma xma ma2 c1 c2 c3) (ma<- x c1 c2 c3) (m <- x ma xma ma2 c1 c2 c3), /// nocapslatent cov(e.y*e.m) qui estat framework, fitted matrix Psi=r(Psi) matrix list Psi scalar rho_dash=(Psi[3,1])/(sqrt(Psi[1,1]*Psi[3,3])) scalar list rho_dash return scalar rho=rho_dash restore end *Trial whether the specification was succesful. sens_rho
For estimating the confidence intervals for the sensitivity estimate 'rho', you can derive those this way:
bootstrap rho_dash=r(rho), reps(500) saving(sens_rho,replace):sens_rho estat bootstrap, all
Turning to the sequentially ordered mediators, the code is much more complex. This is due to the high number of indirect effects that need to be estimated, and also because the results are being conditional on the sensitivity parameter. As a result of this complexity, I will not go into detail regarding the code, I uploaded the whole do file instead. Feel free to ask for further clarifications should you get stuck.
The data is available in this dropbox folder. Because the post-treatment confounding approach can be a time-consuming exercise, I also added the log of the results of this analysis. Finally, you can find the do file for the sequentially ordered mediators in the corresponding folder.
UPDATE 18.06.2019: I have been informed that the link did not work, I have hopefully fixed the issue now.