Mediation Analysis Examples
Hairong Song
Oct/16/2024
The code below demonstrates the uses of two R packages for running mediation analysis.
Install Packages
install.packages("mediation",repos = "http://cran.us.r-project.org")## package 'mediation' successfully unpacked and MD5 sums checked
##
## The downloaded binary packages are in
## C:\Users\song2170\AppData\Local\Temp\RtmpuQxP3O\downloaded_packagesinstall.packages("lavaan",repos = "http://cran.us.r-project.org")## package 'lavaan' successfully unpacked and MD5 sums checked
##
## The downloaded binary packages are in
## C:\Users\song2170\AppData\Local\Temp\RtmpuQxP3O\downloaded_packagesLoad Packages
library(mediation)
library(lavaan)1.Mediation Analysis Using “Mediation” Package
set.seed (2024)
data("framing", package="mediation")
# mediation model with covariates
med.fit = lm(emo ~ treat + age + educ + gender + income, data = framing)
# outcome model with covariates
out.fit = glm(cong_mesg ~ emo + treat + age + educ + gender + income,
data = framing, family = binomial("probit"))
med.out1 = mediate(med.fit, out.fit, treat = "treat", mediator = "emo",
robustSE = TRUE, sims = 1000)
summary(med.out1)##
## Causal Mediation Analysis
##
## Quasi-Bayesian Confidence Intervals
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME (control) 0.0817 0.0307 0.14 0.002 **
## ACME (treated) 0.0821 0.0307 0.14 0.002 **
## ADE (control) 0.0120 -0.1028 0.13 0.840
## ADE (treated) 0.0123 -0.1114 0.13 0.840
## Total Effect 0.0941 -0.0319 0.21 0.148
## Prop. Mediated (control) 0.7881 -3.8596 6.53 0.150
## Prop. Mediated (treated) 0.8050 -3.4385 6.06 0.150
## ACME (average) 0.0819 0.0302 0.14 0.002 **
## ADE (average) 0.0122 -0.1062 0.13 0.840
## Prop. Mediated (average) 0.7965 -3.6084 6.26 0.150
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sample Size Used: 265
##
##
## Simulations: 1000plot(med.out1)
Bootstrapping method for indirect effect
med.out2 = mediate(med.fit, out.fit, boot = TRUE, treat = "treat", mediator = "emo", sims = 1000) ## Running nonparametric bootstrapsummary(med.out2)##
## Causal Mediation Analysis
##
## Nonparametric Bootstrap Confidence Intervals with the Percentile Method
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME (control) 0.0828 0.0358 0.14 <2e-16 ***
## ACME (treated) 0.0837 0.0360 0.14 <2e-16 ***
## ADE (control) 0.0115 -0.0969 0.13 0.85
## ADE (treated) 0.0123 -0.1097 0.14 0.85
## Total Effect 0.0951 -0.0306 0.22 0.12
## Prop. Mediated (control) 0.8707 -5.9094 6.71 0.12
## Prop. Mediated (treated) 0.8796 -5.2294 6.22 0.12
## ACME (average) 0.0833 0.0357 0.14 <2e-16 ***
## ADE (average) 0.0119 -0.1039 0.13 0.85
## Prop. Mediated (average) 0.8751 -5.5696 6.43 0.12
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sample Size Used: 265
##
##
## Simulations: 10002. Mediation Analysis Using “lavaan” Package
model = 'cong_mesg ~ c*treat # direct effect
emo ~ a*treat # mediator
cong_mesg ~ b*emo
ab := a*b # indirect effect (a*b)
total := c + (a*b) # total effect
'
fit = sem(model, data = framing)
summary(fit, fit.measures=TRUE, standardized=TRUE,rsquare=TRUE) ## lavaan 0.6-19 ended normally after 1 iteration
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 5
##
## Number of observations 265
##
## Model Test User Model:
##
## Test statistic 0.000
## Degrees of freedom 0
##
## Model Test Baseline Model:
##
## Test statistic 54.185
## Degrees of freedom 3
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -795.365
## Loglikelihood unrestricted model (H1) -795.365
##
## Akaike (AIC) 1600.731
## Bayesian (BIC) 1618.629
## Sample-size adjusted Bayesian (SABIC) 1602.776
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.000
## P-value H_0: RMSEA <= 0.050 NA
## P-value H_0: RMSEA >= 0.080 NA
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## cong_mesg ~
## treat (c) 0.015 0.063 0.231 0.818 0.015 0.014
## emo ~
## treat (a) 1.480 0.379 3.906 0.000 1.480 0.233
## cong_mesg ~
## emo (b) 0.063 0.010 6.276 0.000 0.063 0.368
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .cong_mesg 0.191 0.017 11.511 0.000 0.191 0.862
## .emo 7.253 0.630 11.511 0.000 7.253 0.946
##
## R-Square:
## Estimate
## cong_mesg 0.138
## emo 0.054
##
## Defined Parameters:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## ab 0.093 0.028 3.316 0.001 0.093 0.086
## total 0.107 0.066 1.626 0.104 0.107 0.099Bootstrapping method for indirect effects
fit2 = sem(model, data = framing,se ="bootstrap",bootstrap=5000)
parameterEstimates(fit2, ci=TRUE, level=0.95, boot.ci.type="perc", standardized=TRUE) ## lhs op rhs label est se z pvalue ci.lower ci.upper
## 1 cong_mesg ~ treat c 0.015 0.064 0.228 0.820 -0.110 0.142
## 2 emo ~ treat a 1.480 0.377 3.921 0.000 0.718 2.214
## 3 cong_mesg ~ emo b 0.063 0.010 6.465 0.000 0.044 0.082
## 4 cong_mesg ~~ cong_mesg 0.191 0.012 16.588 0.000 0.166 0.210
## 5 emo ~~ emo 7.253 0.486 14.917 0.000 6.256 8.118
## 6 treat ~~ treat 0.191 0.000 NA NA 0.191 0.191
## 7 ab := a*b ab 0.093 0.028 3.279 0.001 0.041 0.153
## 8 total := c+(a*b) total 0.107 0.069 1.560 0.119 -0.023 0.242
## std.lv std.all std.nox
## 1 0.015 0.014 0.031
## 2 1.480 0.233 0.534
## 3 0.063 0.368 0.368
## 4 0.191 0.862 0.862
## 5 7.253 0.946 0.946
## 6 0.191 1.000 0.191
## 7 0.093 0.086 0.197
## 8 0.107 0.099 0.228