logo
生物科学研究所 井口研究室
Laboratory of Biology, Okaya, Nagano, Japan
Home

Cox 比例ハザードモデルを多重代入法で計算: 欠損値処理

井口豊（生物科学研究所，長野県岡谷市）
最終更新： 2023 年 3 月 11 日

1. はじめに

Cox 比例ハザードモデル（Cox 回帰）のデータに欠損値がある場合も少なくない。そのような場合の計算手手順として， R の mice パッケージを利用し，多重代入法（multiple imputation）で欠損補完して計算してみる。

2. 多重代入法から Cox 回帰へ

データには， survival パッケージ付属の mgus2 から，一部の変数を抜き出して利用した。

以下が， R スクリプトである。 


#############
library(mice)
library(survival)

# survival パッケージ付属データ mgus2
dat<- mgus2[,
   c("age", "sex", "hgb", "futime", "death")
]

# 欠損確認
dat[!complete.cases(dat),]

# 欠損数確認
sum(is.na(dat))

# 補完
dat.imp<- mice(
  dat,
  m = 100,
  method = "pmm",
  seed = 10,
  printFlag = FALSE
)

# 欠損数確認
sum(is.na(complete(dat.imp)))

# 補完データで Cox 回帰
cox.imp<- with(
   dat.imp,
   coxph(Surv(futime, death) ~ sex + age + hgb)
)

# 結果
options(digits=3)
summary(pool(cox.imp))

# ハザード比，信頼区間も算出
summary(
  pool(cox.imp),
  conf.int = TRUE,
  exponentiate = TRUE
)

################

結果は以下の通りである。最後の summary 部分を抜粋した。 


###########################
> # 結果
> options(digits=3)
> summary(pool(cox.imp))
  term estimate std.error statistic  df  p.value
1 sexM   0.4889   0.06735      7.26 957 8.04e-13
2  age   0.0554   0.00339     16.33 957 0.00e+00
3  hgb  -0.1483   0.01733     -8.56 950 0.00e+00
> 
> # ハザード比，信頼区間も算出
> summary(
+   pool(cox.imp),
+   conf.int = TRUE,
+   exponentiate = TRUE
+ )
  term estimate std.error statistic  df  p.value 2.5 % 97.5 %
1 sexM    1.631   0.06735      7.26 957 8.04e-13 1.429  1.861
2  age    1.057   0.00339     16.33 957 0.00e+00 1.050  1.064
3  hgb    0.862   0.01733     -8.56 950 0.00e+00 0.833  0.892
##############################################################

デフォルトで，ハザード比と信頼区間も算出してくれたほうが便利だと思うが，なぜかそうなっていない。

関連ページ

Home