我是个半路出家的半吊子统计学爱好者,没学过线性代数,关于数学的基础知识薄弱,关注论坛好久了,看到许多帖子讨论的内容还是非常有深度,论坛也是藏龙卧虎,学到了非常多的统计学思想和知识。
最近在看关于倾向性得分相关的文献,由此引申了对混杂因素以及对混杂控制的一些思考,网络上的教程多是关于方法如何使用,少有提及在什么场景下使用,解决什么问题的。(我觉得这也是很多统计方法被误用、乱用的原因)
关于应用场景有几个问题实在是没有想明白,还向各位大佬请教。
正文:
通常对混杂因素U的定义为:U对核心解释变量A有影响,U对结局B有影响,且U不为A→B的中间变量,那么我们称U为A对B影响的混杂因素。
一个经典的例子,如:探讨高血压对死亡的影响,糖尿病能导致高血压的发生,同时糖尿病也能导致结局的发生,如果不考虑糖尿病则会过高估计高血压对死亡的影响。
如何去控制混杂因素,无论是教材还是互联网,给出的回答是:①分层分析;②多元回归,将混杂因素纳入进回归模型;③倾向性得分法。分层分析在大型数据/混杂因素较多的情况下不便于实施,因此在此不作讨论,接下来主要讨论回归与PS法。
首先我模拟了一个这样的框架,A是我们研究的某种干预(二分类)对结局B的影响,其中U1,U2是两个混杂因素且U1对U2有影响。
rm(list = ls()); gc()
set.seed(123)
n <- 100000
U1 <- rnorm(n, mean = 0, sd = 1)
U2 <- 0.7 * U1 + rnorm(n, mean = 0, sd = 1)
A <- 0.5 * U1 + 0.6 * U2 + rnorm(n, mean = 0, sd = 1)
A <- ifelse(A >0, 1, 0)
B <- 0.7 * A + 0.3 * U1 + 0.8 * U2 + rnorm(n, mean = 0, sd = 1)
data <- data.frame(U1, U2, A, B)
首先,我们单独分析了干预A对结局B的影响
f0 <- lm(B ~ as.factor(A))
summary(f0)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.677771 0.006181 -109.6 <2e-16 ***
as.factor(A)1 2.063214 0.008731 236.3 <2e-16 ***
回归系数远大于我们假定的0.7,因为有其他的混杂因素的效应被错误的归因到A的效应上。
第一种,最常见的控制混杂因素的方法-多元回归,将U1,U2作为控制变量纳入回归。
lm <- lm(data = data, B ~ as.factor(A) + U1 + U2)
summary(lm)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0004866 0.0050435 0.096 0.923
as.factor(A)1 0.7002599 0.0078488 89.218 <2e-16 ***
U1 0.3081652 0.0039986 77.068 <2e-16 ***
U2 0.7903921 0.0033975 232.639 <2e-16 ***
混杂因素的效应被控制住。
第二种,使用倾向性得分加权
weight <- glm(data = data, as.factor(A) ~ U1 + U2, family = binomial(link = "logit"))
ps <- weight[["fitted.values"]]
ps <- unname(ps)
data <- cbind(data, ps)
library(dplyr)
data <- data %>%
mutate( weight = ifelse(A == 1, 1/ps, 1/(1-ps)))
library(survey)
design <- svydesign(ids = ~1,
weights = ~weight,
data = data)
f1 <- svyglm(B ~ as.factor(A),
design = design,
family = gaussian(link = "identity"))
summary(f1)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.06745 0.01184 -5.696 1.23e-08 ***
as.factor(A)1 0.84162 0.01731 48.607 < 2e-16 ***
使用双重稳健的倾向性得分加权
f2 <- svyglm(B ~ as.factor(A) + U1 + U2,
design = design,
family = gaussian(link = "identity"))
summary(f2)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.004957 0.006603 0.751 0.453
as.factor(A)1 0.695921 0.009639 72.202 <2e-16 ***
U1 0.310406 0.005776 53.737 <2e-16 ***
U2 0.792583 0.005413 146.427 <2e-16 ***
那么:
①在第一个回归中,U1,U2,A之间存在复杂的相关关系,破坏了回归分析自变量相互独立的假设,为什么不会对估计结果产生影响?为什么不需要纳入U1 U2 的交互项?
②既然将混杂因素纳入回归就可以进行较好的控制,为什么还要作倾向性得分匹配/加权?他们的应用场景到底是什么样的?
还望各位老师不吝赐教!