본문 바로가기

깜신의 통계 이야기

2018.5.23일 10차 통계워크샵 - 생존분석 강의록

R에서 생존분석 따라하기 코드입니다.


#### Survival Analysis #####

require(survival)

data("colon")

View(colon)

colon <- na.omit(colon)

rm(colon1)

str(colon)

colon$TS <- Surv(colon$time,colon$status==1)

fit = survfit(TS ~ rx, data = colon)

plot(fit)

plot(fit,col=1:3, lty=1:3)

legend("topright", legend=levels(colon$rx), col=1:3, lty = 1:3)


###### Cumulative hazard Start #####

plot(fit, col=1:3, lty=1:3, fun ="cumhaz", mark.time=FALSE, ylab="Cumulative hazard")

legend("topleft", legend = levels(colon$rx), col=1:3, lty =1:3)

###### Cumulative hazard End #######


#### Log-rank test ####

survdiff(Surv(time, status ==1)~rx, data=colon)


#### Cox Regression #####

out = coxph(Surv(time, status ==1)~rx, data=colon)

summary(out)


######  Hazard ratios of all individual variables #####

colon$TS <- Surv(colon$time,colon$status==1)

out = coxph(colon$TS~rx, data=colon )

install.packages("moonBook")

require(moonBook)

attach(colon)

out = mycph(TS~.-id-study-time-status-etype, data=colon)

out2 =coxph(TS~. -id-study-time-status-etype, data=colon)

final =step(out2, direction = "backward")

HRplot(out, type =2, show.CI = TRUE)