IPW+Cox回帰

ATEとしてハザード比を求める必要があったので、IPWを使った。 ダミーデータで備忘録を残しておく。

library(survminer)
library(survival)
library(cobalt)
library(tidyverse)

# データの準備
diabetic_mod <- diabetic %>% 
  filter(eye == "left") %>% 
  mutate(flag_highrisk = if_else(risk >= 10, 1, 0))

# KM曲線の描画
surv <- survfit(Surv(time, status) ~ trt, data = diabetic_mod)
ggsurvplot(surv)

# Cox回帰
fit <- coxph(Surv(time, status) ~ trt, data = diabetic_mod)
summary(fit)

# IPW + Cox回帰
# weightsにIPWの重みを与える
library(WeightIt)
model_ipw <- weightit(trt ~ flag_highrisk,
                  data = diabetic_mod,
                  estimand = "ATE", 
                  method = "ps"
)

fit_ipw <- coxph(Surv(time, status) ~ trt, 
                 data = diabetic_mod, 
                 weights = get.w(model_ipw))
summary(fit_ipw)

まあ、専門外の人に説明するときは傾向スコアマッチングの方が説明しやすいんだけど、それだとATTとなってしまうので...

重み付けCox回帰だとcoxphwパッケージを用いる方が適切かもしれないが、今回はナイーブにcoxph()で実施している。 結果はcoxphw()のtemplate="PH"とした場合とほぼ一致する。 stats.stackexchange.com