rggplot2linear-regression

Added variable / partial regression plots for groups in an interaction?

I'm looking to make added variable / partial regression plots for a regression with an interaction, which we can do with car.

# packages
library(car)
# linear model with 2-way interaction
lm1 <- lm(data = iris, Petal.Length ~ Sepal.Width + Sepal.Length*Species)
# av plot data
out <- avPlots(lm1)
# av plot of the interaction
avPlots(lm1, terms = "Sepal.Length:Species")

However, I'd like to plot the partial regression coefficient for each "Species" so that each has a separate regression line. To illustrate the idea, something like:

# ggplot idea to illustrate
library(ggplot2)
ggplot(data = iris, mapping = aes(x = Sepal.Width, y = Petal.Length,
 color = as.factor(Species))) + 
geom_smooth(method = "lm") + 
  labs(y = "Residuals from \nPetal Length ~ Sepal.Length + Species +
       Sepal.Length:Species",
       x = "Residuals from \nSepal Width ~
       Sepal.Length + Species + Sepal.Length:Species")

Is there a way to do this either through a package or by hand?

Solution

  • You could simply add the Species information to the partial plots output out returned by the call to avPlots(lm1) and run your ggplot() code on that data frame.

    For example, the partial plots for Sepal.Width by Species would be obtained by the following piece of code.

    library(ggplot2)

# Add the Species information to the appropriate `out` data returned by `avPlots(lm1)`
# Note: This assumes that the order of the samples in the original dataset
# is not modified by `avPlots()`, as is expected to be the case.
# Otherwise, use `merge()` to add the `Species` information to data frame `out`.
toplot = as.data.frame(out$Sepal.Width)
toplot$Species = iris$Species

# Create the partial plots by Species
ggplot(data = toplot, mapping = aes(x = Sepal.Width, y = Petal.Length,
color = as.factor(Species))) + 
geom_smooth(method="lm") +
  labs(y = "Residuals from \nPetal Length ~ Sepal.Length + Species +
       Sepal.Length:Species",
       x = "Residuals from \nSepal Width ~
       Sepal.Length + Species + Sepal.Length:Species")

    This is the output:

    Partial plots by Species