I would like to plot the marginal effect of a standard deviation change for an interaction term. I am trying to plot a conditional marginal effect of an interaction term, the partial derivative (slope) of the model, but due to the nature and variation of my data, a unit change could be a bit misleading. Thus, I would like to transform the figure into a standard deviation change instead of a unit change, so I can get the substantive effect for better interpretation. I am currently using the marginaleffects package to plot it. So, here is an example.
Suppose I have the following model, and I would like to see the marginal effect of hp on mpg across values of disp (just as example, I am not sure what these variables mean):
data(mtcars)
mod <- lm(mpg ~ hp * disp + factor(am), data = mtcars)
summary(mod)
Call:
lm(formula = mpg ~ hp * disp + factor(am), data = mtcars)
Residuals:
Min 1Q Median 3Q Max
-3.3599 -1.7385 -0.2522 0.7915 5.2211
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.600e+01 3.678e+00 9.786 2.25e-10 ***
hp -9.196e-02 2.441e-02 -3.767 0.000816 ***
disp -5.418e-02 1.858e-02 -2.915 0.007062 **
factor(am)1 2.289e+00 1.455e+00 1.573 0.127297
hp:disp 2.269e-04 9.377e-05 2.419 0.022563 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.623 on 27 degrees of freedom
Multiple R-squared: 0.835, Adjusted R-squared: 0.8105
F-statistic: 34.15 on 4 and 27 DF, p-value: 3.355e-10
To see the conditional marginal effect, I use the marginaleffects package to plot the figure:
library(marginaleffects)
plot1 <- plot_slopes(mod, variables = "hp", condition = "disp") +
geom_hline(yintercept=0,linetype=2) +
geom_rug(aes(x=disp), data=mtcars, sides="b") +
theme_grey()
plot1
Which result in the following figure:
Which displays the marginal effect (slope) of a unit change of hp across disp. Instead of a unit change, I would like to get the marginal effect of a standard deviation change of hp. I am open to use other packages and alternative ways to plot the marginal effect. Not sure how malleable the marginaleffects package is. Let me know if you need further clarification.
The obvious solution is to scale
the variable of interest (in this case hp
), which is equivalent to (hp - mean(hp))/sd(hp)
. That way, a 1-unit change will be a 1-sd change.
library(marginaleffects)
library(ggplot2)
mtcars2 <- within(mtcars, hp <- scale(hp)[,1])
mod <- lm(mpg ~ hp * disp + factor(am), data = mtcars2)
plot_slopes(mod, variables = "hp", condition = "disp") +
geom_hline(yintercept = 0, linetype = 2) +
geom_rug(aes(x = disp), data = mtcars2, sides = "b") +
theme_grey()
Created on 2023-09-28 with reprex v2.0.2