アウトカムが割り付けに影響する場合の因果推論

10880 ワード
R 差分の差 傾向スコア 因果推論 DID R テキストリンク

傾向スコアによるIPW推定量などの前提は

p(z,y|x)=p(z|x) p(y|x)

と共変量を条件づけたときに割り付けとアウトカムが独立になる事だが、現実的に完璧に満たされる事も少ない気がするのでどうなるか実験。
ちなみに傾向スコアだとピンク本には加法処置モデルというので対処できるが、イマイチ理解できなく情報が殆ど無かったのでそれはパス。


library(ggplot2)
library(tidyverse)
library(cobalt)

N=5000

x=rnorm(N,3,2)
rel_drop=rnorm(N,0,2) #dropにかかわる変数

y0=rnorm(N,0,1)

#割り付けの仮定
#q=2*x+1
#割り付けの仮定,開始前の状態y0がy1に影響
#傾向スコアを使う場合、割り付けが欠測値に依存するのでランダムでない欠測
#x,アウトカムの初期値y0が高いほど割り付ける確率も高くなるとする
q=2*y0+x-1

p=plogis(q)
z=rbernoulli(N,p=p)
hist(p)

#アウトカムの仮定,開始前の状態y0と共変量xが高いほどy1も高くする。
#介入の差は2とする
y1=2*z+1.2*x+y0+rnorm(N)

#単純平均
df %>% group_by(z) %>% summarise(mean=mean(dy))

#  z  mean
# <dbl> <dbl>
#   1     0  1.92
# 2     1  6.26

#差は 4.3
#2よりだいぶ過剰。これは割り付けられた人は元からやる気的な共変量が高いので、そのため。

df=tibble(x,z=as.numeric(z),y0,y1,rel_drop,dy=y1-y0,id=1:N,ps)

#ランダム割り付け
z=rbernoulli(N,p=0.8)
y1=2*z+1.2*x+y0+rnorm(N)
# z  mean
# <dbl> <dbl>
#   1     0  3.68
# 2     1  5.63
#効果は2.05と推定

#重回帰
res.lm=lm(dy~x+z+y0,df)
summary(res.lm)
#2.05と推定

library(DRDID)

out <- drdid(yname = "y", tname = "time", idname = "id", dname = "z",
             xformla= ~rel_drop+x,
             data = tdf, panel = TRUE)

summary(out)
#ATE is 2.01と推定


#IPW
library(ltmle)

#df2=tibble(x=df$x,y0=df$y0,drop=df$rel_drop,A=as.numeric(df$z),Y=df$dy)
df2=tibble(x=df$x,A=as.numeric(df$z),Y=df$dy)
#df2=data.frame(Y=y,A=as.numeric(z),W=df$x)

#df2=data
result1 <- ltmle(df2, Anodes="A", Ynodes="Y",abar=list(1,0))

summary(result1)
#ATEは2.01と推定、意外とできている

金メダル:DL,DID
銀メダル:重回帰,RCT
予選落ち:単純平均

という結果に。
IPWはランダムでない欠測があり条件を満たしきっておらず、DIDは並行トレンドの仮定が満たさているという優位があるにもかかわらず同程度。岩波DSでRCTが殴られていたが、この実験は不服従がないので、金メダルでないのは私のやり方が間違っている可能性も。