if (FALSE) { # require("mediation") && require("brms") && require("rstanarm")
# \donttest{
library(mediation)
library(brms)
library(rstanarm)
# load sample data
data(jobs)
set.seed(123)
# linear models, for mediation analysis
b1 <- lm(job_seek ~ treat + econ_hard + sex + age, data = jobs)
b2 <- lm(depress2 ~ treat + job_seek + econ_hard + sex + age, data = jobs)
# mediation analysis, for comparison with Stan models
m1 <- mediate(b1, b2, sims = 1000, treat = "treat", mediator = "job_seek")
# Fit Bayesian mediation model in brms
f1 <- bf(job_seek ~ treat + econ_hard + sex + age)
f2 <- bf(depress2 ~ treat + job_seek + econ_hard + sex + age)
m2 <- brm(f1 + f2 + set_rescor(FALSE), data = jobs, refresh = 0)
# Fit Bayesian mediation model in rstanarm
m3 <- suppressWarnings(stan_mvmer(
list(
job_seek ~ treat + econ_hard + sex + age + (1 | occp),
depress2 ~ treat + job_seek + econ_hard + sex + age + (1 | occp)
),
data = jobs,
refresh = 0
))
summary(m1)
mediation(m2, centrality = "mean", ci = 0.95)
mediation(m3, centrality = "mean", ci = 0.95)
# }
}
Run the code above in your browser using DataLab