中介分析(mediation analysis)是分析因果关系的一种特殊方法。
image.png
- 模拟数据
set.seed(888)
treat_flg <- rbinom(1000,1,0.3) # binary treatment variable with binomial distribution
crp <- round(abs(40*treat_flag+rnorm(1000,100,30)),1) # regress crp on treatment
lp <- 10*treat_flag+0.02*crp+(rnorm(1000,-5,2)) # linear prediction
link_lp <- exp(lp)/(1+exp(lp)) # link function
mort <- (runif(1000) < link_lp) #logic value
df <- data.frame(crp=crp, treat_flg = treat_flg,mort = mort)
head(df)
image.png
- 中介分析
#install.packages("mediation")
library(mediation)
model.m <- lm(crp~treat_flg,data=df[1:100,])
model.y <- glm(mort~treat_flg+crp,family = binomial, data=df[1:100,])
med.out <- mediate(model.m,model.y,treat="treat_flg",
mediator="crp",robustSE=TRUE,sims =100)
summary(med.out)
image.png
plot(med.out)
image.png
敏感性分析(sensitivity analysis)
mediation analysis的前提是没有confounder，其对中间变量和结局变量都有影响。因此应进行敏感性分析。
image.png
# 接上面的结果
probit.y <- glm(mort~treat_flg+crp,family = binomial(probit),data=df[1:100,])
med.out1 <- mediate(model.m,probit.y,treat="treat_flg",mediator = "crp",robustSE=TRUE,sims=100)
sens.out <- medsens(med.out1,rho.by = 0.1,effect.type = "indirect",sims = 100) #会有warning
par(mfrow = c(2,1))
plot(sens.out,sens.par = "rho", main = "Mortality")
image.png
结果解读：
image.png
参考资料
文中代码及截图均来自章仲恒教授的丁香园公开课：中介分析
网友评论