交互作用中的生存分析参考

Posted 2023-02-19

【中文标题】交互作用中的生存分析参考

【英文标题】：survival analysis reference in interaction effect

【发布时间】：2021-11-24 12:55:04

【问题描述】：

我无法更改生存分析中的参考水平。您可以在下面看到我的数据集的一部分，让我从示例中进行。

我希望参考水平为 H1-M1 但 R 将其作为 H3-M3，我该如何更改它？

Ps.：relevel 不起作用

The red circle is the reference part and that will be 1 V1. Thanks


library(survival)
set.seed(123)
st<-sample(0:1,100, replace = T)
tm<-rnorm(100,20,9)
fm<-sample(c("M","F"),100,replace = T)
age<-rnorm(100,75,10)
level1<-sample(c("M1","M2", "M3"),100,replace = T)
level2<-sample(c("V1","V2", "V3"),100,replace = T)
 
dt<-data.frame(status=st,
               Time= tm,
               Sex=   fm,
               Age=age,
               Level1=level1,
               Level2=level2
               )

dt %>%
  coxph(Surv(Time, status) ~Sex + Age + Level1:Level2 , data = .) %>%
  summary()

红圈部分将改为M1 V1。谢谢。

【问题讨论】：

这与您的previous question 有何不同？

status 似乎没有出现在我之前的问题中。

我已经更新了data.frame部分

【参考方案1】：

如果您在公式中使用interaction 而不是使用“:”，您将获得您期望的显示：

dt %>%
     coxph(Surv(Time, status) ~Sex + Age + interaction(Level1 , Level2) , data = .) %>%
      summary()
#-------------------------------
Call:
coxph(formula = Surv(Time, status) ~ Sex + Age + interaction(Level1, 
    Level2), data = .)

  n= 100, number of events= 43 

                                     coef exp(coef) se(coef)      z Pr(>|z|)  
SexM                             -0.33279   0.71692  0.34166 -0.974   0.3300  
Age                              -0.01150   0.98857  0.02012 -0.571   0.5677  
interaction(Level1, Level2)M2.V1  0.16097   1.17465  0.68718  0.234   0.8148  
interaction(Level1, Level2)M3.V1 -0.18250   0.83318  0.68056 -0.268   0.7886  
interaction(Level1, Level2)M1.V2  1.34832   3.85096  0.78176  1.725   0.0846 .
interaction(Level1, Level2)M2.V2  0.06087   1.06276  0.68969  0.088   0.9297  
interaction(Level1, Level2)M3.V2  0.97353   2.64726  0.65962  1.476   0.1400  
interaction(Level1, Level2)M1.V3 -0.07056   0.93187  0.88420 -0.080   0.9364  
interaction(Level1, Level2)M2.V3  0.25160   1.28608  0.63625  0.395   0.6925  
interaction(Level1, Level2)M3.V3  0.69068   1.99506  0.74816  0.923   0.3559  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

                                 exp(coef) exp(-coef) lower .95 upper .95
SexM                                0.7169     1.3949    0.3670     1.401
Age                                 0.9886     1.0116    0.9503     1.028
interaction(Level1, Level2)M2.V1    1.1746     0.8513    0.3055     4.517
interaction(Level1, Level2)M3.V1    0.8332     1.2002    0.2195     3.163
interaction(Level1, Level2)M1.V2    3.8510     0.2597    0.8320    17.824
interaction(Level1, Level2)M2.V2    1.0628     0.9409    0.2750     4.107
interaction(Level1, Level2)M3.V2    2.6473     0.3777    0.7267     9.644
interaction(Level1, Level2)M1.V3    0.9319     1.0731    0.1647     5.272
interaction(Level1, Level2)M2.V3    1.2861     0.7776    0.3696     4.476
interaction(Level1, Level2)M3.V3    1.9951     0.5012    0.4604     8.646

Concordance= 0.6  (se = 0.06 )
Likelihood ratio test= 8.13  on 10 df,   p=0.6
Wald test            = 8.63  on 10 df,   p=0.6
Score (logrank) test = 9.28  on 10 df,   p=0.5

我个人会选择使用交叉公式运算符*。这样你就可以分别得到“主效应”和交互系数。默认情况下，Level1 和 Level2 参考级别也是您所期望的。

dt %>%
     coxph(Surv(Time, status) ~Sex + Age + Level1 * Level2 , data = .) %>%
     summary()
#--------------
Call:
coxph(formula = Surv(Time, status) ~ Sex + Age + Level1 * Level2, 
    data = .)

  n= 100, number of events= 43 

                      coef exp(coef) se(coef)      z Pr(>|z|)  
SexM              -0.33279   0.71692  0.34166 -0.974   0.3300  
Age               -0.01150   0.98857  0.02012 -0.571   0.5677  
Level1M2           0.16097   1.17465  0.68718  0.234   0.8148  
Level1M3          -0.18250   0.83318  0.68056 -0.268   0.7886  
Level2V2           1.34832   3.85096  0.78176  1.725   0.0846 .
Level2V3          -0.07056   0.93187  0.88420 -0.080   0.9364  
Level1M2:Level2V2 -1.44842   0.23494  1.10566 -1.310   0.1902  
Level1M3:Level2V2 -0.19229   0.82506  0.96641 -0.199   0.8423  
Level1M2:Level2V3  0.16120   1.17491  1.07798  0.150   0.8811  
Level1M3:Level2V3  0.94374   2.56957  1.10329  0.855   0.3923  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

                  exp(coef) exp(-coef) lower .95 upper .95
SexM                 0.7169     1.3949    0.3670     1.401
Age                  0.9886     1.0116    0.9503     1.028
Level1M2             1.1746     0.8513    0.3055     4.517
Level1M3             0.8332     1.2002    0.2195     3.163
Level2V2             3.8510     0.2597    0.8320    17.824
Level2V3             0.9319     1.0731    0.1647     5.272
Level1M2:Level2V2    0.2349     4.2564    0.0269     2.052
Level1M3:Level2V2    0.8251     1.2120    0.1241     5.484
Level1M2:Level2V3    1.1749     0.8511    0.1420     9.718
Level1M3:Level2V3    2.5696     0.3892    0.2956    22.335

Concordance= 0.6  (se = 0.06 )
Likelihood ratio test= 8.13  on 10 df,   p=0.6
Wald test            = 8.63  on 10 df,   p=0.6
Score (logrank) test = 9.28  on 10 df,   p=0.5

这意味着对于特定级别值之间的效果比较将使用predict 函数完成。

【讨论】：

感谢@IRTFM .. 这也有效。 dt$Level2