sjp.int: Plot interaction effects of (generalized) linear (mixed) models

Description

Plot regression (predicted values) or probability lines (predicted probabilities) of significant interaction terms to better understand effects of moderations in regression models. This function accepts following fitted model classes:

linear models (lm)
generalized linear models (glm)
linear mixed effects models (lmer)
generalized linear mixed effects models (glmer)
non-linear mixed effects models (nlmer)
linear mixed effects models (lme, but only fortype = "eff")
generalized least squares models (gls, but only fortype = "eff")
panel data estimators (plm)

Note that beside interaction terms, also the single predictors of each interaction (main effects) must be included in the fitted model as well. Thus, lm(dep ~ pred1 * pred2) will work, but lm(dep ~ pred1:pred2) won't!

Usage

sjp.int(fit, type = c("eff", "cond", "emm"), int.term = NULL,
  int.plot.index = NULL, mdrt.values = c("minmax", "meansd", "zeromax",
  "quart", "all"), swap.pred = FALSE, plevel = 0.1, diff = FALSE,
  title = NULL, axis.title = NULL, axis.labels = NULL,
  legend.title = NULL, legend.labels = NULL, wrap.title = 50,
  wrap.legend.labels = 20, wrap.legend.title = 20, geom.colors = "Set1",
  geom.size = NULL, fill.color = "grey", fill.alpha = 0.3,
  show.values = FALSE, show.ci = FALSE, p.kr = TRUE, grid.breaks = NULL,
  xlim = NULL, ylim = NULL, y.offset = 0.07, digits = 2,
  facet.grid = FALSE, prnt.plot = TRUE, ...)

Arguments

Value

(Insisibily) returns the ggplot-objects with the complete plot-list (plot.list) as well as the data frames that were used for setting up the ggplot-objects (data.list).

Details

[object Object],[object Object],[object Object] The argument int.term only applies to type = "eff" and can be used to select a specific interaction term of the model that should be plotted. The function then calls effect(int.term, fit) to compute effects for this specific interaction term only. This approach is recommended, when the fitted model contains many observations and/or variables, which may slow down the effect-computation dramatically. In such cases, consider computing effects for selected interaction terms only with int.terms. See 'Examples'.

References

Aiken and West (1991). Multiple Regression: Testing and Interpreting Interactions.
Brambor T, Clark WR and Golder M (2006) Understanding Interaction Models: Improving Empirical Analyses. Political Analysis 14: 63-82https://files.nyu.edu/mrg217/public/pa_final.pdf{download}
Esarey J, Sumner JL (2015) Marginal Effects in Interaction Models: Determining and Controlling the False Positive Rate.http://jee3.web.rice.edu/interaction-overconfidence.pdf{download}
Fox J (2003) Effect displays in R for generalised linear models. Journal of Statistical Software 8:15, 1–27,http://www.jstatsoft.org/v08/i15/{}
Hayes AF (2012) PROCESS: A versatile computational tool for observed variable mediation, moderation, and conditional process modeling [White paper]http://imaging.mrc-cbu.cam.ac.uk/statswiki/FAQ/SobelTest?action=AttachFile&do=get&target=process.pdf{download}
http://www.theanalysisfactor.com/interpreting-interactions-in-regression/{Grace-Martin K: Interpreting Interactions in Regression}
http://www.theanalysisfactor.com/clarifications-on-interpreting-interactions-in-regression/{Grace-Martin K: Clarifications on Interpreting Interactions in Regression}

Examples

Run this code

# Note that the data sets used in this example may not be perfectly suitable for
# fitting linear models. I just used them because they are part of the R-software.

# fit "dummy" model. Note that moderator should enter
# first the model, followed by predictor. Else, use
# argument "swap.pred" to change predictor on 
# x-axis with moderator
fit <- lm(weight ~ Diet * Time, data = ChickWeight)

# show summary to see significant interactions
summary(fit)

# plot regression line of interaction terms, including value labels
sjp.int(fit, type = "eff", show.values = TRUE)


# load sample data set
library(sjmisc)
data(efc)
# create data frame with variables that should be included
# in the model
mydf <- data.frame(usage = efc$tot_sc_e,
                   sex = efc$c161sex,
                   education = efc$c172code,
                   burden = efc$neg_c_7,
                   dependency = efc$e42dep)
# convert gender predictor to factor
mydf$sex <- relevel(factor(mydf$sex), ref = "2")
# fit "dummy" model
fit <- lm(usage ~ .*., data = mydf)
summary(fit)

# plot interactions. note that type = "cond" only considers 
# significant interactions by default. use "plevel" to 
# adjust p-level sensivity
sjp.int(fit, type = "cond")

# plot only selected interaction term for
# type = "eff"
sjp.int(fit, type = "eff", int.term = "sex*education")

# plot interactions, using mean and sd as moderator
# values to calculate interaction effect
sjp.int(fit, type = "eff", mdrt.values = "meansd")
sjp.int(fit, type = "cond", mdrt.values = "meansd")

# plot interactions, including those with p-value up to 0.1
sjp.int(fit, type = "cond", plevel = 0.1)

# -------------------------------
# Predictors for negative impact of care.
# Data from the EUROFAMCARE sample dataset
# -------------------------------
library(sjmisc)
data(efc)
# create binary response
y <- ifelse(efc$neg_c_7 < median(stats::na.omit(efc$neg_c_7)), 0, 1)
# create data frame for fitted model
mydf <- data.frame(y = as.factor(y),
                   sex = as.factor(efc$c161sex),
                   barthel = as.numeric(efc$barthtot))
# fit model
fit <- glm(y ~ sex * barthel, data = mydf, family = binomial(link = "logit"))
# plot interaction, increase p-level sensivity
sjp.int(fit, type = "eff", legend.labels = get_labels(efc$c161sex), plevel = 0.1)
sjp.int(fit, type = "cond", legend.labels = get_labels(efc$c161sex), plevel = 0.1)
        
# -------------------------------
# Plot estimated marginal means
# -------------------------------
# load sample data set
library(sjmisc)
data(efc)
# create data frame with variables that should be included
# in the model
mydf <- data.frame(burden = efc$neg_c_7,
                   sex = efc$c161sex,
                   education = efc$c172code)
# convert gender predictor to factor
mydf$sex <- factor(mydf$sex)
mydf$education <- factor(mydf$education)
# name factor levels and dependent variable
levels(mydf$sex) <- c("female", "male")
levels(mydf$education) <- c("low", "mid", "high")
mydf$burden <- set_label(mydf$burden, "care burden")
# fit "dummy" model
fit <- lm(burden ~ .*., data = mydf)
summary(fit)

# plot marginal means of interactions, no interaction found
sjp.int(fit, type = "emm")
# plot marginal means of interactions, including those with p-value up to 1
sjp.int(fit, type = "emm", plevel = 1)
# swap predictors
sjp.int(fit, type = "emm", plevel = 1, swap.pred = TRUE)

# -------------------------------
# Plot effects
# -------------------------------
# add continuous variable
mydf$barthel <- efc$barthtot
# re-fit model with continuous variable
fit <- lm(burden ~ .*., data = mydf)

# plot effects
sjp.int(fit, type = "eff", show.ci = TRUE)

# plot effects, faceted
sjp.int(fit, type = "eff", int.plot.index = 3, show.ci = TRUE, facet.grid = TRUE)

Run the code above in your browser using DataLab