sjp.lmer {sjPlot}

This document shows examples for using the sjp.lmer function of the sjPlot package.


(back to table of content)

Data initialization

Please refer to this document.

sjp.setTheme(theme = "forestgrey", 
             geom.label.size = 3.5, 
             axis.textsize = .9, 
             axis.title.size = .9)

Customizing plot appearance

Please refer to this document

Plotting random effects of linear mixed effects models

The sjp.lmer function plots effects of merMod objects, which were fitted using the lmer function of the lme4 package.

First, we need fit a sample model.

# fit model
fit <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)

By default, random effects are plotted. In this example, the random effects of random intercept and random coefficient(s) are plotted as an integrated (faceted plot.) Note that the y.offset argument is used to adjust the value label position. Depending on text size and screen resolution, the default position of text labels may vary.

sjp.lmer(fit, y.offset = .4)

With the sort.est argument you can specify the random intercept or any random coefficient in order to sort the effects accordingly.

# sort by predictor Days
sjp.lmer(fit, sort.est = "Days", y.offset = .4)

Since the grouping levels in a faceted plot define the x-axis for each facet, sorting can only be applied to one coefficient or the intercept. If all random intercepts and random coefficients should be sorted, turn faceting off and use sort.est = "sort.all".

# sort all predictors
         facet.grid = FALSE,
         sort.est = "sort.all",
         y.offset = .4)

Plotting fixed effects of linear mixed effects models

For the next examples, we fit another model.

# for sample data and tool functions
# prepare group variable
efc$grp = as.factor(efc$e15relat)
levels(x = efc$grp) <- get_labels(efc$e15relat)
# data frame for fitted model
mydf <- data.frame(neg_c_7 = efc$neg_c_7,
                   sex = to_factor(efc$c161sex),
                   c12hour = efc$c12hour,
                   barthel = efc$barthtot,
                   grp = efc$grp)
# fit glmer
fit2 <- lmer(neg_c_7 ~ sex + c12hour + barthel + (1 | grp),
             data = mydf)

With the type argument, you can decide which effects to plot. type = "fe" plots the fixed effects. In this case, the value label for the negative estimate is not completely visible in the output, due to tight axis limits. Use axis.lim to manually change the axis limits.

# plot fixed effects
sjp.lmer(fit2, type = "fe", axis.lim = c(-5, 15))

Plotting standardized coefficients of fixed effects

Using type = "fe.std", we can plot the standardized estimates of the fixed effects.

# plot standardized fixed effects
sjp.lmer(fit2, type = "fe.std")

Plotting slopes of fixed effects

With type = "fe.slope", slopes for each fixed coefficient (i.e. each predictor regressed on the response) are plotted. Specific predictors can be selected with vars or removed from the plot with the remove.estimates arguments.

# plot fixed effects slopes
sjp.lmer(fit2, type = "fe.slope", vars = c("c12hour", "barthel"))

Effect plots

type = "eff" produces effect plots of the mixed model. The Effect function from the effects package is used to compute the data.

# plot effects
sjp.lmer(fit2, type = "eff")

Predicted Values

type = "pred" or type = "pred.fe" plot predicted values for response, conditional on fixed effects only or on random intercept. It’s calling predict(fit, type = "response", re.form = NA) resp. predict(fit, type = "response", re.form = NULL) to compute the values. This plot type requires the vars argument to select specific terms that should be used for the x-axis and - optional - as grouping factor. Hence, vars must be a character vector with the names of one or two model predictors.

# plot effects
sjp.lmer(fit2, type = "pred", vars = "c12hour")

Or grouped by gender:

# plot effects
sjp.lmer(fit2, type = "pred", vars = c("c12hour", "sex"))

You can also plot predictions by groups in a non-facted plot:

# plot effects
sjp.lmer(fit2, type = "pred", facet.grid = FALSE,
         vars = c("c12hour", "sex"))

Plotting fixed effects slopes for each random intercept (group levels)

To get a better picture of the linear relationship between fixed effects and response depending on the grouping levels (random intercepts), you can plot straight slope lines (ablines) for each coefficient with varying random intercepts.

Basically, the formula is b0 + b0[r1-rn] + bi * xi (where xi is the estimate of fixed effects, b0 is the intercept of the fixed effects and b0[r1-rn] are all random intercepts).

Use type = "ri.slope" for this kind of plots.

# random intercepts
## $grp
##                         (Intercept)
## spouse/partner           0.62318882
## child                    0.42390010
## sibling                 -0.05435509
## daughter or son -in-law  0.05759517
## ancle/aunt              -0.10111187
## nephew/niece            -0.55060420
## cousin                  -0.11383605
## other, specify          -0.28477688
# fixed effects
##  (Intercept)         sex2      c12hour      barthel 
## 14.135947301  0.478500190  0.003365667 -0.047946653
# plot fixed effects depending on group levels
sjp.lmer(fit2, type = "ri.slope")

In non-faceted plots, grouping levels might be difficult to distinguish. However, you can emphasize specific groups with emph.grp. Note that emph.grp only works in non-faceted plots (i.e. facet.grid = FALSE)! Remaining (non-emphasized) groups have a light grey color.

# plot fixed effects depending on group levels
# emphasize group levels 1, 2 and 5
         type = "ri.slope", 
         emph.grp = c(1, 2, 5), 
         facet.grid = FALSE)

Grouping levels can also be specified via their names.

# plot fixed effects depending on group levels
# emphasize groups child and cousin
# only for predictor barthel
         type = "ri.slope", 
         emph.grp = c("child", "cousin"), 
         facet.grid = FALSE,
         vars = "barthel")

Plotting random slopes depending on random intercepts

Use this plot type to visualize the random parts of random slope-intercept (or repeated measure) models. When having too many groups, use the sample.n argument to randomly select a specific amount of subjects.

# plot random-slope-intercept
sjp.lmer(fit, type = "rs.ri", sample.n = 15)

If sample.n is a vector of length greater than one, the specifc “subjects” indicated by sample.n are plotted.

# plot random-slope-intercept, plot subjects 1, 5 and 7.
sjp.lmer(fit, type = "rs.ri", 
         sample.n = c(1, 5, 7),
         show.legend = TRUE)

Plotting selected fixed effects slopes for each random intercept

If you want to plot only selected coefficients for each random intercept, use the vars argument and specify the coefficient’s name.

# plot fixed effect "c12hour" only,
# depending on group levels
         type = "ri.slope",
         vars = "c12hour")

Diagnostic plots

Currently, there are two type options to plot diagnostic plots.

Correlation matrix of fixed effects

To plot a correlation matrix of the fixed effects, use type = "fe.cor"

# plot fixed effects correlation matrix
sjp.lmer(fit2, type = "fe.cor")

qq-plot of random effects

Another diagnostic plot is the qq-plot for random effects. Use type = "re.qq" to plot random against standard quantiles. The dots should be plotted along the line.

# plot qq-plot of random effects
sjp.lmer(fit2, type = "re.qq")

If you have other random effects, like random coefficients, qq-plots for these effects are plotted as well. We refer to the first model to demonstrate this.

# plot qq-plot of random effects
sjp.lmer(fit, type = "re.qq")