方法一

基于R包epiR: Tools for the Analysis of Epidemiological Data
https://cran.r-project.org/web/packages/epiR/index.html

一、R包加载，以及示例数据构建

rm(list = ls())
#install.packages("epiR")
library(epiR) 

## EXAMPLE 1:
## Data from Rothman and Keller (1972) evaluating the effect of joint exposure
## to alcohol and tabacco on 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<- data.frame(alc, smk, can)

#因子化前后，回归分析结果一致
dat$smk <- factor(dat$smk)
dat$alc <- factor(dat$alc)
dat$can <- factor(dat$can)
summary(dat)

二、相乘交互作用

fit <- glm(can ~ alc + smk + alc:smk, family = binomial, data = dat)
summary(fit)$coefficients
coef <- summary(fit)$coefficients[,1]
se <- summary(fit)$coefficients[,2]
#CI <-exp(confint(fit))#另外一种形式显示可信区间
Results <- cbind(exp(coef),exp(coef-1.96*se),exp(coef+1.96*se))

P <- summary(fit)$coefficients[4,4]#0.9219744
dimnames(Results)[[2]] <- c("OR", "lower","upper")
Results

#                   OR      lower      upper
#(Intercept) 0.1500000 0.04457283  0.5047918
#alc1        3.3333333 0.70058649 15.8597280
#smk1        2.9629630 0.68002749 12.9099922
#alc1:smk1   0.9149096 0.15435605  5.4229143
#饮酒与吸烟的乘法交互效应为0.91<1，交互P值为0.92，不显著

三、相加交互作用

## 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)
## 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.
beta1 <- as.numeric(dat.glm01$coefficients[2])
beta2 <- as.numeric(dat.glm01$coefficients[3])
beta3 <- as.numeric(dat.glm01$coefficients[4])
exp(beta3) / (exp(beta1) * exp(beta2))
## Measure of interaction on the multiplicative scale: 0.093.

四、哑变量实现交互作用

## 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:
## 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).

image.png
image.png

五、结果可视化

加法交互效应，堆叠图

ORalc <- 3.3333333
ORsmk <- 2.9629630
RERI <- 3.739848
bar_d <- matrix(c(1, 1, 1, 1,
                  ORalc-1, 0, ORalc-1, 0,
                  ORsmk-1 ,ORsmk-1, 0, 0,
                  RERI, 0, 0, 0),
                c(4,4), byrow = T,
                dimnames = list(c('U','alc','smk1','alc1 & smk1'),c("OR_A1B1","OR_A1B0","OR_A0B1","OR_A0B0")))

plot <- barplot(bar_d, legend=rownames(bar_d),
                args.legend=c(x=4,y=10)#图示的位置
                )

加法交互效应，堆叠图

https://blog.csdn.net/qq_22253901/article/details/118981902

方法二

基于interactionR
https://github.com/epi-zen/interactionR

rm(list = ls())
#devtools::install_github("epi-zen/interactionR")
library(interactionR)
#Example: The joint effect of alcohol and smoking on oral cancer.
data (OCdata)

## fit the interaction model
model.glm <- glm(oc ~ alc*smk, family = binomial(link = "logit"), 
                 data = OCdata)
table_object = interactionR(model.glm, 
                            exposure_names = c("alc", "smk"), 
                            ci.type = "mover", ci.level = 0.95, 
                            em = F, recode = F)
interactionR_table(table_object)#可以输出成Word格式

输出结果

#Example 2: Effect of sports participation and smoking on herniated lumbar disc.
data(HDiscdata)
m2 = glm(h ~ ns*smk, family = binomial(link = 'logit'), 
         data = HDiscdata)
table_object2 = interactionR_boot(m2, ci.level = 0.95, 
                                  em = F, recode = F, 
                                  seed = 1234, s = 1000)
interactionR_table(table_object2)
hist(table_object2$bootstrap)

输出结果
image.png

背景知识学习
https://max.book118.com/html/2018/0607/171193334.shtm
https://www.docin.com/p-467347939.html
https://zhuanlan.zhihu.com/p/138018634
https://max.book118.com/html/2016/1031/60871118.shtm
https://max.book118.com/html/2017/0216/92080186.shtm
https://www.sohu.com/a/155275085_170798