Does an additional hour of studying have the same association with test scores for every study method? An interaction lets that slope differ between groups; here we fit and interpret one using only base R.
These are synthetic scores for 120 students. We deliberately build in a steeper study-hours slope for the workshop group, then add random noise.
Each row of students represents one student. runif() draws study hours between one and nine, while rep() assigns 60 rows to each method. factor() tells R that the method is a category, not a numeric measurement. set.seed() makes the random draws reproducible when you rerun the example.
set.seed(2026)
students <- data.frame(
hours = runif(120, min = 1, max = 9),
method = factor(
rep(c("Self-study", "Workshop"), each = 60),
levels = c("Self-study", "Workshop")
)
)
students$hours_centered <- students$hours - 5
workshop <- as.numeric(students$method == "Workshop")
students$score <- 65 +
3 * students$hours_centered +
4 * workshop +
2 * students$hours_centered * workshop +
rnorm(120, mean = 0, sd = 3)With R’s default treatment coding, the first factor level, Self-study, is the reference group. Subtracting five makes hours_centered = 0 mean five hours of studying, not zero hours.
For example, two hours becomes -3 and eight hours becomes 3. Centering changes where we interpret the intercept and group difference; it does not change the fitted predictions when both predictors and their interaction remain in the model.
The workshop variable is 0 for self-study and 1 for workshop students. That makes the score-generation rule easier to follow:
rnorm() adds normally distributed noise with mean 0 and standard deviation 3, so students do not all fall exactly on these lines. The fitted coefficients will therefore differ somewhat from the values we used to generate the data.In an R model formula, hours_centered * method includes both predictors and their interaction: hours_centered + method + hours_centered:method. There is no need to create a product column yourself.
The ~ separates the outcome, score, from the predictors. Using only hours_centered + method would allow different group intercepts but require the same hourly slope in both groups. Using * allows the slopes to differ too.
interaction_model <- lm(score ~ hours_centered * method, data = students)
coefficients <- coef(interaction_model)
round(coefficients, 2) (Intercept) hours_centered
64.54 2.69
methodWorkshop hours_centered:methodWorkshop
4.95 2.08 lm() fits the model, coef() extracts its named coefficient estimates, and round() changes only how we display them. Read these estimates as follows:
| Coefficient | Interpretation |
|---|---|
(Intercept) |
Predicted self-study score at five hours: 64.54 points. |
hours_centered |
Self-study slope: 2.69 additional predicted points per extra hour. |
methodWorkshop |
Workshop minus self-study difference at five hours: 4.95 points. |
hours_centered:methodWorkshop |
Workshop minus self-study slope difference: 2.08 points per hour. |
The interaction coefficient is not the workshop slope. Add it to the self-study slope:
The [["..."]] expressions below retrieve individual coefficients by name. The interaction estimate is the amount we add to the reference group’s slope, not a replacement for that slope.
self_study_slope <- coefficients[["hours_centered"]]
workshop_slope <- self_study_slope +
coefficients[["hours_centered:methodWorkshop"]]
round(c("Self-study" = self_study_slope, "Workshop" = workshop_slope), 2)Self-study Workshop
2.69 4.77 An extra hour is associated with about 2.69 predicted points for self-study and 4.77 for workshop students. Neither coefficient is an overall slope averaged across both groups.
Use predict() to compare groups at specific hours. New data must use the same centering value as the training data.
Here, expand.grid() creates all six combinations of three study durations and two methods. These rows are scenarios, not additional observed students. predict() uses the fitted model to estimate a mean score for each scenario; it does not generate new random scores.
predictions <- expand.grid(
hours = c(2, 5, 8),
method = levels(students$method)
)
predictions$hours_centered <- predictions$hours - 5
predictions$predicted_score <- predict(interaction_model, newdata = predictions)
data.frame(
hours = predictions$hours,
method = predictions$method,
predicted_score = round(predictions$predicted_score, 1)
) hours method predicted_score
1 2 Self-study 56.5
2 5 Self-study 64.5
3 8 Self-study 72.6
4 2 Workshop 55.2
5 5 Workshop 69.5
6 8 Workshop 83.8Compare rows with the same hours: the workshop minus self-study difference changes as hours increase. The methodWorkshop coefficient describes that difference only at five hours, not at every study duration.
For example, at eight hours the predicted scores are 72.6 for self-study and 83.8 for the workshop. Their difference is about 11.2 points: the five-hour group difference plus three hours of the additional workshop slope.
To see that changing gap, subset() separates the prediction rows by method. plot() draws the self-study line, and lines() adds the workshop line to the same axes. Both lines come from our single interaction model, not separate regressions.
self_study_predictions <- subset(predictions, method == "Self-study")
workshop_predictions <- subset(predictions, method == "Workshop")
plot(predicted_score ~ hours, data = self_study_predictions,
type = "l", col = "steelblue4", lwd = 2,
ylim = range(predictions$predicted_score),
xlab = "Study hours", ylab = "Predicted score")
lines(predicted_score ~ hours, data = workshop_predictions,
col = "darkorange3", lwd = 2, lty = 2)
legend("topleft", legend = c("Self-study", "Workshop"),
col = c("steelblue4", "darkorange3"), lwd = 2, lty = c(1, 2), bty = "n")
The different slopes show the fitted interaction. These are point estimates, not guaranteed individual scores or evidence that either method causes improvement. For uncertainty, confint(interaction_model) gives coefficient confidence intervals; the interaction interval concerns the difference in slopes.
For syntax, see the R documentation for lm(); for interpretation, see UCLA’s interaction guide.
The companion app makes the page’s central comparison adjustable: include the interaction to allow different fitted slopes, or remove it to require parallel slopes. You can also change the centering value, compare the methods at a selected number of study hours, and reveal the synthetic observations behind the lines.
This shows the app’s default state. The controls are intentionally disabled because the workshop page is not running a Shiny server.
Self-study slope 2.69 points per hour
Workshop slope 4.77 points per hour
At eight hours, the fitted scores are 72.6 for self-study and 83.8 for the workshop.
Download the complete Shiny app
After extracting the ZIP, install shiny if needed, open R in the extracted example-linear-regression-interactions folder, and run:
shiny::runApp()The app runs locally and does not use an API or upload data. Like the example above, it uses synthetic observations and illustrates fitted associations rather than causal effects. The app intentionally omits confidence intervals so the interaction and centering interpretations remain the focus.