# Lavaan mediation + moderation + 2 X's

I'm trying to build a SEM that looks like the picture shown below: Where:

• $$X$$ and $$X^2$$ are the independent variables
• $$Y$$ is the dependent variable
• $$M_1$$ and $$M_2$$ are the mediators
• $$W$$ is a moderator

None of the default mediation packages that I've tried so far support such a structure so I'm using lavaan SEM instead. I've been reading up on Hayes and have figured out the model without the moderator:

  # direct effect
Y ~ cprime1*X + cprime2*X2 + controls

# mediator 1 paths
M1 ~ a1*X+a2*X2  + controls
Y ~ b1*M1

# mediator 2 paths
M2 ~ a1prime*X + a2prime*X2 + d*M1 + controls
Y ~ b2*M2

# indirect effect (a*b) using the Hayes2010 formula (Table 1)
ab :=  (a1 + 2*a2*X)*b1 + (a1prime + 2*a2prime + (a1+2*a2*X)*d ) * b2

# total effect
cprime := (cprime1+2*cprime2)
total := cprime + (ab)


I'm aware of how to moderate with 1 input variable, for instance for M1 you'd have something like:

M1 ~ a1*X1 + a2*W + a3*X1:W + controls


But it's unclear to me how one would add the moderators for X2 in these and the other parts of the model.

• I am quite unsure about your indirect effect. You are looking to 2 indirect effects and the equation only compute a single one.
– POC
Apr 2, 2021 at 16:51
• I'm following Hayes to compute the total indirect effect as the sum of the individual indirect effects.
– ciri
Apr 2, 2021 at 20:04

I simulated some irrelevant data:

library(lavaan)

dat <- as.data.frame(mvrnorm(5e2, rep(0, 6), matrix(.25, 6, 6) + .75 * diag(6)))
colnames(dat) <- c("X", "X2", "M1", "M2", "W", "Y")
dat$$X2 <- dat$$X ^ 2
#            X         X2         M1         M2          W          Y
# 1 -1.2556613 1.57668526  0.9558917 -0.6155703 -0.4540778 -0.3747208
# 2 -0.8463471 0.71630340  0.2066086 -0.1680590 -0.9181445 -0.2819832
# 3  0.1081230 0.01169059 -0.1934604  1.2662057  0.4817797 -0.5342619
# 4  0.6180336 0.38196555  0.3301925 -0.1277026  1.3222274  0.3626635


Given these data, this is the model code you want:

summary(sem(
"M1 ~ a11 * X + a12 * X2 + w10 * W + w11 * X:W + w12 * X2:W
M2 ~ a21 * X + a22 * X2 + w20 * W + w21 * X:W + w22 * X2:W + d * M1
Y ~ c1 * X + c2 * X2 + w30 * W + w31 * X:W + w32 * X2:W + b1 * M1 + b2 * M2
a11db2 := a11 * d * b2
a12db2 := a12 * d * b2
a11b1 := a11 * b1
a12b1 := a12 * b1
a21b2 := a21 * b2
a22b2 := a22 * b2",
dat), rsquare = TRUE)


With this, you can work out the exact combinations that get you the total indirect, direct and total effects. The nature of W makes things a bit difficult. Right now, the computed effects at the bottom are for when W = 0. You might want to play around with the computed effects depending on the type of variable W is.

• Thank you, this is very helpful! How would one get the effect size for one standard deviation above/below the value of a continuous W?
– ciri
Apr 9, 2021 at 1:08
• @ciri standardize W before adding it to the model. Then add its interaction coefficient wherever it might be of interest. Apr 9, 2021 at 2:32
• I think library(MASS) is required for running mvrnorm. Oct 18, 2022 at 21:39
• I have requested the developers of semopy to reproduce this example. Oct 18, 2022 at 23:47

Moderation is basically the product of the IVs ($$X_i$$) with the moderator ($$W$$). Since you have two IVs, you have two products $$X_1W$$ and $$W_2W$$, you just have to do it once for each predictor $$X_i$$ on each outcome ($$M_1,M_2,Y$$), since you want the moderator on each of them.

I used the notation in the graph for the code,

Y ~ cprime1*X + cprime2*X2 + w3_1*X1:W + w3_2*X2:W + controls
M2 ~ aprime1*X1 + aprime2*W + w2_1*X1:W + w2_2*X2:W + + d*M1 + controls
M1 ~ a1*X1 + a2*W + w2_1*X1:W + w2_2*X2:W + controls


You can also consider a triple interaction X1:X2:W. You would have to add the double interaction of $$X_1X_2$$ and the triple interaction $$X_1X_2W$$. At this point you would have an IV and two moderators (conceptually speaking; statistically, it is arbitrary because there is no difference between the variables). I don't think it is what you want though.