Back to Index

R mediation package

• The dataset is simulated as having 1/3 of mediation effect
• ACME: average causal mediation effects is the total effect minus the direct effect
• ADE: average direct effect

library(mediation)
dt <- data.table(x = rnorm(100))
dt <- dt[, noise1 := rnorm(100) / 10
         ][, noise2 := rnorm(100) / 10
           ][, noise3 := rnorm(100) / 10
             ][, mediator := (x + noise1) * (1 + noise2) + noise3
               ][, y := mediator * (1 + noise3) + (x + noise2) * (2 + noise1)]
dt <- dt[x > -2][x < 2]

m <- lm(mediator ~ x, data = dt)
o <- lm(y ~ x + mediator, data = dt)

md <- mediate(m, o, sims = 1000
            , treat = "x"
            , mediator = "mediator"
            , boot = TRUE
            , control.value = -1
            , treat.value = 1
              )


mds <- summary(md)


extract_mediation_summary <- function (x) { 

  clp <- 100 * x$conf.level
  isLinear.y <- ((class(x$model.y)[1] %in% c("lm", "rq")) || 
                   (inherits(x$model.y, "glm") && x$model.y$family$family == 
                      "gaussian" && x$model.y$family$link == "identity") || 
                   (inherits(x$model.y, "survreg") && x$model.y$dist == 
                      "gaussian"))

  printone <- !x$INT && isLinear.y

  if (printone) {

    smat <- c(x$d1, x$d1.ci, x$d1.p)
    smat <- rbind(smat, c(x$z0, x$z0.ci, x$z0.p))
    smat <- rbind(smat, c(x$tau.coef, x$tau.ci, x$tau.p))
    smat <- rbind(smat, c(x$n0, x$n0.ci, x$n0.p))

    rownames(smat) <- c("ACME", "ADE", "Total Effect", "Prop. Mediated")

  } else {
    smat <- c(x$d0, x$d0.ci, x$d0.p)
    smat <- rbind(smat, c(x$d1, x$d1.ci, x$d1.p))
    smat <- rbind(smat, c(x$z0, x$z0.ci, x$z0.p))
    smat <- rbind(smat, c(x$z1, x$z1.ci, x$z1.p))
    smat <- rbind(smat, c(x$tau.coef, x$tau.ci, x$tau.p))
    smat <- rbind(smat, c(x$n0, x$n0.ci, x$n0.p))
    smat <- rbind(smat, c(x$n1, x$n1.ci, x$n1.p))
    smat <- rbind(smat, c(x$d.avg, x$d.avg.ci, x$d.avg.p))
    smat <- rbind(smat, c(x$z.avg, x$z.avg.ci, x$z.avg.p))
    smat <- rbind(smat, c(x$n.avg, x$n.avg.ci, x$n.avg.p))

    rownames(smat) <- c("ACME (control)", "ACME (treated)", 
                        "ADE (control)", "ADE (treated)", "Total Effect", 
                        "Prop. Mediated (control)", "Prop. Mediated (treated)", 
                        "ACME (average)", "ADE (average)", "Prop. Mediated (average)")

  }

  colnames(smat) <- c("Estimate", paste(clp, "% CI Lower", sep = ""), 
                      paste(clp, "% CI Upper", sep = ""), "p-value")
  smat

}

prt_med <- extract_mediation_summary(summary(md))


prt_med %>% prt()

R sessionInfo

sessionInfo()

R version 4.2.0 (2022-04-22) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 20.04.3 LTS

Matrix products: default BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0 LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0

locale: [1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
[4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
[7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
[10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C

attached base packages: [1] stats graphics grDevices utils datasets methods base

other attached packages: [1] mediation_4.5.0 sandwich_3.0-2 mvtnorm_1.1-3
[4] MASS_7.3-54 Wu_0.0.0.9000 flexdashboard_0.5.2 [7] lme4_1.1-29 Matrix_1.4-0 mgcv_1.8-38
[10] nlme_3.1-152 png_0.1-7 scales_1.2.0
[13] nnet_7.3-16 labelled_2.9.1 kableExtra_1.3.4
[16] plotly_4.10.0 gridExtra_2.3 ggplot2_3.3.6
[19] DT_0.23 tableone_0.13.2 magrittr_2.0.3
[22] lubridate_1.8.0 dplyr_1.0.9 plyr_1.8.7
[25] data.table_1.14.2 rmdformats_1.0.4 knitr_1.39

loaded via a namespace (and not attached): [1] webshot_0.5.3 RColorBrewer_1.1-3 httr_1.4.3
[4] backports_1.4.1 tools_4.2.0 bslib_0.3.1
[7] utf8_1.2.2 R6_2.5.1 rpart_4.1-15
[10] Hmisc_4.7-0 DBI_1.1.2 lazyeval_0.2.2
[13] colorspace_2.0-3 withr_2.5.0 tidyselect_1.1.2
[16] compiler_4.2.0 cli_3.3.0 rvest_1.0.2
[19] htmlTable_2.4.0 xml2_1.3.3 bookdown_0.27
[22] sass_0.4.1 checkmate_2.1.0 systemfonts_1.0.4
[25] stringr_1.4.0 digest_0.6.29 foreign_0.8-81
[28] minqa_1.2.4 rmarkdown_2.14 svglite_2.1.0
[31] base64enc_0.1-3 jpeg_0.1-9 pkgconfig_2.0.3
[34] htmltools_0.5.3 highr_0.9 fastmap_1.1.0
[37] htmlwidgets_1.5.4 rlang_1.0.4 rstudioapi_0.13
[40] jquerylib_0.1.4 generics_0.1.3 zoo_1.8-10
[43] jsonlite_1.8.0 Formula_1.2-4 Rcpp_1.0.9
[46] munsell_0.5.0 fansi_1.0.3 lifecycle_1.0.1
[49] stringi_1.7.8 yaml_2.3.5 grid_4.2.0
[52] forcats_0.5.1 lattice_0.20-45 haven_2.5.0
[55] splines_4.2.0 hms_1.1.1 klippy_0.0.0.9500
[58] pillar_1.8.0 boot_1.3-28 lpSolve_5.6.15
[61] glue_1.6.2 evaluate_0.15 mitools_2.4
[64] latticeExtra_0.6-29 vctrs_0.4.1 nloptr_2.0.3
[67] gtable_0.3.0 purrr_0.3.4 tidyr_1.2.0
[70] assertthat_0.2.1 xfun_0.31 survey_4.1-1
[73] survival_3.2-13 viridisLite_0.4.0 tibble_3.1.8
[76] cluster_2.1.2 ellipsis_0.3.2