Drawing a regression line with interaction in ggplot2

To draw a regression line with interaction in ggplot2, you can use the ggpredict() function from the ggeffects package.

library(tidyverse)
library(ggeffects)

First, you need to generate an example dataset for reproducibility:

set.seed(42)
example_data = 
  data.frame(treatment   = as.factor(c(rep("t1", 50), rep("t2", 50))),
             time        = rep(1:50, 2), 
             error       = rnorm(100, 0, 1))                      %>% 
  mutate(    slope       = ifelse(treatment == "t1", -0.2, -0.3)) %>% 
  mutate(    measurement = 2.5 + time * slope + error)          

Next, you need to fit a linear model to your data:

model = lm(measurement ~ time:treatment, data = example_data)
summary(model)

Once you have fitted your model, you can generate predictions using ggpredict():

g <- ggpredict(model, terms=c("time", "treatment")) 

plot(g)

The plot() function will produce a plot with superposed lines, showing the predicted values for each treatment group.

You can also add a facet to the plot to show the results for each treatment group separately:

plot(g) + facet_wrap(~treatment)

Here is a base R approach if you want to see what is going on under the hood with the prediction lines:

#### Get Prediction Data ####
nd1 <- data.frame(
  time = seq(
    min(example_data$time),
    max(example_data$time),
    length.out=200
  ),
  treatment = "t1"
) # for treatment group 1

nd2 <- data.frame(
  time = seq(
    min(example_data$time),
    max(example_data$time),
    length.out=200
  ),
  treatment = "t2"
) # for treatment group 2

After creating both dataframes, one for the first group and the other for the second group, you can create predictions based off those values:

pd1 <- predict(model,newdata=nd1)
pd2 <- predict(model,newdata=nd2)

Then you just have to plot the points based on group:

#### Plot Points ####
plot(example_data$time,
     example_data$measurement,
     xlab="Time",
     ylab="Measurement",
     main="Interaction Plot",
     pch=23,
     bg=ifelse(example_data$treatment == "t1", "steelblue",'gray'))

And finally overlay the lines from the predictions:

#### Draw Lines ####
lines(nd1$time,pd1,col="steelblue",lwd=2)
lines(nd2$time,pd2,col="gray",lwd=2)