Suppose I have a dataset with outcome y, and covariates x and z, for example:
N = 100
x = rnorm(N)
z = rnorm(N)
y = 3*x + 2*z + rnorm(N)
dataset = data.frame(x=x, z=z, y=y)
For a univariate regression of y on x, I can obtain a plot with confidence interval, as follows:
ggplot(dataset) +
geom_point(aes(x=x, y=y)) +
stat_smooth(method='lm', formula = y~x)
QUESTION: How could I get the same plot for a multivariate regression of y on x and z, where the line corresponds to a specific value of z (say, z=0.42)?
I can draw the line as follows:
model <- lm(data=dataset, formula = y~x+z)
special_z = 0.42
ggplot() +
geom_point(data=dataset, aes(x=x, y=y)) +
geom_abline(
slope = coef(model)["x"],
intercept = coef(model)["(Intercept)"] + special_z*coef(model)["z"],
color = "blue")
However, how could I add the accompanying confidence interval to this line?
emmeans
orggeffects
packages to compute the predicted values and CIs you need might be the easiest way to get there