Conditional Moderated Mediation in R
中山大學管理學院 卓雍然 2025-04-13
References:
Loading data & packages
Rows: 260
Columns: 7
$ PI <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
$ SD <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ SS <dbl> 1.3571, 1.8036, 3.0357, 4.4643, 2.5714, 2.7143, 2.2143, 3.0893, 2.0…
$ SA <dbl> 2.4286, 5.2857, 2.8571, 5.1429, 2.8571, 2.2857, 2.1429, 3.1429, 2.8…
$ L <dbl> 1.25, 6.25, 1.50, 3.00, 3.25, 2.00, 2.00, 3.75, 1.75, 4.75, 6.00, 1…
$ D <dbl> 1.0000, 1.0000, 2.0000, 3.0000, 1.3333, 2.0000, 1.0000, 2.6667, 3.0…
$ P <dbl> 1.0000, 1.6667, 4.3333, 3.0000, 2.0000, 2.6667, 2.0000, 2.3333, 6.3…【A】Path Analysis with lavvan
A1: Specify the Path Model with a Script
Model scripts of lavvan is very similar to those of
mplus.
For Example, when we have a mediation model with:
- Independent Var.:PI & SD
- Mediator:SS;
- Dependent Var.:D
We can write the model script in character string as below …
Dm1 = '
SS ~ a*PI + c*SD
D ~ d*PI + f*SS + g*SD
directPI := d
indirectPI := a*f
totalPI := d + (a*f)
directSD := g
indirectSD := c*f
totalSD := g + (c*f)
'A2: Fillting the model
Like
mplus, there're many optional argument inlavven<br> For detail, simply dohelp(“sem”)`
A3: Examine the Results
The result of the fitting is kept in the data object
fit. To examine the result, we simply call the
summary() function
lavaan 0.6-19 ended normally after 1 iteration
Estimator ML
Optimization method NLMINB
Number of model parameters 7
Number of observations 260
Model Test User Model:
Test statistic 0.000
Degrees of freedom 0
Model Test Baseline Model:
Test statistic 275.123
Degrees of freedom 5
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) -750.925
Loglikelihood unrestricted model (H1) -750.925
Akaike (AIC) 1515.849
Bayesian (BIC) 1540.774
Sample-size adjusted Bayesian (SABIC) 1518.581
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
SS ~
PI (a) 0.533 0.137 3.879 0.000 0.533 0.205
SD (c) 1.259 0.137 9.167 0.000 1.259 0.484
D ~
PI (d) 0.079 0.121 0.654 0.513 0.079 0.029
SS (f) 0.722 0.053 13.582 0.000 0.722 0.685
SD (g) 0.160 0.136 1.178 0.239 0.160 0.058
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.SS 1.226 0.108 11.402 0.000 1.226 0.724
.D 0.902 0.079 11.402 0.000 0.902 0.479
R-Square:
Estimate
SS 0.276
D 0.521
Defined Parameters:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
directPI 0.079 0.121 0.654 0.513 0.079 0.029
indirectPI 0.385 0.103 3.730 0.000 0.385 0.140
totalPI 0.464 0.154 3.013 0.003 0.464 0.169
directSD 0.160 0.136 1.178 0.239 0.160 0.058
indirectSD 0.910 0.120 7.598 0.000 0.910 0.332
totalSD 1.069 0.154 6.942 0.000 1.069 0.390We can extract the estimated coefficients & their CI95 by call
the parameterEstimates() function.
lhs op rhs label est se z pvalue ci.lower ci.upper
1 SS ~ PI a 0.533 0.137 3.879 0.000 0.264 0.802
2 SS ~ SD c 1.259 0.137 9.167 0.000 0.990 1.528
3 D ~ PI d 0.079 0.121 0.654 0.513 -0.158 0.317
4 D ~ SS f 0.722 0.053 13.582 0.000 0.618 0.827
5 D ~ SD g 0.160 0.136 1.178 0.239 -0.106 0.425
6 SS ~~ SS 1.226 0.108 11.402 0.000 1.015 1.437
7 D ~~ D 0.902 0.079 11.402 0.000 0.747 1.057
8 PI ~~ PI 0.250 0.000 NA NA 0.250 0.250
9 PI ~~ SD 0.000 0.000 NA NA 0.000 0.000
10 SD ~~ SD 0.250 0.000 NA NA 0.250 0.250
11 directPI := d directPI 0.079 0.121 0.654 0.513 -0.158 0.317
12 indirectPI := a*f indirectPI 0.385 0.103 3.730 0.000 0.183 0.587
13 totalPI := d+(a*f) totalPI 0.464 0.154 3.013 0.003 0.162 0.766
14 directSD := g directSD 0.160 0.136 1.178 0.239 -0.106 0.425
15 indirectSD := c*f indirectSD 0.910 0.120 7.598 0.000 0.675 1.144
16 totalSD := g+(c*f) totalSD 1.069 0.154 6.942 0.000 0.767 1.371A4: Re-Estimate the model with Bootstraping
Because the products of coefficients do not followed any known distributions, usually we need to re-estimate the CI’s by bootstrapping. It may takes a couple of minutes to do a bootstrap of 1000 repetitions.
The result of the bootstrap is stored in fit2.
lhs op rhs label est se z pvalue ci.lower ci.upper
1 SS ~ PI a 0.533 0.140 3.803 0.000 0.254 0.792
2 SS ~ SD c 1.259 0.138 9.094 0.000 0.982 1.520
3 D ~ PI d 0.079 0.120 0.661 0.509 -0.162 0.303
4 D ~ SS f 0.722 0.056 12.908 0.000 0.605 0.833
5 D ~ SD g 0.160 0.135 1.186 0.236 -0.111 0.413
6 SS ~~ SS 1.226 0.092 13.283 0.000 1.039 1.399
7 D ~~ D 0.902 0.081 11.203 0.000 0.734 1.052
8 PI ~~ PI 0.250 0.000 NA NA 0.250 0.250
9 PI ~~ SD 0.000 0.000 NA NA 0.000 0.000
10 SD ~~ SD 0.250 0.000 NA NA 0.250 0.250
11 directPI := d directPI 0.079 0.120 0.660 0.509 -0.162 0.303
12 indirectPI := a*f indirectPI 0.385 0.108 3.551 0.000 0.180 0.599
13 totalPI := d+(a*f) totalPI 0.464 0.152 3.057 0.002 0.171 0.756
14 directSD := g directSD 0.160 0.135 1.185 0.236 -0.111 0.413
15 indirectSD := c*f indirectSD 0.910 0.120 7.551 0.000 0.683 1.152
16 totalSD := g+(c*f) totalSD 1.069 0.148 7.239 0.000 0.792 1.367A5: Complare the CIs’
When we compare the CI’s in fit & fit2,
we’d see the CIs’ with bootstrap is usually narrower.
label est se z pvalue ci.lower ci.upper
11 directPI 0.079 0.121 0.654 0.513 -0.158 0.317
12 indirectPI 0.385 0.103 3.730 0.000 0.183 0.587
14 directSD 0.160 0.136 1.178 0.239 -0.106 0.425
15 indirectSD 0.910 0.120 7.598 0.000 0.675 1.144
111 directPI 0.079 0.120 0.660 0.509 -0.162 0.303
121 indirectPI 0.385 0.108 3.551 0.000 0.180 0.599
141 directSD 0.160 0.135 1.185 0.236 -0.111 0.413
151 indirectSD 0.910 0.120 7.551 0.000 0.683 1.152【B】Conditional Moderated Mediation
When the Mediation Path is Moderated by an Moderator
(SD),
- All of the direct, indirect and total effect coefficients depend on the value of the moderator.
- It takes a pretty length job script to specify an job of
Conditional Moderated Mediation. - R help to make the script easier and help to avoid mistakes.
B1: Make a Script
In the case of Conditional Moderated Mediation, the effects of the mediation paths are affected by a moderator. Usually moderation is implemented by introducing interaction terms into the model spec. For example, when we have a path diagram as below, …
Both the direct and indirect effect of PT to
D are moderated by SD via the interaction term
PI*SD.
- DIRECT Effect:
D ~ (d + e*SD)*PI + … - INDIRECT Effect:
D ~ (a + b*SD)*f*PI + …
A0$PISD = A0$PI * A0$SD # define the interaction term
Dm2 = '
SS ~ a*PI + b*PISD + c*SD
D ~ d*PI + e*PISD + f*SS + g*SD
'
Z = '
dirPI.D.SD% := d + %*e
indPI.D.SD% := a*f + %*b*f
totPI.D.SD% := d + %*e + a*f + %*b*f
'B2: Append the Conditional Lines
In order to evaluate how PI’s direct/indirect effects
varies with SD, we need to append the base model spec with
different value of SD. The code chuck below appends the
dir, ind & tot lines to the
script for each of
SD = 0, 0.125, 0.25, 0.375, 0.5, 0.625, 0.75, 0.875, 1
The Script becomes …
SS ~ a*PI + b*PISD + c*SD
D ~ d*PI + e*PISD + f*SS + g*SD
dirPI.D.SD0 := d + 0*e
indPI.D.SD0 := a*f + 0*b*f
totPI.D.SD0 := d + 0*e + a*f + 0*b*f
dirPI.D.SD0.125 := d + 0.125*e
indPI.D.SD0.125 := a*f + 0.125*b*f
totPI.D.SD0.125 := d + 0.125*e + a*f + 0.125*b*f
dirPI.D.SD0.25 := d + 0.25*e
indPI.D.SD0.25 := a*f + 0.25*b*f
totPI.D.SD0.25 := d + 0.25*e + a*f + 0.25*b*f
dirPI.D.SD0.375 := d + 0.375*e
indPI.D.SD0.375 := a*f + 0.375*b*f
totPI.D.SD0.375 := d + 0.375*e + a*f + 0.375*b*f
dirPI.D.SD0.5 := d + 0.5*e
indPI.D.SD0.5 := a*f + 0.5*b*f
totPI.D.SD0.5 := d + 0.5*e + a*f + 0.5*b*f
dirPI.D.SD0.625 := d + 0.625*e
indPI.D.SD0.625 := a*f + 0.625*b*f
totPI.D.SD0.625 := d + 0.625*e + a*f + 0.625*b*f
dirPI.D.SD0.75 := d + 0.75*e
indPI.D.SD0.75 := a*f + 0.75*b*f
totPI.D.SD0.75 := d + 0.75*e + a*f + 0.75*b*f
dirPI.D.SD0.875 := d + 0.875*e
indPI.D.SD0.875 := a*f + 0.875*b*f
totPI.D.SD0.875 := d + 0.875*e + a*f + 0.875*b*f
dirPI.D.SD1 := d + 1*e
indPI.D.SD1 := a*f + 1*b*f
totPI.D.SD1 := d + 1*e + a*f + 1*b*fas shown above, now we can evaluate the direct, indirect and total
effect of PI to D when
SD = 0, 0.125, 0.25, 0.375, 0.5, 0.625, 0.75, 0.875, 1
B3: Fitting the Model with Bootstrapping
The bootstrap process may takes a few minutes …
B4: Extract and Visualize the Result
Extract the Coefficient Estimation and CI95 ranges
x = summary(fitDm2,fit.measures=TRUE,standardized=TRUE,rsquare=TRUE)
p = parameterEstimates(fitDm2); p lhs op rhs label est se
1 SS ~ PI a 0.830 0.203
2 SS ~ PISD b -0.594 0.282
3 SS ~ SD c 1.556 0.196
4 D ~ PI d 0.110 0.143
5 D ~ PISD e -0.059 0.239
6 D ~ SS f 0.721 0.056
7 D ~ SD g 0.191 0.166
8 SS ~~ SS 1.204 0.097
9 D ~~ D 0.902 0.081
10 PI ~~ PI 0.250 0.000
11 PI ~~ PISD 0.125 0.000
12 PI ~~ SD 0.000 0.000
13 PISD ~~ PISD 0.188 0.000
14 PISD ~~ SD 0.125 0.000
15 SD ~~ SD 0.250 0.000
16 dirPI.D.SD0 := d+0*e dirPI.D.SD0 0.110 0.143
17 indPI.D.SD0 := a*f+0*b*f indPI.D.SD0 0.598 0.157
18 totPI.D.SD0 := d+0*e+a*f+0*b*f totPI.D.SD0 0.708 0.174
19 dirPI.D.SD0.125 := d+0.125*e dirPI.D.SD0.125 0.102 0.128
20 indPI.D.SD0.125 := a*f+0.125*b*f indPI.D.SD0.125 0.544 0.140
21 totPI.D.SD0.125 := d+0.125*e+a*f+0.125*b*f totPI.D.SD0.125 0.647 0.155
22 dirPI.D.SD0.25 := d+0.25*e dirPI.D.SD0.25 0.095 0.118
23 indPI.D.SD0.25 := a*f+0.25*b*f indPI.D.SD0.25 0.491 0.125
24 totPI.D.SD0.25 := d+0.25*e+a*f+0.25*b*f totPI.D.SD0.25 0.586 0.144
25 dirPI.D.SD0.375 := d+0.375*e dirPI.D.SD0.375 0.088 0.115
26 indPI.D.SD0.375 := a*f+0.375*b*f indPI.D.SD0.375 0.437 0.115
27 totPI.D.SD0.375 := d+0.375*e+a*f+0.375*b*f totPI.D.SD0.375 0.525 0.143
28 dirPI.D.SD0.5 := d+0.5*e dirPI.D.SD0.5 0.080 0.120
29 indPI.D.SD0.5 := a*f+0.5*b*f indPI.D.SD0.5 0.384 0.109
30 totPI.D.SD0.5 := d+0.5*e+a*f+0.5*b*f totPI.D.SD0.5 0.464 0.152
31 dirPI.D.SD0.625 := d+0.625*e dirPI.D.SD0.625 0.073 0.131
32 indPI.D.SD0.625 := a*f+0.625*b*f indPI.D.SD0.625 0.330 0.110
33 totPI.D.SD0.625 := d+0.625*e+a*f+0.625*b*f totPI.D.SD0.625 0.403 0.170
34 dirPI.D.SD0.75 := d+0.75*e dirPI.D.SD0.75 0.065 0.148
35 indPI.D.SD0.75 := a*f+0.75*b*f indPI.D.SD0.75 0.277 0.116
36 totPI.D.SD0.75 := d+0.75*e+a*f+0.75*b*f totPI.D.SD0.75 0.342 0.194
37 dirPI.D.SD0.875 := d+0.875*e dirPI.D.SD0.875 0.058 0.168
38 indPI.D.SD0.875 := a*f+0.875*b*f indPI.D.SD0.875 0.223 0.127
39 totPI.D.SD0.875 := d+0.875*e+a*f+0.875*b*f totPI.D.SD0.875 0.281 0.222
40 dirPI.D.SD1 := d+1*e dirPI.D.SD1 0.051 0.191
41 indPI.D.SD1 := a*f+1*b*f indPI.D.SD1 0.170 0.142
42 totPI.D.SD1 := d+1*e+a*f+1*b*f totPI.D.SD1 0.221 0.254
z pvalue ci.lower ci.upper
1 4.082 0.000 0.435 1.225
2 -2.105 0.035 -1.163 -0.055
3 7.919 0.000 1.180 1.943
4 0.766 0.444 -0.171 0.385
5 -0.248 0.804 -0.538 0.407
6 12.916 0.000 0.609 0.829
7 1.150 0.250 -0.136 0.512
8 12.374 0.000 1.003 1.388
9 11.191 0.000 0.734 1.048
10 NA NA 0.250 0.250
11 NA NA 0.125 0.125
12 NA NA 0.000 0.000
13 NA NA 0.188 0.188
14 NA NA 0.125 0.125
15 NA NA 0.250 0.250
16 0.766 0.444 -0.171 0.385
17 3.799 0.000 0.306 0.927
18 4.062 0.000 0.361 1.052
19 0.802 0.423 -0.152 0.346
20 3.889 0.000 0.290 0.837
21 4.177 0.000 0.348 0.960
22 0.807 0.420 -0.140 0.325
23 3.915 0.000 0.258 0.745
24 4.080 0.000 0.299 0.878
25 0.763 0.445 -0.137 0.321
26 3.811 0.000 0.222 0.668
27 3.683 0.000 0.244 0.815
28 0.671 0.502 -0.152 0.320
29 3.513 0.000 0.174 0.602
30 3.057 0.002 0.165 0.756
31 0.556 0.578 -0.192 0.329
32 3.013 0.003 0.121 0.546
33 2.375 0.018 0.066 0.719
34 0.443 0.658 -0.231 0.354
35 2.390 0.017 0.056 0.509
36 1.765 0.078 -0.051 0.710
37 0.345 0.730 -0.282 0.390
38 1.759 0.079 -0.024 0.475
39 1.265 0.206 -0.175 0.705
40 0.265 0.791 -0.329 0.430
41 1.196 0.232 -0.112 0.457
42 0.869 0.385 -0.288 0.703Arrange the data for visualiztion
df = filter(p, str_detect(label,"^ind")) %>% transmute(
SD = as.numeric(str_remove(lhs, "^.*SD")),
est,ci.lower,ci.upper); df SD est ci.lower ci.upper
1 0.000 0.598 0.306 0.927
2 0.125 0.544 0.290 0.837
3 0.250 0.491 0.258 0.745
4 0.375 0.437 0.222 0.668
5 0.500 0.384 0.174 0.602
6 0.625 0.330 0.121 0.546
7 0.750 0.277 0.056 0.509
8 0.875 0.223 -0.024 0.475
9 1.000 0.170 -0.112 0.457Visualize the result in an Interactive Chart
gg = gather(df, "CI95", "indirect_effect", 2:4) %>%
ggplot(aes(x=SD, y=indirect_effect, col=CI95)) +
geom_smooth(linewidth=1,se=F,method="loess") +
theme_bw() + geom_hline(yintercept=0) +
ggtitle("Conditional Indirect Effect of PI on D",
"mediated by SS which is moderated by SD") +
labs(x="SD")
ggplotly(gg) %>% layout(legend=list(
bgcolor='rgba(0,0,0,0)'),
font=list(size=0.5)
)🌻 With the interactive chart above, we can conclude that in the
confidence level of 95%, the PI’s mediation effect via
SS to D is significant in the condition that
SD < 0.85.
【C】 Multiple Conditional Moderated Mediation
Sometimes we need to examine multiple SEM model at once. For example
In the path diagram above
- we have two dependent variables (
D&P), and PI&SD’s indirect path to thes DVs are moderated by each other.
To cope with this complicated situation, we can write parameterized scripts in R. By feeding in different parameters, we can create multiple scripts at once.
C1. PROGRAMMABLE model script
MED = "SS ~ a*PI + b*PISD + c*SD\n"
OUT = paste(c('D','P'),"~ d*PI + e*PISD + f*SS + g*SD")
PATH = c("
dirPI.D.SD% := d + %*e
indPI.D.SD% := a*f + %*b*f
totPI.D.SD% := d + %*e + a*f + %*b*f
","
dirSD.D.PI% := g + %*e
indSD.D.PI% := c*f + %*b*f
totSD.D.PI% := g + %*e + c*f + %*b*f
")C2. Turn on Parallel Processing
The doParallel package let us making multiple scripts
and executing multiple scripts in parallel.
C3. Make & Execute Scripts in Parallel
In this code chuck we create fours scripts
- Independent Var.:PI; Mediator:SS; Moderator:SD; Dependent
Var.:D
- Independent Var.:SD; Mediator:SS; Moderator:PI; Dependent
Var.:D
- Independent Var.:PI; Mediator:SS; Moderator:SD; Dependent
Var.:P
- Independent Var.:SD; Mediator:SS; Moderator:PI; Dependent Var.:P
and bootstrap the four models in parallel. The resulting models will
be keep in the list X.
X = foreach(i = 1:2) %:% foreach(j = 1:2) %dopar% {
library(lavaan)
mod = paste(MED, OUT[i])
for(k in seq(0,1,0.125)) mod = paste0(mod, gsub("%",k,PATH[j]))
set.seed(123)
fit = sem(mod,A0,se="bootstrap",bootstrap=2000) }Stop parallel processing before proceeding.
C4. Extract the Coefficients and CI’s
Then we extract the Coefficients & CI’s from X and
keep them in a data frame df for visulazation.
C5. Multi-Panel Visulization
Now we can visualize all of the four conditional moderated mediation effects within a multi-panel interactive chart.
gg = gather(df, ".95CI", "indirect_effect", 4:6) %>%
ggplot(aes(x=moderator, y=indirect_effect, col=`.95CI`)) +
geom_smooth(linewidth=1,se=F,method="loess") +
theme_bw() + geom_hline(yintercept=0) +
facet_grid(output~input) +
ggtitle("PI/SD's Conditional Indirect Effect on D/P",
"mediated by SS which is modertaed by SD/PI") +
xlab("SD <Midiation Moderator> PI")
ggplotly(gg) %>% layout(legend=list(
bgcolor='rgba(0,0,0,0)'),
font=list(size=0.5))🌻 With the multi-panel interactive chart above, we can conclude that, in the confidence level of 95%,
PI’s indirect effect toD&Pare significant, only when the value ofSDis lower than~0.85.- Regardless the value of
PI,SD’s indirect effect toD&Pare always significant.