# NOT RUN {
## From the mediation tutorial:
## http://lavaan.ugent.be/tutorial/mediation.html
set.seed(1234)
X <- rnorm(100)
M <- 0.5*X + rnorm(100)
Y <- 0.7*M + rnorm(100)
dat <- data.frame(X = X, Y = Y, M = M)
mod <- ' # direct effect
Y ~ c*X
# mediator
M ~ a*X
Y ~ b*M
# indirect effect (a*b)
ind := a*b
# total effect
total := ind + c
'
fit <- sem(mod, data = dat)
summary(fit, ci = TRUE) # print delta-method CIs
## Automatically extract information from lavaan object
set.seed(1234)
monteCarloCI(fit) # CIs more robust than delta method in smaller samples
## save samples to calculate more precise intervals:
# }
# NOT RUN {
set.seed(1234)
foo <- monteCarloCI(fit, append.samples = TRUE)
library(HDInterval)
hdi(fit$Samples)
# }
# NOT RUN {
## Parameters can also be obtained from an external analysis
myParams <- c("a","b","c")
(coefs <- coef(fit)[myParams]) # names must match those in the "expression"
## Asymptotic covariance matrix from an external analysis
(AsyCovMat <- vcov(fit)[myParams, myParams])
## Compute CI, include a plot
set.seed(1234)
monteCarloCI(expr = c(ind = 'a*b', total = 'ind + c',
## other arbitrary functions are also possible
meaningless = 'sqrt(a)^b / log(abs(c))'),
coefs = coefs, ACM = AsyCovMat,
plot = TRUE, ask = TRUE) # print a plot for each
# }
Run the code above in your browser using DataLab