## EXAMPLE 1:
## Data from Rothman and Keller (1972) evaluating the effect of joint exposure
## to smoking and alcohol use on the risk of cancer of the mouth and pharynx
## (cited in Hosmer and Lemeshow, 1992):
can <- c(rep(1, times = 231), rep(0, times = 178), rep(1, times = 11),
rep(0, times = 38))
smk <- c(rep(1, times = 225), rep(0, times = 6), rep(1, times = 166),
rep(0, times = 12), rep(1, times = 8), rep(0, times = 3), rep(1, times = 18),
rep(0, times = 20))
alc <- c(rep(1, times = 409), rep(0, times = 49))
dat.df01 <- data.frame(alc, smk, can)
## Table 2 of Hosmer and Lemeshow (1992):
dat.glm01 <- glm(can ~ alc + smk + alc:smk, family = binomial, data = dat.df01)
summary(dat.glm01)
## What is the measure of effect modification on the additive scale?
epi.interaction(model = dat.glm01, param = "product", coef = c(2,3,4),
conf.level = 0.95)$reri
## Measure of interaction on the additive scale: RERI 3.73
## (95% CI -1.84 to 9.32), page 453 of Hosmer and Lemeshow (1992).
## What is the measure of effect modification on the multiplicative scale?
## See VanderWeele and Knol (2014) page 36 and Knol and Vanderweele (2012)
## for details.
epi.interaction(model = dat.glm01, param = "product", coef = c(2,3,4),
conf.level = 0.95)$multiplicative
## Measure of interaction on the multiplicative scale: 0.091 (95% CI 0.14 to
## 5.3).
## EXAMPLE 2:
## Rothman defines an alternative coding scheme to be employed for
## parameterising an interaction term. Using this approach, instead of using
## two risk factors and one product term to represent the interaction (as
## above) the risk factors are combined into one variable comprised of
## (in this case) four levels. Dummy variables are added to the data set using
## the following code:
## a.neg b.neg: 0 0 0
## a.pos b.neg: 1 0 0
## a.neg b.pos: 0 1 0
## a.pos b.pos: 0 0 1
dat.df01$d <- rep(NA, times = nrow(dat.df01))
dat.df01$d[dat.df01$alc == 0 & dat.df01$smk == 0] <- 0
dat.df01$d[dat.df01$alc == 1 & dat.df01$smk == 0] <- 1
dat.df01$d[dat.df01$alc == 0 & dat.df01$smk == 1] <- 2
dat.df01$d[dat.df01$alc == 1 & dat.df01$smk == 1] <- 3
dat.df01$d <- factor(dat.df01$d)
## Table 3 of Hosmer and Lemeshow (1992):
dat.glm02 <- glm(can ~ d, family = binomial, data = dat.df01)
summary(dat.glm02)
## What is the measure of effect modification on the additive scale?
epi.interaction(model = dat.glm02, param = "dummy", coef = c(2,3,4),
conf.level = 0.95)
## Measure of interaction on the additive scale: RERI 3.74
## (95% CI -1.84 to 9.32), page 455 of Hosmer and Lemeshow (1992).
## EXAMPLE 3:
## Here we demonstrate the use of epi.interaction when you're working with
## multilevel data. Imagine each of the study subjects listed in data frame
## dat.df01 are aggregated into clusters (e.g., community health centres).
## Assuming there are five clusters, assign each subject to a cluster:
if (FALSE) {
set.seed(1234)
dat.df01$inst <- round(runif(n = nrow(dat.df01), min = 1, max = 5), digits = 0)
table(dat.df01$inst)
## Fit a generalised linear mixed-effects model using function glmer in the
## lme4 package, with variable inst as a random intercept term:
dat.glmer01 <- glmer(can ~ alc + smk + alc:smk + (1 | inst), family = binomial,
data = dat.df01)
summary(dat.glmer01)
## What is the measure of effect modification on the additive scale?
epi.interaction(model = dat.glmer01, param = "product", coef = c(2,3,4),
conf.level = 0.95)
## Measure of interaction on the additive scale: RERI 3.74
## (95% CI -1.84 to 9.32), identical to that produced above largely because
## there's no strong institution-level effect due to the contrived way we've
## created the multilevel data.
}
Run the code above in your browser using DataLab