Example 1: The Simple Mediation Model

In this example, the Monte Carlo method is used to generate confidence intervals for the indirect effect in a simple mediation model where variable X has an effect on variable Y, through a mediating variable M. This mediating or indirect effect is a product of path coefficients from the fitted model.

The Simple Mediation Model

library(semmcci)
library(lavaan)

Data

n <- 1000
a <- 0.50
b <- 0.50
cp <- 0.25
s2_em <- 1 - a^2
s2_ey <- 1 - cp^2 - a^2 * b^2 - b^2 * s2_em - 2 * cp * a * b
em <- rnorm(n = n, mean = 0, sd = sqrt(s2_em))
ey <- rnorm(n = n, mean = 0, sd = sqrt(s2_ey))
X <- rnorm(n = n)
M <- a * X + em
Y <- cp * X + b * M + ey
df <- data.frame(X, M, Y)

Model Specification

The indirect effect is defined by the product of the slopes of paths X to M labeled as a and M to Y labeled as b. In this example, we are interested in the confidence intervals of indirect defined as the product of a and b using the := operator in the lavaan model syntax.

model <- "
  Y ~ cp * X + b * M
  M ~ a * X
  indirect := a * b
  direct := cp
  total := cp + (a * b)
"

Model Fitting

We can now fit the model using the sem() function from lavaan.

fit <- sem(data = df, model = model)

Monte Carlo Confidence Intervals

The fit lavaan object can then be passed to the MC() function to generate Monte Carlo confidence intervals.

MC(fit, R = 20000L, alpha = c(0.001, 0.01, 0.05))
#> Monte Carlo Confidence Intervals
#>             est     se     R  0.05%   0.5%   2.5%  97.5%  99.5% 99.95%
#> cp       0.2333 0.0264 20000 0.1465 0.1651 0.1817 0.2846 0.3012 0.3165
#> b        0.5082 0.0272 20000 0.4161 0.4385 0.4551 0.5611 0.5786 0.5988
#> a        0.4820 0.0264 20000 0.3959 0.4140 0.4299 0.5337 0.5495 0.5642
#> Y~~Y     0.5462 0.0244 20000 0.4660 0.4832 0.4979 0.5944 0.6095 0.6263
#> M~~M     0.7527 0.0339 20000 0.6441 0.6655 0.6858 0.8187 0.8387 0.8620
#> indirect 0.2449 0.0188 20000 0.1885 0.1983 0.2094 0.2831 0.2964 0.3105
#> direct   0.2333 0.0264 20000 0.1465 0.1651 0.1817 0.2846 0.3012 0.3165
#> total    0.4782 0.0264 20000 0.3911 0.4102 0.4262 0.5295 0.5457 0.5629

Standardized Monte Carlo Confidence Intervals

Standardized Monte Carlo Confidence intervals can be generated by passing the result of the MC() function to the MCStd() function.

Note: We recommend setting fixed.x = FALSE when generating standardized estimates and confidence intervals to model the variances and covariances of the exogenous observed variables if they are assumed to be random. If fixed.x = TRUE, which is the default setting in lavaan, MC() will fix the variances and the covariances of the exogenous observed variables to the sample values.

fit <- sem(data = df, model = model, fixed.x = FALSE)
unstd <- MC(fit, R = 20000L, alpha = c(0.001, 0.01, 0.05))

MCStd(unstd)
#> Standardized Monte Carlo Confidence Intervals
#>             est     se     R  0.05%   0.5%   2.5%  97.5%  99.5% 99.95%
#> cp       0.2422 0.0268 20000 0.1535 0.1722 0.1893 0.2945 0.3101 0.3294
#> b        0.5123 0.0246 20000 0.4289 0.4494 0.4639 0.5607 0.5764 0.5937
#> a        0.4963 0.0239 20000 0.4157 0.4338 0.4488 0.5427 0.5575 0.5754
#> Y~~Y     0.5558 0.0235 20000 0.4805 0.4951 0.5095 0.6018 0.6164 0.6334
#> M~~M     0.7537 0.0237 20000 0.6690 0.6892 0.7055 0.7986 0.8118 0.8272
#> X~~X     1.0000 0.0000 20000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000
#> indirect 0.2542 0.0176 20000 0.1974 0.2101 0.2206 0.2892 0.3005 0.3154
#> direct   0.2422 0.0268 20000 0.1535 0.1722 0.1893 0.2945 0.3101 0.3294
#> total    0.4964 0.0239 20000 0.4159 0.4342 0.4485 0.5418 0.5564 0.5712

Standardized Monte Carlo Confidence Intervals - An Alternative Approach

In this example, confidence intervals for the standardized indirect effect are generated by specifying the standardized indirect effect as a derived parameter using the := operator. The standardized indirect effect in a simple mediation model involves paths \(a\) and \(b\), and the standard deviations of \(X\) and \(Y\). It is given by

\[\begin{equation} a b \frac{s_X}{s_Y} \end{equation}\]

where

\[\begin{equation} s_X = \sqrt{ s_{X}^{2} } \end{equation}\]

and

\[\begin{equation} s_Y = \sqrt{ c^{\prime 2} s_{X}^{2} + a^2 b^2 s_{X}^{2} + b^2 s_{e_{M}}^{2} + 2 c^{\prime} b a s_{X}^{2} + s_{e_{Y}}^{2} } \end{equation}\]

where

\(s_{e_{Y}}^{2}\) and \(s_{e_{M}}^{2}\) are the residual variances in the regression equations.

The standardized indirect effect can be defined using the := operator and the named parameters in the model.

model <- "
  Y ~ cp * X + b * M
  M ~ a * X
  X ~~ s2_X * X
  M ~~ s2_em * M
  Y ~~ s2_ey * Y
  indirect_std := a * b * (
    sqrt(s2_X) / sqrt(
      (
        cp^2 * s2_X + a^2 * b^2 * s2_X
      ) + (
        b^2 * s2_em
      ) + (
        2 * cp * b * a * s2_X
      ) + s2_ey
    )
  )
"

fit <- sem(data = df, model = model, fixed.x = FALSE)

The row indirect_std corresponds to the confidence intervals for the standardized indirect effect.

MC(fit, R = 20000L, alpha = c(0.001, 0.01, 0.05))
#> Monte Carlo Confidence Intervals
#>                 est     se     R  0.05%   0.5%   2.5%  97.5%  99.5% 99.95%
#> cp           0.2333 0.0263 20000 0.1480 0.1648 0.1810 0.2838 0.3016 0.3204
#> b            0.5082 0.0270 20000 0.4154 0.4387 0.4553 0.5613 0.5769 0.5920
#> a            0.4820 0.0264 20000 0.3974 0.4155 0.4309 0.5338 0.5499 0.5733
#> s2_X         1.0590 0.0472 20000 0.9108 0.9390 0.9675 1.1531 1.1797 1.2133
#> s2_em        0.7527 0.0339 20000 0.6404 0.6655 0.6865 0.8192 0.8401 0.8665
#> s2_ey        0.5462 0.0243 20000 0.4661 0.4840 0.4987 0.5940 0.6086 0.6258
#> indirect_std 0.2542 0.0174 20000 0.2011 0.2104 0.2204 0.2884 0.3000 0.3128