Translate

2020년 3월 4일 수요일

생존분석- 콕스 비례위험모형 (6)





 콕스 비례위험 모형 역시 여러 가지 모델을 지닐 수 있습니다. 이 가운데 가장 좋은 모형이 어떤 것인지를 알아내는 과정은 쉽지 않지만, 앞서 다른 분석 방법에서 본 것처럼 몇 가지 판단할 수 있는 방법이 있습니다. 콕스 비례위험 모형 역시 AIC (Akaike information criterion)를 기준으로 가장 적합한 모형을 찾는 방법이 있습니다. kidney 데이터에서 이를 찾아보겠습니다. 




library(survival)
library(MASS)

kd3<-coxph data="kidney)</span" disease="" factor="" sex="" status="=1)~age" ties="efron" time="" urv="">
summary(kd3)

model_selection<- kd3="" model="" selection="" span="" stepaic="" stepwise="" trace="TRUE)">
model_selection


 stepAIC를 통해 stepwise model selection이 가능한데, 당연히 AIC가 가장 작은 모델을 선택합니다. 만약 과정 없이 최종 모델만 보고 싶다면 trace + FALSE로 지정하면 됩니다. 


> kd3<-coxph data="kidney)</span" disease="" factor="" sex="" status="=1)~age" ties="efron" time="" urv="">
> summary(kd3)
Call:
coxph(formula = Surv(time, status == 1) ~ age + factor(sex) + 
    factor(disease), data = kidney, ties = "efron")

  n= 76, number of events= 58 

                        coef exp(coef)  se(coef)      z Pr(>|z|)    
age                 0.003181  1.003186  0.011146  0.285   0.7754    
factor(sex)2       -1.483137  0.226925  0.358230 -4.140 3.47e-05 ***
factor(disease)GN   0.087957  1.091941  0.406369  0.216   0.8286    
factor(disease)AN   0.350794  1.420195  0.399717  0.878   0.3802    
factor(disease)PKD -1.431108  0.239044  0.631109 -2.268   0.0234 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

                   exp(coef) exp(-coef) lower .95 upper .95
age                   1.0032     0.9968   0.98151    1.0253
factor(sex)2          0.2269     4.4067   0.11245    0.4579
factor(disease)GN     1.0919     0.9158   0.49238    2.4216
factor(disease)AN     1.4202     0.7041   0.64880    3.1088
factor(disease)PKD    0.2390     4.1833   0.06939    0.8235

Concordance= 0.697  (se = 0.041 )
Rsquare= 0.207   (max possible= 0.993 )
Likelihood ratio test= 17.65  on 5 df,   p=0.003
Wald test            = 19.9  on 5 df,   p=0.001
Score (logrank) test = 20.13  on 5 df,   p=0.001

> model_selection<- kd3="" model="" selection="" span="" stepaic="" stepwise="" trace="TRUE)">
Start:  AIC=368.16
Surv(time, status == 1) ~ age + factor(sex) + factor(disease)

                  Df    AIC
- age              1 366.24
               368.16
- factor(disease)  3 372.69
- factor(sex)      1 381.12

Step:  AIC=366.24
Surv(time, status == 1) ~ factor(sex) + factor(disease)

                  Df    AIC
               366.24
- factor(disease)  3 370.74
- factor(sex)      1 379.15
> model_selection
Call:
coxph(formula = Surv(time, status == 1) ~ factor(sex) + factor(disease), 
    data = kidney, ties = "efron")

                      coef exp(coef) se(coef)      z        p
factor(sex)2       -1.4774    0.2282   0.3569 -4.140 3.48e-05
factor(disease)GN   0.1392    1.1494   0.3635  0.383   0.7017
factor(disease)AN   0.4132    1.5116   0.3360  1.230   0.2188
factor(disease)PKD -1.3671    0.2549   0.5889 -2.321   0.0203

Likelihood ratio test=17.56  on 4 df, p=0.001501
n= 76, number of events= 58 


 여기서는 sex와 disease만 넣은 변수가 가장 좋은 모델로 나왔습니다. 하지만 만약 우리가 진짜 알려는 독립 변수는 age이고 sex와 disease는 보정변수이기만 하다면 age를 뺄 수 없습니다. 그리고 아무래도 age가 중요한 영향인자라면 주된 독립 변수가 아니라도 뺄 수 없습니다. 선형 회귀나 로지스틱 회귀 모형과 마찬가지로 참고일 뿐이지 모델을 결정하는 절대적 방법은 아니지만, 여러 개의 변수가 존재할 때 참고할 수 있는 방법임은 분명합니다. 


 콕스 비례위험 모형은 이름처럼 각 변수의 위험도가 시간에 따라 비례해 증가하거나 감소한다고 가정하고 있습니다. 그러나 몇몇 변수들은 실제로는 그렇지 않을 가능성이 있습니다. 이럴 때 가장 간단한 해결책은 비례 가정을 만족시키지 않는 변수를 뺀 모형을 만드는 것이지만, 앞서 설명한 것처럼 그 변수가 빼기 어려운 중요한 보정 변수인 경우 고민에 빠지게 됩니다. 이 경우 층화(strata)가 한 가지 방법이 될 수 있습니다. 기본 개념은 비례 정도가 다른 집단에서 (예를 들어 남자와 여자, 노인과 청년 등)을 여러 개의 층으로 만들어 하나의 HR 값을 찾는 것입니다. 


 여기서는 age가 보려는 변수인데 주요 보정 변수인 sex가 비례 위험 가정을 만족시키지 않는 것처럼 보인다고 가정하겠습니다. 이 경우 sex를 모델에서 빼기 어렵기 때문에 성별에 따른 각각의 값을 구하거나 혹은 층화 콕스 비례 위험모형을 적용해야 합니다. strata()를 통해 적용할 수 있습니다. 


 #strata

kd3<-coxph age="" data="kidney)</span" factor="" sex="" status="=1)~" ties="efron" time="" urv="">
summary(kd3)

kd4<-coxph age="" data="kidney)</span" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
summary(kd4)

> kd3<-coxph age="" data="kidney)</span" factor="" sex="" status="=1)~" ties="efron" time="" urv="">
> summary(kd3)
Call:
coxph(formula = Surv(time, status == 1) ~ age + factor(sex), 
    data = kidney, ties = "efron")

  n= 76, number of events= 58 

                  coef exp(coef)  se(coef)      z Pr(>|z|)   
age           0.002032  1.002034  0.009246  0.220  0.82607   
factor(sex)2 -0.829314  0.436349  0.298955 -2.774  0.00554 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

             exp(coef) exp(-coef) lower .95 upper .95
age             1.0020      0.998    0.9840     1.020
factor(sex)2    0.4363      2.292    0.2429     0.784

Concordance= 0.662  (se = 0.045 )
Rsquare= 0.089   (max possible= 0.993 )
Likelihood ratio test= 7.12  on 2 df,   p=0.03
Wald test            = 8.02  on 2 df,   p=0.02
Score (logrank) test = 8.45  on 2 df,   p=0.01

> kd4<-coxph age="" data="kidney)</span" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
> summary(kd4)
Call:
coxph(formula = Surv(time, status == 1) ~ age + strata(sex), 
    data = kidney, ties = "efron")

  n= 76, number of events= 58 

        coef exp(coef) se(coef)     z Pr(>|z|)
age 0.008007  1.008039 0.009253 0.865    0.387

    exp(coef) exp(-coef) lower .95 upper .95
age     1.008      0.992    0.9899     1.026

Concordance= 0.547  (se = 0.048 )
Rsquare= 0.01   (max possible= 0.982 )
Likelihood ratio test= 0.77  on 1 df,   p=0.4
Wald test            = 0.75  on 1 df,   p=0.4
Score (logrank) test = 0.75  on 1 df,   p=0.4


 결과를 보면 층화를 한 경우와 아닌 경우 모두 통계적 유의성이 없지만, 비례 계수가 1.002에서 1.008로 바뀌는 것을 볼 수 있습니다. 만약 sex가 비례 위험 가정을 만족하지 않는다면 층화를 통해 이 문제를 극복하고 더 정확한 값을 구할 수 있습니다. strata()를 적용한 변수에 대해서는 HR 값이 얻어지지 않는다는 점도 주목해야 합니다. 당연히 값이 나올 수 없습니다. 다만 이 결과를 믿기 위해서는 sex와 age간의 중요한 교호작용(상호작용)이 있는지를 검증해야 합니다. 


kd5<-coxph age="" data="kidney)</span" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
summary(kd5)

> kd5<-coxph age="" data="kidney)</span" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
> summary(kd5)
Call:
coxph(formula = Surv(time, status == 1) ~ age * strata(sex), 
    data = kidney, ties = "efron")

  n= 76, number of events= 58 

                         coef exp(coef) se(coef)      z Pr(>|z|)
age                  -0.01339   0.98670  0.01692 -0.791    0.429
age:strata(sex)sex=2  0.02893   1.02935  0.02009  1.440    0.150

                     exp(coef) exp(-coef) lower .95 upper .95
age                     0.9867     1.0135    0.9545     1.020
age:strata(sex)sex=2    1.0294     0.9715    0.9896     1.071

Concordance= 0.559  (se = 0.047 )
Rsquare= 0.035   (max possible= 0.982 )
Likelihood ratio test= 2.74  on 2 df,   p=0.3
Wald test            = 2.68  on 2 df,   p=0.3
Score (logrank) test = 2.71  on 2 df,   p=0.3


 P = 0.15로 유의한 교호작용은 없으며 층화 모형이 의미 있다고 할 수 있습니다. 만약 다른 변수도 층화를 통해 구해야 한다면 같은 과정을 반복하면 됩니다. 


kd3<-coxph age="" data="kidney)</span" disease="" factor="" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
summary(kd3)

kd4<-coxph age="" data="kidney)</span" disease="" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
summary(kd4)

kd5<-coxph age="" data="kidney)</span" disease="" sex="" status="=1)~" strata="" ties="efron" time="" urv="">
summary(kd5)


 그리고 만약 연속변수가 비례 위험 모형을 만족시키지 않는 경우 변수를 몇 개의 그룹으로 쪼개 문제를 피할 수도 있습니다. 예를 들어 나이 60세를 기준으로 나누는 등입니다. 그래도 비례 위험 가정이 만족되지 않으면 역시 층화를 시도하거나 그룹을 나눠 분석할 수 있습니다. 

댓글 없음:

댓글 쓰기