如何在多元 Coxph 中获得特定变量的 wald 检验？

Question

我使用 R 生存包拟合了一个多元 Cox 模型，如下所示：

library(survival)
data(lung)
res.cox1 <- coxph(Surv(time, status) ~  sex + ph.karno + wt.loss, data =  lung)
res.cox1
Call:
coxph(formula = Surv(time, status) ~ sex + ph.karno + wt.loss, 
    data = lung)

              coef exp(coef)  se(coef)      z       p
sex      -0.521839  0.593428  0.174454 -2.991 0.00278
ph.karno -0.015243  0.984873  0.005988 -2.546 0.01091
wt.loss  -0.002523  0.997480  0.006233 -0.405 0.68558

Likelihood ratio test=16.42  on 3 df, p=0.0009298
n= 214, number of events= 152 
   (14 observations deleted due to missingness)

如何获得多元 Cox 模型 (sex + ph.karno + wt.loss) 中每个变量（性别、ph.karno 和 wt.loss）的 Wald 检验的 3 个值？

我试图查看coxph的结构和coxph对象的摘要，发现只有一个wald测试的单一值$wald.test : num 16.5，$ waldtest : Named num [1:3] 1.65e+01 3.00 8.81e-04..- attr(*, "names")= chr [1:3] "test" "df" "pvalue"！

这个测试值对应什么？ 如何获得性别、ph.karno 和 wt.loss 的 Wald 测试的 3 个值？

str(res.cox1)
List of 20
 $ coefficients     : Named num [1:3] -0.52184 -0.01524 -0.00252
  ..- attr(*, "names")= chr [1:3] "sex" "ph.karno" "wt.loss"
 $ var              : num [1:3, 1:3] 3.04e-02 -6.78e-05 2.77e-05 -6.78e-05 3.59e-05 ...
 $ loglik           : num [1:2] -680 -672
 $ score            : num 16.9
 $ iter             : int 4
 $ linear.predictors: num [1:214] 0.0756 0.0756 0.0857 -0.039 0.7232 ...
 $ residuals        : Named num [1:214] -0.147 -2.93 0.58 -1.613 -5.599 ...
  ..- attr(*, "names")= chr [1:214] "2" "3" "4" "5" ...
 $ means            : Named num [1:3] 1.4 82.06 9.83
  ..- attr(*, "names")= chr [1:3] "sex" "ph.karno" "wt.loss"
 $ method           : chr "efron"
 $ n                : int 214
 $ nevent           : num 152
 $ terms            :Classes 'terms', 'formula'  language Surv(time, status) ~ sex + ph.karno + wt.loss
  .. ..- attr(*, "variables")= language list(Surv(time, status), sex, ph.karno, wt.loss)
  .. ..- attr(*, "factors")= int [1:4, 1:3] 0 1 0 0 0 0 1 0 0 0 ...
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:4] "Surv(time, status)" "sex" "ph.karno" "wt.loss"
  .. .. .. ..$ : chr [1:3] "sex" "ph.karno" "wt.loss"
  .. ..- attr(*, "term.labels")= chr [1:3] "sex" "ph.karno" "wt.loss"
  .. ..- attr(*, "specials")=Dotted pair list of 2
  .. .. ..$ strata: NULL
  .. .. ..$ tt    : NULL
  .. ..- attr(*, "order")= int [1:3] 1 1 1
  .. ..- attr(*, "intercept")= num 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. ..- attr(*, "predvars")= language list(Surv(time, status), sex, ph.karno, wt.loss)
  .. ..- attr(*, "dataClasses")= Named chr [1:4] "nmatrix.2" "numeric" "numeric" "numeric"
  .. .. ..- attr(*, "names")= chr [1:4] "Surv(time, status)" "sex" "ph.karno" "wt.loss"
 $ assign           :List of 3
  ..$ sex     : int 1
  ..$ ph.karno: int 2
  ..$ wt.loss : int 3
 $ wald.test        : num 16.5
 $ concordance      : Named num [1:7] 11071 6046 96 22 0 ...
  ..- attr(*, "names")= chr [1:7] "concordant" "discordant" "tied.x" "tied.y" ...
 $ na.action        : 'omit' Named int [1:14] 1 20 36 44 56 63 108 138 178 183 ...
  ..- attr(*, "names")= chr [1:14] "1" "20" "36" "44" ...
 $ y                : 'Surv' num [1:214, 1:2]  455  1010+  210   883  1022+  310   361   218   166   170  ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:214] "2" "3" "4" "5" ...
  .. ..$ : chr [1:2] "time" "status"
  ..- attr(*, "type")= chr "right"
 $ timefix          : logi TRUE
 $ formula          :Class 'formula'  language Surv(time, status) ~ sex + ph.karno + wt.loss
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
 $ call             : language coxph(formula = Surv(time, status) ~ sex + ph.karno +      wt.loss, data = lung)
 - attr(*, "class")= chr "coxph"


str(summary(res.cox1))
List of 14
 $ call        : language coxph(formula = Surv(time, status) ~ sex + ph.karno +      wt.loss, data = lung)
 $ fail        : NULL
 $ na.action   : 'omit' Named int [1:14] 1 20 36 44 56 63 108 138 178 183 ...
  ..- attr(*, "names")= chr [1:14] "1" "20" "36" "44" ...
 $ n           : int 214
 $ loglik      : num [1:2] -680 -672
 $ nevent      : num 152
 $ coefficients: num [1:3, 1:5] -0.52184 -0.01524 -0.00252 0.59343 0.98487 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "sex" "ph.karno" "wt.loss"
  .. ..$ : chr [1:5] "coef" "exp(coef)" "se(coef)" "z" ...
 $ conf.int    : num [1:3, 1:4] 0.593 0.985 0.997 1.685 1.015 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "sex" "ph.karno" "wt.loss"
  .. ..$ : chr [1:4] "exp(coef)" "exp(-coef)" "lower .95" "upper .95"
 $ logtest     : Named num [1:3] 16.42029 3 0.00093
  ..- attr(*, "names")= chr [1:3] "test" "df" "pvalue"
 $ sctest      : Named num [1:3] 1.69e+01 3.00 7.52e-04
  ..- attr(*, "names")= chr [1:3] "test" "df" "pvalue"
 $ rsq         : Named num [1:2] 0.0739 0.9983
  ..- attr(*, "names")= chr [1:2] "rsq" "maxrsq"
 $ waldtest    : Named num [1:3] 1.65e+01 3.00 8.81e-04
  ..- attr(*, "names")= chr [1:3] "test" "df" "pvalue"
 $ used.robust : logi FALSE
 $ concordance : Named num [1:2] 0.646 0.0274
  ..- attr(*, "names")= chr [1:2] "C" "se(C)"
 - attr(*, "class")= chr "summary.coxph"

谢谢！

Answer 1

http://www.sthda.com/english/wiki/cox-proportional-hazards-model

“z”列与多变量 Cox 模型中每个协变量的 Wald 检验统计量相同。

您也可以这样调用 Cox 模型统计数据：

summary(res.cox1)

Answer 2

“Wald 检验”基于回归过程中的参数值呈正态分布的假设。您检查系数估计值（“coef”）除以估计值的标准误差（“coef(se)”）的比率，并查看该比率的 95% 置信区间是否包含零值。操作上说：取 coef +/- 1.96*se(coef) 并查看区间是否包括零。或者等效地，您可以采用比率：coef/se(coef)，并查看它的绝对值是否大于 1.96。当我说“测试”是一个是/否的结果时，也许我是在迂腐，回答“比率值是否位于临界区间内”的问题，而“测试统计量”，如 z 值, 是一个纯数。

实际上在您构建的摘要中报告了 4 个 Wald 测试。其中三个用于单个系数，一个用于整体模型，即名为“wald”的模型。但是您不想要整体模型 Wald 测试。您想要来自summary() 处理结果的“系数”矩阵的结果（而不是来自coxph() 结果的“系数”值。）当您采用这样的比率时，它被分析为 z 检验，所以您这样做不要平方统计（当然，除非您想使用卡方表，此时 Z^2 将用于评估。）

summ.coef <- summary(res.cox1)$coefficients

( Wald.ratios <- summ.coef[,"coef"]/summ.coef[,"se(coef)"] )
       sex   ph.karno    wt.loss 
-2.9912645 -2.5456273 -0.4048609 
identical(Wald.ratios, summ.coef[, "z"])
#[1] TRUE

如果您想按名称关注单个变量：

 summ.coef["sex", "coef"]/summ.coef["sex", "se(coef)"]