Two-Way-Interactions

Note: To better understand the principle of plotting interaction terms, it might be helpful to read the vignette on marginal effects first.

To plot marginal effects of interaction terms, at least two model terms need to be specified (the terms that define the interaction) in the terms-argument, for which the effects are computed. To plot marginal effects for three-way-interactions, all three terms need to be specified in terms.

A convenient way to automatically plot interactions is type = "int", which scans the model formula for interaction terms and then uses these as terms-argument.

library(sjPlot)
library(sjmisc)
library(ggplot2)
data(efc)
theme_set(theme_sjplot())

# make categorical
efc$c161sex <- to_factor(efc$c161sex)

# fit model with interaction
fit <- lm(neg_c_7 ~ c12hour + barthtot * c161sex, data = efc)

plot_model(fit, type = "pred", terms = c("barthtot", "c161sex"))

For type = "int", no terms need to be specified. Note that this plot type automatically uses the first interaction term in the formula for the x-axis, while the second term is used as grouping factor. Furthermore, if continuous variables are used as second term, you can specify preset-values for this term with the mdrt.values-argument, which are then used as grouping levels.

In this example, the second term is a factor with two levels (male/female), so there is no need for choosing specific values for the moderator.

plot_model(fit, type = "int")

To switch the terms, in this example barthtot and c161sex, simply switch the order of these terms on the terms-argument and use type = "pred".

plot_model(fit, type = "pred", terms = c("c161sex", "barthtot [0, 100]"))

To switch the terms for plot-type type = "int", you need to re-fit the model and change the formula accordingly, i.e. using c161sex as first term in the interaction.

# fit model with interaction, switching terms in formula
fit <- lm(neg_c_7 ~ c12hour + c161sex * barthtot, data = efc)
plot_model(fit, type = "int")

By default, for continuous variables, the minimum and maximum values are chosen as grouping levels, which are 0 and 100 - that’s why the previous two plots are identical. You have other options as well, e.g. the mean-value and +/- 1 standard deviation (as suggested by Cohen and Cohen for continuous variables and popularized by Aiken and West 1991), which can be specified using mdrt.values.

plot_model(fit, type = "int", mdrt.values = "meansd")

Three-Way-Interactions

Since the terms-argument accepts up to three model terms, you can also compute marginal effects for a 3-way-interaction.

# fit model with 3-way-interaction
fit <- lm(neg_c_7 ~ c12hour * barthtot * c161sex, data = efc)

# select only levels 30, 50 and 70 from continuous variable Barthel-Index
plot_model(fit, type = "pred", terms = c("c12hour", "barthtot [30,50,70]", "c161sex"))

Again, type = "int" will automatically plot the interaction terms, however, using mdrt.values = "minmax" as default - in this case, the “levels” 0 and 100 from continuous variable barthtot are chosen by default.

plot_model(fit, type = "int")
#> [[1]]

#> 
#> [[2]]

#> 
#> [[3]]

#> 
#> [[4]]

References

Aiken and West (1991). Multiple Regression: Testing and Interpreting Interactions.