R codes for Kaplan - Meier analysis in my youtube channel
# Raw data
# 6, 6, 6, 6*, 7, 9*, 10, 10*, 11*, 13, 16, 17*, 19*, 20*, 22, 23, 25*, 32*, 32*, 34*, 35*
library(survival); library(ggfortify); library(survminer)
## Loading required package: ggplot2
## Loading required package: ggpubr
time = c(6, 6, 6, 6, 7, 9, 10, 10, 11, 13, 16, 17, 19, 20, 22, 23, 25, 32, 32, 34, 35)
status = c(1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0)
dat = data.frame(time, status)
# Obtaining survival object with Surv function
sm = Surv(time, status)
# KM analysis - two functions are the same
km = survfit(sm ~ 1)
km = survfit(Surv(time, status) ~ 1, data=dat)
summary(km)
## Call: survfit(formula = Surv(time, status) ~ 1, data = dat)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 6 21 3 0.857 0.0764 0.720 1.000
## 7 17 1 0.807 0.0869 0.653 0.996
## 10 15 1 0.753 0.0963 0.586 0.968
## 13 12 1 0.690 0.1068 0.510 0.935
## 16 11 1 0.627 0.1141 0.439 0.896
## 22 7 1 0.538 0.1282 0.337 0.858
## 23 6 1 0.448 0.1346 0.249 0.807
# Using survminer package
ggsurvplot(km, risk.table=T, color = "#2E9FDF")
## Warning: Now, to change color palette, use the argument palette= '#2E9FDF'
## instead of color = '#2E9FDF'
# Plot the probability of event
ggsurvplot(km, risk.table=T, fun="event", color = "#2E9FDF")
## Warning: Now, to change color palette, use the argument palette= '#2E9FDF'
## instead of color = '#2E9FDF'
Test for difference between 2 survival curves
# Survival time (yrs) for 2 groups, each with 25 patients
time1 = c(12.3, 5.4, 8.2, 12.2, 11.7, 10.0, 5.7, 9.8, 2.6, 11.0, 9.2, 12.1, 6.6, 2.2, 1.8, 10.2, 10.7, 11.1, 5.3, 3.5, 9.2, 2.5, 8.7, 3.8, 3.0)
status1 = c(0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1)
time2 = c(5.8, 2.9, 8.4, 8.3, 9.1, 4.2, 4.1, 1.8, 3.1, 11.4, 2.4, 1.4, 5.9, 1.6, 2.8, 4.9, 3.5, 6.5, 9.9, 3.6, 5.2, 8.8, 7.8, 4.7, 3.9)
status2 = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1)
time = c(time1, time2)
status = c(status1, status2)
group = c(rep("1", 25), rep("2", 25))
dat = data.frame(time, status, group)
km = survfit(Surv(time, status) ~ group, data=dat)
km
## Call: survfit(formula = Surv(time, status) ~ group, data = dat)
##
## n events median 0.95LCL 0.95UCL
## group=1 25 22 8.7 5.4 10.7
## group=2 25 25 4.7 3.6 7.8
summary(km)
## Call: survfit(formula = Surv(time, status) ~ group, data = dat)
##
## group=1
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 1.8 25 1 0.96 0.0392 0.8862 1.000
## 2.2 24 1 0.92 0.0543 0.8196 1.000
## 2.5 23 1 0.88 0.0650 0.7614 1.000
## 2.6 22 1 0.84 0.0733 0.7079 0.997
## 3.0 21 1 0.80 0.0800 0.6576 0.973
## 3.5 20 1 0.76 0.0854 0.6097 0.947
## 3.8 19 1 0.72 0.0898 0.5639 0.919
## 5.3 18 1 0.68 0.0933 0.5197 0.890
## 5.4 17 1 0.64 0.0960 0.4770 0.859
## 5.7 16 1 0.60 0.0980 0.4357 0.826
## 6.6 15 1 0.56 0.0993 0.3956 0.793
## 8.2 14 1 0.52 0.0999 0.3568 0.758
## 8.7 13 1 0.48 0.0999 0.3192 0.722
## 9.2 12 2 0.40 0.0980 0.2475 0.646
## 9.8 10 1 0.36 0.0960 0.2135 0.607
## 10.0 9 1 0.32 0.0933 0.1807 0.567
## 10.2 8 1 0.28 0.0898 0.1493 0.525
## 10.7 7 1 0.24 0.0854 0.1195 0.482
## 11.0 6 1 0.20 0.0800 0.0913 0.438
## 11.1 5 1 0.16 0.0733 0.0652 0.393
## 11.7 4 1 0.12 0.0650 0.0415 0.347
##
## group=2
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 1.4 25 1 0.96 0.0392 0.88618 1.000
## 1.6 24 1 0.92 0.0543 0.81957 1.000
## 1.8 23 1 0.88 0.0650 0.76141 1.000
## 2.4 22 1 0.84 0.0733 0.70791 0.997
## 2.8 21 1 0.80 0.0800 0.65761 0.973
## 2.9 20 1 0.76 0.0854 0.60974 0.947
## 3.1 19 1 0.72 0.0898 0.56386 0.919
## 3.5 18 1 0.68 0.0933 0.51967 0.890
## 3.6 17 1 0.64 0.0960 0.47698 0.859
## 3.9 16 1 0.60 0.0980 0.43566 0.826
## 4.1 15 1 0.56 0.0993 0.39563 0.793
## 4.2 14 1 0.52 0.0999 0.35681 0.758
## 4.7 13 1 0.48 0.0999 0.31919 0.722
## 4.9 12 1 0.44 0.0993 0.28275 0.685
## 5.2 11 1 0.40 0.0980 0.24749 0.646
## 5.8 10 1 0.36 0.0960 0.21346 0.607
## 5.9 9 1 0.32 0.0933 0.18071 0.567
## 6.5 8 1 0.28 0.0898 0.14934 0.525
## 7.8 7 1 0.24 0.0854 0.11947 0.482
## 8.3 6 1 0.20 0.0800 0.09132 0.438
## 8.4 5 1 0.16 0.0733 0.06517 0.393
## 8.8 4 1 0.12 0.0650 0.04151 0.347
## 9.1 3 1 0.08 0.0543 0.02117 0.302
## 9.9 2 1 0.04 0.0392 0.00586 0.273
## 11.4 1 1 0.00 NaN NA NA
ggsurvplot(km, pval=T, risk.table=T)
## Warning: Vectorized input to `element_text()` is not officially supported.
## Results may be unexpected or may change in future versions of ggplot2.
