By default, modsem() creates product indicators for you
based on the interaction specified in your model. Behind the scenes,
modsem() generates a total of 9 variables (product
indicators) that are used as the indicators for your latent product.
m1 <- '
# Outer Model
X =~ x1 + x2 + x3
Y =~ y1 + y2 + y3
Z =~ z1 + z2 + z3
# Inner model
Y ~ X + Z + X:Z
'
est1 <- modsem(m1, oneInt)
cat(est1$syntax)
#> X =~ x1
#> X =~ x2
#> X =~ x3
#> Y =~ y1
#> Y =~ y2
#> Y =~ y3
#> Z =~ z1
#> Z =~ z2
#> Z =~ z3
#> Y ~ X
#> Y ~ Z
#> Y ~ XZ
#> XZ =~ 1*x1z1
#> XZ =~ x2z1
#> XZ =~ x3z1
#> XZ =~ x1z2
#> XZ =~ x2z2
#> XZ =~ x3z2
#> XZ =~ x1z3
#> XZ =~ x2z3
#> XZ =~ x3z3
#> x1z1 ~~ 0*x2z2
#> x1z1 ~~ 0*x2z3
#> x1z1 ~~ 0*x3z2
#> x1z1 ~~ 0*x3z3
#> x1z2 ~~ 0*x2z1
#> x1z2 ~~ 0*x2z3
#> x1z2 ~~ 0*x3z1
#> x1z2 ~~ 0*x3z3
#> x1z3 ~~ 0*x2z1
#> x1z3 ~~ 0*x2z2
#> x1z3 ~~ 0*x3z1
#> x1z3 ~~ 0*x3z2
#> x2z1 ~~ 0*x3z2
#> x2z1 ~~ 0*x3z3
#> x2z2 ~~ 0*x3z1
#> x2z2 ~~ 0*x3z3
#> x2z3 ~~ 0*x3z1
#> x2z3 ~~ 0*x3z2
#> x1z1 ~~ x1z2
#> x1z1 ~~ x1z3
#> x1z1 ~~ x2z1
#> x1z1 ~~ x3z1
#> x1z2 ~~ x1z3
#> x1z2 ~~ x2z2
#> x1z2 ~~ x3z2
#> x1z3 ~~ x2z3
#> x1z3 ~~ x3z3
#> x2z1 ~~ x2z2
#> x2z1 ~~ x2z3
#> x2z1 ~~ x3z1
#> x2z2 ~~ x2z3
#> x2z2 ~~ x3z2
#> x2z3 ~~ x3z3
#> x3z1 ~~ x3z2
#> x3z1 ~~ x3z3
#> x3z2 ~~ x3z3While this is often sufficient, you might want more control over how
these indicators are created. In general, modsem() offers
two mechanisms for controlling the creation of product indicators: 1. By
specifying the measurement model for your latent product yourself. 2. By
using the mean() and sum() functions,
collectively known as parceling operations.
Specifying the Measurement Model
By default, modsem() creates all possible combinations
of product indicators. However, another common approach is to match the
indicators by order. For example, let’s say you have an interaction
between the latent variables X and Z:
X =~ x1 + x2 and Z =~ z1 + z2. By default, you
would get XZ =~ x1z1 + x1z2 + x2z1 + x2z2. If you prefer to
use the matching approach, you would expect
XZ =~ x1z1 + x2z2 instead. To achieve this, you can use the
match = TRUE argument.
m2 <- '
# Outer Model
X =~ x1 + x2
Y =~ y1 + y2
Z =~ z1 + z2
# Inner model
Y ~ X + Z + X:Z
'
est2 <- modsem(m2, oneInt, match = TRUE)
summary(est2)
#> modsem (version 1.0.4, approach = dblcent):
#> lavaan 0.6-19 ended normally after 41 iterations
#>
#> Estimator ML
#> Optimization method NLMINB
#> Number of model parameters 22
#>
#> Number of observations 2000
#>
#> Model Test User Model:
#>
#> Test statistic 11.355
#> Degrees of freedom 14
#> P-value (Chi-square) 0.658
#>
#> Parameter Estimates:
#>
#> Standard errors Standard
#> Information Expected
#> Information saturated (h1) model Structured
#>
#> Latent Variables:
#> Estimate Std.Err z-value P(>|z|)
#> X =~
#> x1 1.000
#> x2 0.819 0.021 38.127 0.000
#> Y =~
#> y1 1.000
#> y2 0.807 0.010 82.495 0.000
#> Z =~
#> z1 1.000
#> z2 0.836 0.024 35.392 0.000
#> XZ =~
#> x1z1 1.000
#> x2z2 0.645 0.024 26.904 0.000
#>
#> Regressions:
#> Estimate Std.Err z-value P(>|z|)
#> Y ~
#> X 0.688 0.029 23.366 0.000
#> Z 0.576 0.029 20.173 0.000
#> XZ 0.706 0.032 22.405 0.000
#>
#> Covariances:
#> Estimate Std.Err z-value P(>|z|)
#> .x1z1 ~~
#> .x2z2 0.000
#> X ~~
#> Z 0.202 0.025 8.182 0.000
#> XZ 0.003 0.026 0.119 0.905
#> Z ~~
#> XZ 0.042 0.026 1.621 0.105
#>
#> Variances:
#> Estimate Std.Err z-value P(>|z|)
#> .x1 0.179 0.022 8.029 0.000
#> .x2 0.151 0.015 9.956 0.000
#> .y1 0.184 0.021 8.577 0.000
#> .y2 0.136 0.014 9.663 0.000
#> .z1 0.197 0.025 7.802 0.000
#> .z2 0.138 0.018 7.831 0.000
#> .x1z1 0.319 0.035 9.141 0.000
#> .x2z2 0.244 0.016 15.369 0.000
#> X 0.962 0.042 23.120 0.000
#> .Y 0.964 0.042 23.110 0.000
#> Z 0.987 0.044 22.260 0.000
#> XZ 1.041 0.054 19.441 0.000More Complicated Models
If you want even more control, you can use the
get_pi_syntax() and get_pi_data() functions to
extract the modified syntax and data from modsem(),
allowing you to modify them as needed. This can be particularly useful
in cases where you want to estimate a model using a feature in
lavaan that isn’t available in modsem().
For example, the syntax for ordered and multigroup models (as of now)
isn’t as flexible in modsem() as it is in
lavaan. You can modify the auto-generated syntax (along
with the altered dataset) from modsem() to suit your
needs.
m3 <- '
# Outer Model
X =~ x1 + x2
Y =~ y1 + y2
Z =~ z1 + z2
# Inner model
Y ~ X + Z + X:Z
'
syntax <- get_pi_syntax(m3)
cat(syntax)
#> X =~ x1
#> X =~ x2
#> Y =~ y1
#> Y =~ y2
#> Z =~ z1
#> Z =~ z2
#> Y ~ X
#> Y ~ Z
#> Y ~ XZ
#> XZ =~ 1*x1z1
#> XZ =~ x2z1
#> XZ =~ x1z2
#> XZ =~ x2z2
#> x1z1 ~~ 0*x2z2
#> x1z2 ~~ 0*x2z1
#> x1z1 ~~ x1z2
#> x1z1 ~~ x2z1
#> x1z2 ~~ x2z2
#> x2z1 ~~ x2z2
data <- get_pi_data(m3, oneInt)
head(data)
#> x1 x2 y1 y2 z1 z2 x1z1
#> 1 2.4345722 1.3578655 1.4526897 0.9560888 0.8184825 1.60708140 -0.4823019
#> 2 0.2472734 0.2723201 0.5496756 0.7115311 3.6649148 2.60983102 -2.2680403
#> 3 -1.3647759 -0.5628205 -0.9835467 -0.6697747 1.7249386 2.10981827 -1.9137416
#> 4 3.0432836 2.2153763 6.4641465 4.7805981 2.5697116 3.26335379 2.9385205
#> 5 2.8148841 2.7029616 2.2860280 2.1457643 0.3467850 0.07164577 -1.4009548
#> 6 -0.5453450 -0.7530642 1.1294876 1.1998472 -0.2362958 0.60252657 1.7465860
#> x2z1 x1z2 x2z2
#> 1 -0.1884837 0.3929380 -0.0730934
#> 2 -2.6637694 -1.2630544 -1.4547433
#> 3 -1.4299711 -2.3329864 -1.7383407
#> 4 1.3971422 3.9837389 1.9273102
#> 5 -1.1495704 -2.2058995 -1.8169042
#> 6 2.2950753 0.7717365 1.0568143