생존분석을 할 때는 time 변수, status 변수(0,1로 코딩된 변수)를 이용하여 Surv클래스의 변수를 만들고 이를 이용하여 Kaplan-Meier plot을 그릴 수 있다. 전체집단을 특정변수를 이용하여 몇개의 군으로 나누어 생존을 비교할때는 log-rank test 를 실시한다. R 에서는 다음과 같이 한다. survival package에 있는 대장암 데이타에서 치료(rx)별로 생존을 비교할떄는 다음과 같이 할 수 있다.

require(survival)
Loading required package: survival
require(plyr)
Loading required package: plyr
source("ggkm.R")
require(maxstat)
Loading required package: maxstat
data(colon)

fit=survfit(Surv(time,status)~rx,data=colon)
plot(fit)

세 군별로 색깔과 선종류를 달리하고 time mark를 없애고 범례를 표시하려면 다음과 같이 하면 된다.

plot(fit,col=1:3,lty=1:3,mark.time=FALSE)
legend("topright",legend=levels(colon$rx),col=1:3,lty=1:3)

인터넷에 공개되어 있는 ggkm함수를 이용하여 그래프를 그리면 numbers at risk table과 log-rank test의 p값을 같이 그릴 수 있다.

source("ggkm.R")
ggkm(fit,timeby=500)

문제) 그런데 나이에 따른 생존을 보려면 어떻게 할까 ? 나이는 연속형 변수로 다음과 같은 분포를 가진다.

hist(colon$age)

그냥 아까와 같은 방법으로 plot을 하면 다음과 같이 된다.

fit=survfit(Surv(time,status)~age,data=colon)
plot(fit)

ggkm(fit,timeby=500)

즉 모든 나이 별로 다 그래프를 그리기 때문에 알아볼 수가 없다. 우리가 원하는 것은 연속형변수의 생존을 비교할때 생존에 차이가 나는 cutpoint 를 구해 cutpoint를 이용해 두개의 군으로 나누어 비교한 그래프를 그리고 싶다.

먼저 cutpoint를 결정하는 것은 maxstat 패키지의 maxstat.test() 를 이용하거나 party package의 ctree()를 이용할 수 있는데 maxstat.text를 이용하면 다음과 같다.

mstat=maxstat.test(Surv(time,status)~age,data=colon,
                   smethod="LogRank",pmethod="condMC",B=999)
mstat

    Maximally selected LogRank statistics using condMC

data:  Surv(time, status) by age
M = 2.4682, p-value = 0.1331
sample estimates:
estimated cutpoint 
                44 

여기서 cutpoint는 다음과 같이 구할 수 있다.

(cutpoint=mstat$estimate)
estimated cutpoint 
                44 

데이타프레임을 새로 만들어 cutpoint를 기준으로 cutpoint이하는 0으로 cutpoint 초과는 1로 코딩한다.

mydf1=colon
mydf1$temp=ifelse(mydf1$age<=cutpoint,0,1)

새로 survfit을 적용한 후 이를 이용해 그림을 그린다.

fit=survfit(Surv(time,status)~temp,data=mydf1)
labs=paste("age",c("≤",">"),cutpoint)
strataname=paste("age","group")
ggkm(fit,timeby=500,ystratalabs=labs,ystrataname=strataname)

웹에서 하는 R통계(web-r.org)에서는 이와 같은 방법으로 생존분석을 자동화하여 연속형변수는 cutpoint를 구해 자동으로 그래프를 그려준다.

또한 여러개의 변수를 영향변수(독립변수)로 넣어주면 각각의 변수에 대한 Cox비례위험 모형 및 각각의 survival curve를 그려주고 다변량 분석 및 stepwise backward elemination 까지 한번 입력으로 얻을 수 있으며 그 결과를 html또는 pdf로 다운받을 수 있고 그래프도 고해상도로 다운 받을 수 있다.