IPTW カプランマイヤー曲線に Number at Risk を入れる方法

IPTW カプランマイヤー曲線において、任意の時点における Number at Risk を書き入れる方法

R で作成する方法

＞＞もう統計で悩むのを終わりにしませんか？

↑1万人以上の医療従事者が購読中

IPTW カプランマイヤー曲線
IPTW カプランマイヤー曲線の書き方
まとめ

IPTW カプランマイヤー曲線

IPTW は、逆確率重み付けとも呼ばれ、ランダム割り付けしていない群間の交絡因子を調整する方法の一つ

詳しくはこちらも参照

toukeier.hatenablog.com

IPTＷは、EZR を使うと簡単に作成できるが、後述のようなスクリプトでも作成可能

IPTW で調整したカプランマイヤー推定量を計算して、カプランマイヤー曲線を書く方法と、図の下に Number at Risk（何人イベントを起こさずに残っているか）を示す方法のスクリプトを紹介する

IPTW カプランマイヤー曲線の書き方

R の survival パッケージの、survfit を使って推定量を求めるときに、weight で IPTW を指定するという方法を取る

曲線自体を書くのは、survminer パッケージの ggsurvplot を使う

そうすると、Number at Risk を書き入れるのが簡単である

まず、必要なパッケージを読み込む

# 必要なパッケージを読み込み
library(survival)
library(survminer)

＞＞もう統計で悩むのを終わりにしませんか？

↑1万人以上の医療従事者が購読中

標準のカプランマイヤー曲線

ここでは、サンプルデータとして、ovarian データセットを用いる

survfit でカプランマイヤー推定量を求める

まずは、IPTW 調整なしで書いてみる

km_fit0 <- survfit(Surv(futime, fustat) ~ rx, data=ovarian)
plot(km_fit0, col=1:2)

IPTW で調整していない、rx 因子別のカプランマイヤー曲線が書かれる

傾向スコア及び IPTW の作成

# rx = 1,2 を rx1 = 0,1 に変換
rx1 <- ovarian$rx-1

# 傾向スコアの作成（EZR で作成しても可）
propensity <- glm(rx1 ~ age + resid.ds + ecog.ps, 
                  family = binomial(link='logit'),
                  data = ovarian)$fitted.values

# IPTW の作成（EZR で作成しても可）
weights <- ifelse(rx1==1, 1/propensity, 1/(1-propensity))

IPTW カプランマイヤー曲線を書く

# カプランマイヤー曲線の作成
km_fit <- survfit(Surv(futime, fustat) ~ rx, data = ovarian, 
                  weights = weights)

plot(km_fit, col=1:2)

IPTW カプランマイヤー曲線はこのようになる

未調整のカプランマイヤー曲線と比較して、ほんの少し黒と赤が近づいたように見える

Number at Risk 入りのグラフを書く

# グラフを作成
ggsurv <- ggsurvplot(km_fit, data = ovarian, 
                     risk.table = TRUE, 
                     risk.table.col = "strata",
                   ggtheme = theme_minimal())

# Number at risk を表示
ggsurv$table$variables$nrisk <- 
  paste("At risk:", ggsurv$table$variables$n.risk)

# グラフを表示
ggsurv