Interpret Interaction Terms

interaction.r

#generate data
library(MASS)
covvar <- matrix(c(1,.6,.6,.6,1,.6,.6,.6,1),nr = 3)
means <- c(20,25,30)

data <- as.data.frame(mvrnorm(100,means,covvar))

#dichotomize the dependent variable
data$V1 <- ifelse(data$V1 > 20,1,0)

#define the interesting values for variable V3
iv <- seq(27,33,1)

#create the vector to store the betas and standard errors
betas <- rep(0, length(iv))
error <- rep(0, length(iv))

#build the loop to carry out the estimation
for (i in 1:length(iv)){
   data.m <- data
   data.m$V3 <- data$V3 - iv[i]
   model <- glm(V1 ~ V2 + V3 + V2:V3, data = data.m, family = "binomial")
   betas[i] <- coef(model)[2]
   error[i] <- coef(summary(model))[2,2]
}

#put the vectors together in a data frame for ggplot2
dat <- data.frame(iv = iv, middle = betas, upper = betas + error*1.96, lower = betas - error*1.96)

#plot the changes in significance with ggplot2
plot <- ggplot(dat, aes(x = factor(iv)))
plot <- plot + geom_point(aes(y = middle), size = 5)
plot <- plot + geom_point(aes(y = lower, shape = 2), size = 4)
plot <- plot + geom_point(aes(y = upper, shape = 6), size = 4)
plot <- plot + geom_segment(aes(x = factor(iv), y = middle, xend = factor(iv), yend = upper), linetype = 2)
plot <- plot + geom_segment(aes(x = factor(iv), y = middle, xend = factor(iv), yend = lower), linetype = 2)
plot <- plot + geom_hline(yintercept = 0)
plot <- plot + xlab("Interesting Values of V3") + ylab("Logit coef of V2. w/95%ci")
plot <- plot + theme_bw()
plot