생존분석 (survival analysis)는 의학 연구는 물론 기계의 내구성에 대한 판단 등 여러 분야에서 사용되는 통계 분석 방법입니다. 기본적인 개념은 사망이나 질병 발생, 기계 고장처럼 특정 현상이 발생하기까지 걸리는 시간을 분석하는 통계 분석 방법이라고 할 수 있습니다. 대표적인 방법은 카플란-메이어 (Kaplan-Meier) 생존분석과 콕스 비례위험 모형 (Cox proportional hazards model), 로그 순위법 (Log-rank test) 등이 있습니다.
이 방법들은 각각 적용 범위가 다른데 카플란-메이어 생존분석은 특정 집단의 시간에 따른 생존율을 분석하는 것이고 (예를 들어 암 수술후 재발 없이 지내는 시간이나 생존 기간) 로그 순위법은 두 집단의 생존율을 비교하는 것입니다. 콕스 비례위험 모형은 생존율에 영향을 미치는 위험 인자를 분석하는 것으로 예를 들어 허리둘레가 1cm씩 증가할 때 마다 당뇨 위험도가 얼마나 증가하는지를 보는 것입니다.
여기서는 각각의 통계 분석 방법에 대한 상세한 설명은 생략하고 R을 통해 어떻게 분석할 수 있는지를 간단히 설명합니다. 우선 콕스 비례 위험 모형에 대해서 설명해 봅니다. 비례위험모형은 기저 분포에 대한 가정 (예를 들어 정규 분포)이 필요없고 여러 공변량을 모형에 넣어 보정해 사용할 수 있어 다양한 경우에 사용될 수 있기 때문에 생존분석에서 가장 널리 사용되는 방식 중 하나입니다.
R에서 생존분석을 하기 위해 가장 기본적으로 사용되는 패키지는 survival 패키지 입니다. survival 패키지에는 여러 가지 예제 데이터가 있는데 우선 kidney 데이터를 살펴보겠습니다.
require(survival)
data(kidney)
str(kidney)
kidney 데이터에 대한 설명은 위의 링크나 ?kidney를 통해 볼 수 있습니다.
patient: id
time: time
status: event status
age: in years
sex: 1=male, 2=female
disease: disease type (0=GN, 1=AN, 2=PKD, 3=Other)
frail: frailty estimate from original paper
Note
The original paper ignored the issue of tied times and so is not exactly reproduced by the survival package.
Source
CA McGilchrist, CW Aisbett (1991), Regression with frailty in survival analysis. Biometrics 47, 461–66.
이 데이터는 신장염 환자들이 투석을 위해 도뇨관 (catheter)를 삽입한 후 감염이 재발하기까지의 시간을 나타낸 것입니다. 감염 이외에 다른 요인으로 도뇨관을 제거하는 경우에는 중도 절단으로 계산하며 이 경우가 0, 감염이 관찰된 경우 1로 status에 코딩되어 있습니다. id는 환자 아이디로 한 환자에서 두 개의 도뇨관이 삽입되었다는 것을 알 수 있습니다.
모든 분석에서 그렇듯이 예제를 수행하기 전 kidney 데이터를 한 번 보고서 어떤 데이터인지 대략 확인하는 것이 좋습니다. 다만 의학 쪽 전공자가 아니라면 다소 이해가 되지 않을 데이터도 있습니다. 예를 들어 신장염의 종류의 경우 GN은 사구체신염 (glomerulonephritis)라는 뜻입니다. 카테터는 삽입 후 시간이 지나면서 감염이 생기기 때문에 결국 다시 제거해야 합니다. 하지만 그 전에 다른 문제가 생겨 중단될수도 있습니다. 다만 기본적으로 오래전 데이터라 현재의 술기와 좀 차이가 있는 듯 합니다.
이 예제에서는 나이와 성별이 카테터 감염 재발 위험도에 미치는 영향을 알아보기로 하겠습니다. 기본적으로 coxph 기능을 사용하는데, ?coxph를 통해 기본적인 문법을 알 수 있습니다. (coxph(Surv(사건 발생까지 시간 변수, status 는 상태를 코딩한 변수) ~ 알고자 하는 공변량, 방법, 데이터) 생존시간에 동점이 있는 경우 기본적인 처리 방법은 efron입니다.
kd1<-coxph data="kidney)</span" status="=1)~age," ties="efron" time="" urv="">
summary(kd1)
> kd1<-coxph data="kidney)</span" status="=1)~age," ties="efron" time="" urv="">
> summary(kd1)
Call:
coxph(formula = Surv(time, status == 1) ~ age, data = kidney,
ties = "efron")
n= 76, number of events= 58
coef exp(coef) se(coef) z Pr(>|z|)
age 0.004575 1.004586 0.008957 0.511 0.609
exp(coef) exp(-coef) lower .95 upper .95
age 1.005 0.9954 0.9871 1.022
Concordance= 0.503 (se = 0.046 )
Rsquare= 0.003 (max possible= 0.993 )
Likelihood ratio test= 0.26 on 1 df, p=0.6
Wald test = 0.26 on 1 df, p=0.6
Score (logrank) test = 0.26 on 1 df, p=0.6
여기서 핵심적인 결과는 exp(coef) 입니다. HR 값을 의미하기 때문입니다. 이 결과를 보면 나이가 한 살 더 많아지면 위험도는 1.005 (95% CI 0.987 - 1.022)배 증가하는 것으로 보입니다. 통계적 유의성이 없는 것으로 판단할 수 있습니다. 그렇다면 성별은 어떨까요?
kd2<-coxph data="kidney)</span" status="=1)~factor(sex)," ties="efron" time="" urv="">
summary(kd2)
> kd2<-coxph data="kidney)</span" status="=1)~factor(sex)," ties="efron" time="" urv="">
> summary(kd2)
Call:
coxph(formula = Surv(time, status == 1) ~ factor(sex), data = kidney,
ties = "efron")
n= 76, number of events= 58
coef exp(coef) se(coef) z Pr(>|z|)
factor(sex)2 -0.8377 0.4327 0.2966 -2.824 0.00474 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
factor(sex)2 0.4327 2.311 0.2419 0.7738
Concordance= 0.636 (se = 0.031 )
Rsquare= 0.089 (max possible= 0.993 )
Likelihood ratio test= 7.07 on 1 df, p=0.008
Wald test = 7.98 on 1 df, p=0.005
Score (logrank) test = 8.4 on 1 df, p=0.004
결과를 보면 여성인 경우 위험도는 0.4327배 (95% CI 0.2419 - 0.7738)인 것으로 나타납니다. 여성에서 카테터 감염이 생길 비례위험이 훨씬 작다는 의미입니다. P값은 0.00474로 역시 유의합니다. 이번에는 모델에 나이와 성별 모두 넣어 보겠습니다.
kd3<-coxph data="kidney)</span" factor="" sex="" status="=1)~age" ties="efron" time="" urv="">
summary(kd3)
> kd3<-coxph data="kidney)</span" factor="" sex="" status="=1)~age" 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.046 )
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
역시 비슷한 결과이지만, 그래도 나이가 보정된 성별 비례 위험도를 알 수 있습니다. 이런식으로 모델에 여러 변수를 넣어 서로 보정한 결과를 알 수 있습니다. survival 패키지는 여러 가지 기능을 제공하는데, 그 중에서는 각 관측치의 위험률 추정치나 회귀식 추정치 등을 보여주는 predict 기능도 있습니다.
> predict(kd3, type = "risk")
[1] 1.7845712 1.7845712 0.8109912 0.8109912 1.7991344 1.7991344 0.7834563 0.7850498 1.7204816 1.7204816 0.7599381 0.7614838 1.8699495 1.8699495 0.8226085 0.8242817
[17] 0.8463448 0.8463448 1.8699495 1.8737529 0.8044266 0.8044266 0.7882465 0.7882465 0.7898498 0.7898498 0.8011642 0.8011642 0.7614838 0.7614838 1.9044598 1.9044598
[33] 0.8310083 0.8310083 0.8027937 0.8044266 0.8192724 0.8192724 0.8044266 0.8044266 1.8510481 1.8548131 0.7818660 0.7818660 0.8343922 0.8360893 0.8011642 0.8027937
[49] 1.8397991 1.8397991 0.8259582 0.8276382 0.7507297 0.7507297 0.8176094 0.8176094 1.8775640 1.8775640 0.8209388 0.8209388 0.8242817 0.8242817 0.8142936 0.8159498
[65] 0.8259582 0.8259582 0.8044266 0.8060627 0.7692594 0.7692594 0.8011642 0.8011642 0.8176094 0.8176094 1.9044598 1.9044598
76개 관측치 각각의 비례 위험도를 보여주는데, 어떤 사람에서 특히 위험도가 높은지를 추정할 수 있습니다.
댓글
댓글 쓰기