I have cross-sectional data and I am trying to specify a model with multiple mediations.
My independent variable (IV) is measured by a tool with 24 items, which make up 5 subscales (latent variables), which in turn load onto a total "higher-order" factor. I try to estimate this model in two different ways: (1) A five-factor model (without a higher order factor) in which all 5 subscales are allowed to correlate and (2) a higher-order model with a TOTAL latent variable made up of those 5 suscales.The first model has five correlated latents (FNR, FOB...FAA) factors, with variance fixed to 1.
My first model has an IV (FTOTAL), 4 mediators (ER 1, 2, 3 and 4), and one DV (PH). The relationship between ER 4 and the IV is mediated by one of the other mediators (ER3), making it a mediated mediation. I was able to specify the SEM for the first model with the higher order TOTAL factor without a problem, and the code is shown below:
"#Measurements model
FNR =~ FNR1 + FNR2 + FNR3 +FNR4 +FNR5
FOB =~ FOB1 + FOB2 +FOB3 +FOB4
FDS =~ FDS1 +FDS2 +FDS3 + FDS4 + FDS5
FNJ =~ FNJ1 + FNJ2 + FNJ3 +FNJ4 + FNJ5
FAA =~ FAA1 + FAA2 +FAA3 + FAA4 +FAA5
FTOTAL =~ FNR + FOB + FDS + FNJ+ FAA
#Regressions
ER3~ a*FTOTAL
ER4~ b*RSTOTAL +FTOTAL
ER1 ~ u1*FTOTAL
ER2 ~ u2*FTOTAL
PHQTOTAL ~ z1*ER1 + z2*ER2 + d*FTOTAL + c*ER4 + ER3
indirect1 := u1 * z1
indirect2 := u2 * z2
indirect3 := a*b*c
total := d + (u1 * z1) + (u2 * z2) + a*b*c
#Residual correlations
CRTOTAL~~SUPTOTAL+SLEEPTOTAL+RSTOTAL
SUPTOTAL~~SLEEPTOTAL+RSTOTAL
"
fitPHtotal <- sem(model = multipleMediationPH, data = SEMDATA, std.lv=TRUE)
summary(fitPHtotal)
However, I cannot figure out how to specify the model with all the 5 subscales as independent variables in the same model. I tried following the same logic and but including the 5 IVs in the model and it did not work out. Could anyone suggest a solution? or is the only way to just run 5 different models with one subscale as an IV at a time?
Thank you in advance for the help.
Since you did not provide data. I show you how to do it using the HolzingerSwineford1939 which comes with the lavaan package.
First, a mediation using a second-order latent factor (3 first-order factors):
library(lavaan)
#> This is lavaan 0.6-8
#> lavaan is FREE software! Please report any bugs.
model_2L <- "
visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
higher =~ visual + textual + speed
#grade will be your Y
#higher order latent factor will be your X
#agemo will be your M
grade ~ c*higher + b*agemo
agemo ~ a*higher
# indirect effect (a*b)
ab := a*b
# total effect
total := c + (a*b)
"
fit_2L <- sem(model = model_2L, data = HolzingerSwineford1939)
summary(object = fit_2L, std=T)
#> lavaan 0.6-8 ended normally after 48 iterations
#>
#> Estimator ML
#> Optimization method NLMINB
#> Number of model parameters 26
#>
#> Used Total
#> Number of observations 300 301
#>
#> Model Test User Model:
#>
#> Test statistic 116.110
#> Degrees of freedom 40
#> P-value (Chi-square) 0.000
#>
#> Parameter Estimates:
#>
#> Standard errors Standard
#> Information Expected
#> Information saturated (h1) model Structured
#>
#> Latent Variables:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> visual =~
#> x1 1.000 0.849 0.727
#> x2 0.621 0.109 5.680 0.000 0.527 0.448
#> x3 0.824 0.124 6.641 0.000 0.699 0.619
#> textual =~
#> x4 1.000 0.990 0.851
#> x5 1.117 0.066 16.998 0.000 1.106 0.859
#> x6 0.922 0.056 16.563 0.000 0.913 0.834
#> speed =~
#> x7 1.000 0.648 0.595
#> x8 1.130 0.148 7.612 0.000 0.732 0.726
#> x9 1.010 0.135 7.465 0.000 0.655 0.649
#> higher =~
#> visual 1.000 0.673 0.673
#> textual 0.849 0.185 4.586 0.000 0.490 0.490
#> speed 0.810 0.179 4.519 0.000 0.714 0.714
#>
#> Regressions:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> grade ~
#> higher (c) 0.421 0.089 4.730 0.000 0.241 0.482
#> agemo (b) -0.004 0.008 -0.519 0.604 -0.004 -0.029
#> agemo ~
#> higher (a) 0.322 0.469 0.687 0.492 0.184 0.053
#>
#> Variances:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> .x1 0.641 0.110 5.822 0.000 0.641 0.471
#> .x2 1.108 0.102 10.848 0.000 1.108 0.799
#> .x3 0.786 0.094 8.398 0.000 0.786 0.616
#> .x4 0.373 0.048 7.750 0.000 0.373 0.276
#> .x5 0.436 0.058 7.453 0.000 0.436 0.263
#> .x6 0.364 0.044 8.369 0.000 0.364 0.304
#> .x7 0.767 0.080 9.629 0.000 0.767 0.646
#> .x8 0.482 0.070 6.924 0.000 0.482 0.474
#> .x9 0.589 0.068 8.686 0.000 0.589 0.579
#> .grade 0.192 0.020 9.767 0.000 0.192 0.768
#> .agemo 11.881 0.972 12.220 0.000 11.881 0.997
#> .visual 0.394 0.111 3.535 0.000 0.547 0.547
#> .textual 0.745 0.101 7.397 0.000 0.760 0.760
#> .speed 0.206 0.062 3.312 0.001 0.490 0.490
#> higher 0.327 0.097 3.375 0.001 1.000 1.000
#>
#> Defined Parameters:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> ab -0.001 0.004 -0.366 0.715 -0.001 -0.002
#> total 0.420 0.089 4.728 0.000 0.240 0.481
Second, a mediation using a three first-order factors. Three indirect effects and three total effects are estimated:
library(lavaan)
model_1L <- "
visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
#grade will be your Y
#higher order latent factor will be your X
#agemo will be your M
grade ~ c1*visual + c2*textual + c3*speed + b*agemo
agemo ~ a1*visual + a2*textual + a3*speed
# indirect effect (a*b)
a1b := a1*b
a2b := a2*b
a3b := a3*b
# total effect
total1 := c1 + (a1*b)
total2 := c2 + (a2*b)
total3 := c3 + (a3*b)
"
fit_1L <- sem(model = model_1L, data = HolzingerSwineford1939)
summary(object = fit_1L, std=T)
#> lavaan 0.6-8 ended normally after 55 iterations
#>
#> Estimator ML
#> Optimization method NLMINB
#> Number of model parameters 30
#>
#> Used Total
#> Number of observations 300 301
#>
#> Model Test User Model:
#>
#> Test statistic 101.925
#> Degrees of freedom 36
#> P-value (Chi-square) 0.000
#>
#> Parameter Estimates:
#>
#> Standard errors Standard
#> Information Expected
#> Information saturated (h1) model Structured
#>
#> Latent Variables:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> visual =~
#> x1 1.000 0.904 0.775
#> x2 0.555 0.100 5.564 0.000 0.501 0.426
#> x3 0.724 0.109 6.657 0.000 0.655 0.580
#> textual =~
#> x4 1.000 0.993 0.853
#> x5 1.108 0.065 17.017 0.000 1.101 0.855
#> x6 0.921 0.055 16.667 0.000 0.915 0.836
#> speed =~
#> x7 1.000 0.668 0.613
#> x8 1.115 0.142 7.840 0.000 0.744 0.737
#> x9 0.945 0.125 7.540 0.000 0.631 0.625
#>
#> Regressions:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> grade ~
#> visual (c1) 0.012 0.048 0.246 0.806 0.011 0.021
#> textual (c2) 0.048 0.035 1.376 0.169 0.047 0.095
#> speed (c3) 0.295 0.063 4.689 0.000 0.197 0.394
#> agemo (b) -0.003 0.008 -0.361 0.718 -0.003 -0.020
#> agemo ~
#> visual (a1) 0.354 0.355 0.996 0.319 0.320 0.093
#> textual (a2) -0.233 0.256 -0.912 0.362 -0.231 -0.067
#> speed (a3) 0.098 0.421 0.232 0.817 0.065 0.019
#>
#> Covariances:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> visual ~~
#> textual 0.412 0.074 5.565 0.000 0.459 0.459
#> speed 0.265 0.058 4.554 0.000 0.438 0.438
#> textual ~~
#> speed 0.180 0.052 3.448 0.001 0.271 0.271
#>
#> Variances:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> .x1 0.545 0.115 4.747 0.000 0.545 0.400
#> .x2 1.135 0.102 11.115 0.000 1.135 0.819
#> .x3 0.846 0.091 9.322 0.000 0.846 0.664
#> .x4 0.368 0.048 7.698 0.000 0.368 0.272
#> .x5 0.447 0.058 7.657 0.000 0.447 0.270
#> .x6 0.361 0.043 8.343 0.000 0.361 0.301
#> .x7 0.741 0.079 9.422 0.000 0.741 0.624
#> .x8 0.465 0.069 6.724 0.000 0.465 0.456
#> .x9 0.620 0.067 9.217 0.000 0.620 0.609
#> .grade 0.201 0.018 11.307 0.000 0.201 0.806
#> .agemo 11.813 0.969 12.191 0.000 11.813 0.991
#> visual 0.817 0.147 5.564 0.000 1.000 1.000
#> textual 0.986 0.113 8.752 0.000 1.000 1.000
#> speed 0.446 0.091 4.906 0.000 1.000 1.000
#>
#> Defined Parameters:
#> Estimate Std.Err z-value P(>|z|) Std.lv Std.all
#> a1b -0.001 0.003 -0.344 0.731 -0.001 -0.002
#> a2b 0.001 0.002 0.335 0.738 0.001 0.001
#> a3b -0.000 0.002 -0.183 0.855 -0.000 -0.000
#> total1 0.011 0.048 0.226 0.821 0.010 0.020
#> total2 0.048 0.035 1.399 0.162 0.048 0.096
#> total3 0.295 0.063 4.685 0.000 0.197 0.394
Created on 2021-03-30 by the reprex package (v1.0.0)