Let's say I have a 2x2 design where participants are either in condition A or condition B and, within each condition, either get exposed to exposure C or exposure D.
First, I want to test whether exposure D > exposure C, within each condition. I can do one of two things.
First, I can subset the data to the appropriate condition and check the coefficient on a dummy variable set to 1 for exposure_d.
lm(outcome ~ exposure_d, data = subset(data, condition == "A") lm(outcome ~ exposure_d, data = subset(data, condition == "B")
Or I can specify a full model and run a contrast
library(emmeans) model <- lm(outcome ~ exposure_d * condition_a, data = data) means = emmeans(model, ~ condition_a | exposure_d) contrast_table = contrast(means, method = "pairwise")
Are these methods formally equivalent? I tested the two on simulated data and they seem to produce identical coefficients and standard errors, but I worry things might change when I start adding in more control variables (especially since the
contrast() function throws the warning that the cell means and standard errors may not be accurate in the presence of interactions). Which of these methods is preferred?