# Holding covariates constant to plot MLR model on 2d scatterplot in R

I am currently working with a MLR model comprising 1 numeric/continuous predictor variable (x1), several nominal categorical variables (x2 ... xi), and an interaction term between the continuous variable and one of the categorical variables (x1*x2). I wish to plot the relationship between y and x1 on a 2d scatterplot (ideally in ggplot2) with a line of best fit and confidence intervals. I believe that this is possible in theory, so long as the categorical variables are fixed at a pre-selected/reference level for plotting purposes. My understanding is that the model terms associated with the different levels of the categorical variables will just shift the intercept up/down the y-axis, and not fundamentally change the nature of the relationship between y and x1 (i.e., the slope). However, I have not been able to work out how to do generate such a plot in practice.

A reprex of a similar (toy) model is provided below:

# Generate data.frame
df<-
data.frame(
"y"=c(32,27,29,41,26,23,35,36,35,32,29,30,40,27,38,21,31,26,26,34,41,29,26,24),
"x1"=c(28,32,36,40,44,48,52,56,60,64,68,72,72,68,64,60,56,52,48,44,40,36,32,28),
"x2"=c("M","F"),
"x3"=c("A","B","C"),
"x4"=c("I","II","III","IV")
)

df$$x2<- as.factor(df$$x2)

df$$x3<- as.factor(df$$x3)

df$$x4<- as.factor(df$$x4)

# Generate MLR model
lm<-
lm(
y~
x1*x2+
x3+
x4,
data=df
)

summary(lm)


Imagine that I wish to plot the relationship between y and x1 on a 2d scatterplot, and use x2=M, x3=A, and x4=III as the reference levels at which I wish to fix these covariates. How would I do so? I have tried manually calculating the predicted values for each of the data points if they were associated with these reference levels, and plotting them all, like so:

# Manually generate predictions
df$$fixed<- 32.424825+ # intercept df$$x1*-0.003497+ # term for x1
1*-4.525641+ # term for x2=M
1*0+ # term for x3=A
1*-3.666667+ # term for x4=III
1*(df\$x1*0.153846) # x1*x2 interaction term when x2=M

# Plot df$$fixed vs df$$x1
library(ggplot2)

p<-
ggplot(
data=df,
mapping=aes(
x=x1,
y=fixed
)
)+
geom_point()+
geom_smooth(
method="lm",
se=T
)

p


This approach has not worked for me, specifically as i) I wish to show the confidence intervals around the line of best fit, and ii) I have a very large dataset (~500k observations). In essence, I think that what I am doing here is passing ~500k fitted values to ggplot2, and then asking it to plot the line of best fit and associated confidence intervals. Unsurprisingly, there is essentially no uncertainty--given that the values are fitted, and that there are so many of them--so the confidence intervals are arbitrarily small. This would thus not appear to be the correct approach.

Is anyone aware of a method/package/function where I can plot y~x1 (ideally in ggplot2 graphics) with x2 .... xi held constant, in a way that still shows the uncertainty in the data (i.e., with confidence intervals)?

Thank you very much.

• You might be able to do this with the ggeffects package. See its vignettes. Commented Apr 7, 2023 at 16:12
• @dipetkov thank you very much. I was able to use ggeffects to make these plots. Specifically, for anyone else who stumbles across this in future, I used ggpredict() and was able to plot using the ggpredict object directly through ggplot(), which proved much more powerful/flexible than the ggeffects::plot() function (which is nevertheless still helpful). What I actually wanted to plot was Marginal Effects--specifically, Marginal Effects at Representative Values. I did not know the terminology, hence why I struggled to find the correct tool for the job. Commented Apr 14, 2023 at 8:19
• Thanks for the update. Aside: It's okay to answer one's own question and accept the answer; see here. Your question seemed to emphasize the plotting aspects; it wasn't clear that you have questions about the theory too. On that topic, ggeffects is a visualization interface to other packages that do the hard lifting of computing marginal effects such as emmeans. Admittedly, emmeans has a bit of a learning curve but its documentation and vignettes are very helpful. See here. Commented Apr 14, 2023 at 8:31
• @dipetkov thank you again. I think that I understand the theory behind the plots that I am generating, but wasn't clear on the terminology. Other people might be in a similar situation in future if I was. Hopefully this post will help them. I have answered the question myself as suggested. Commented Apr 14, 2023 at 8:43