Указание пользовательских временных точек для графика выживания - PullRequest
3 голосов
/ 22 апреля 2020

Я работаю над созданием графика выживания / кумулятивного события, используя функцию ggsurvplot из пакета survminer. Я хочу указать собственные временные точки для своего графика, но не могу понять, как это сделать. Параметры xlim и break.x.by помогают, но они создают равномерно распределенные временные точки и больше временных точек, чем я хочу. В частности, мне нужны только временные точки в 0, 30 и 365 дней. Возможно ли это сделать с помощью пакета survminer? Если возможно почти повторить вывод из ggsurvplot, используя другой пакет, это тоже сработает - я действительно хочу получить информацию о числе риска и совокупном количестве событий. Большое спасибо!

Пример данных:

data <- structure(list(status = c(1, 1, 1, 1, 1, 0, 0, 1, 0, 1, 0, 1, 
0, 0, 1, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 
1, 1, 0, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 0, 1, 
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0, 
0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 
0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 
0, 0, 1, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 
0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 
0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 0, 
0, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 
1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 
0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 
0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 
0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 
0, 0, 0, 0, 0, 0), time_to_trt = c(48, 202, 83, 235, 195, 158, 
81, 251, 237, 408, 794, 82, 734, 697, 208, 385, 96, 384, 126, 
277, 102, 95, 651, 475, 522, 119, 51, 248, 202, 76, 606, 129, 
675, 399, 26, 969, 76, 491, 91, 68, 261, 185, 164, 395, 93, 176, 
389, 84, 67, 245, 299, 500, 487, 347, 325, 125, 102, 62, 403, 
496, 298, 264, 249, 167, 146, 67, 311, 175, 283, 654, 599, 371, 
172, 367, 151, 277, 234, 81, 54, 419, 405, 244, 152, 948, 332, 
174, 451, 347, 872, 436, 738, 699, 449, 578, 41, 517, 332, 56, 
214, 214, 172, 48, 343, 23, 157, 874, 668, 375, 201, 55, 633, 
112, 43, 245, 23, 720, 60, 405, 799, 173, 446, 363, 333, 393, 
335, 294, 46, 48, 776, 237, 435, 132, 216, 58, 266, 50, 333, 
171, 816, 844, 201, 135, 24, 833, 398, 354, 269, 573, 51, 391, 
48, 713, 206, 148, 538, 52, 294, 139, 647, 655, 8, 25, 486, 458, 
370, 314, 224, 841, 933, 191, 601, 648, 755, 727, 717, 194, 684, 
37, 577, 580, 580, 399, 462, 378, 148, 71, 273, 265, 190, 160, 
46, 58, 323, 370, 181, 377, 291, 254, 535, 113, 129, 440, 202, 
523, 155, 556, 284, 266, 945, 281, 914, 293, 175, 805, 78, 327, 
82, 77, 606, 376, 292, 168, 110, 124, 98, 170, 83, 25, 18, 54, 
26, 561, 106, 45, 528, 41, 341, 259, 102, 277, 591, 256, 165, 
354, 53, 356, 391, 221, 127, 444, 69, 188, 377, 54, 874, 851, 
252, 42, 762, 76, 54, 96, 315, 347, 313, 109, 74, 231, 283, 223, 
237, 194, 172, 321, 217, 384, 486, 466, 111, 105, 378, 140, 129, 
207, 110, 518, 60, 278, 252, 363, 213, 356, 167, 45, 438, 308, 
95, 374, 363, 252, 46, 172, 139, 52, 92, 34, 158, 87, 87, 207, 
95, 175), trtsdt = structure(c(17605, 17885, 18074, 18065, 18025, 
NA, NA, 18156, NA, 17906, NA, 17688, NA, NA, 17883, 18135, 17932, 
NA, 18191, NA, NA, NA, NA, 18281, NA, 18106, 18171, NA, NA, NA, 
NA, 18169, 18225, 18058, 18171, NA, 17702, 18164, 17780, 17814, 
18102, NA, NA, 18085, 17963, 18326, 17820, 17634, 17687, 17834, 
17974, NA, NA, 18282, NA, NA, NA, NA, 18282, NA, NA, NA, NA, 
NA, NA, NA, 18103, 17800, 17935, NA, NA, 18291, 18128, NA, 18246, 
NA, NA, NA, NA, NA, NA, NA, 17744, NA, 17819, 17826, NA, NA, 
NA, 18003, NA, NA, 18131, NA, 17842, NA, NA, 18187, NA, NA, NA, 
NA, 18116, 17871, NA, NA, NA, NA, 18248, NA, 18018, 17603, 17617, 
17927, 17626, NA, 17940, NA, 18220, 17624, 17989, 18038, 18058, 
18184, NA, NA, NA, NA, NA, 17856, 18304, 18045, NA, 17562, 17910, 
17779, NA, 18351, NA, NA, 18247, 17598, 17507, NA, 18351, NA, 
NA, 18071, 17619, NA, NA, NA, 18001, 18227, 18274, 17980, NA, 
NA, 18225, NA, 17757, 17854, NA, NA, NA, NA, 18303, 18200, NA, 
17745, 18172, 18233, NA, NA, NA, 17883, NA, 17771, NA, NA, NA, 
18284, NA, NA, 18137, 18156, NA, NA, NA, NA, NA, 17506, 18066, 
NA, 18311, NA, NA, 18059, NA, 18122, NA, NA, NA, 18010, 17723, 
NA, NA, 17682, NA, 17736, NA, 17827, 17729, NA, 17729, 18009, 
17834, 17841, NA, 18198, 18205, 18109, 18114, 18198, 18219, 18319, 
NA, 18317, 17576, 17625, 17681, 18199, 17879, 17914, NA, 17983, 
18332, NA, 17673, 18037, NA, 18065, 18121, NA, 18121, 18088, 
NA, 17652, 17638, NA, 17660, 18080, NA, NA, NA, NA, 17806, 17652, 
NA, 17759, 17869, 17925, 18233, NA, NA, NA, NA, 18305, NA, NA, 
NA, NA, NA, 18060, 18025, 18235, NA, NA, 18200, NA, NA, 18109, 
18249, NA, 17932, NA, 18121, NA, NA, NA, NA, NA, NA, NA, NA, 
18226, NA, NA, NA, 18296, 18093, NA, NA, NA, 18199, 18241, NA, 
NA, NA, NA, NA, NA), label = "Index Procedure Date", class = "Date", format.sas = "DATE"), 
    randdt = structure(c(17557, 17683, 17991, 17830, 17830, 18214, 
    18291, 17905, 18135, 17498, 17578, 17606, 17638, 17675, 17675, 
    17750, 17836, 17988, 18065, 18095, 18270, 18277, 17721, 17806, 
    17850, 17987, 18120, 18124, 18170, 18296, 17766, 18040, 17550, 
    17659, 18145, 17403, 17626, 17673, 17689, 17746, 17841, 18187, 
    18208, 17690, 17870, 18150, 17431, 17550, 17620, 17589, 17675, 
    17872, 17885, 17935, 18047, 18247, 18270, 18310, 17879, 17876, 
    18074, 18108, 18123, 18205, 18226, 18305, 17792, 17625, 17652, 
    17718, 17773, 17920, 17956, 18005, 18095, 18095, 18138, 18291, 
    18318, 17953, 17967, 18128, 17592, 17424, 17487, 17652, 17921, 
    18025, 17500, 17567, 17634, 17673, 17682, 17794, 17801, 17855, 
    18040, 18131, 18158, 18158, 18200, 18324, 17773, 17848, 18215, 
    17498, 17704, 17997, 18047, 18317, 17385, 17491, 17574, 17682, 
    17603, 17652, 17880, 17967, 17421, 17451, 17543, 17675, 17725, 
    17791, 18037, 18078, 18326, 18324, 17596, 17619, 17869, 17913, 
    18156, 17504, 17644, 17729, 18039, 18180, 17556, 17528, 18046, 
    17463, 17483, 17539, 17953, 18018, 18103, 17498, 17568, 17981, 
    18324, 17659, 17795, 18079, 17736, 17928, 18078, 18233, 17578, 
    17717, 17749, 17829, 17886, 17914, 18002, 18058, 18079, 17359, 
    17439, 17554, 17571, 17585, 17617, 17645, 17655, 17689, 17688, 
    17734, 17795, 17792, 17792, 17885, 17910, 17994, 17989, 18085, 
    18099, 18107, 18182, 18212, 18326, 17448, 17743, 18002, 18130, 
    17995, 18081, 17805, 17837, 18009, 18243, 17932, 18170, 17487, 
    17568, 17816, 18088, 17416, 17427, 17455, 17458, 17534, 17554, 
    17567, 17651, 17682, 17752, 17764, 17766, 17822, 17913, 17941, 
    18004, 18074, 18121, 18149, 18289, 18292, 17558, 17571, 17655, 
    17638, 17773, 17869, 17844, 17942, 17991, 18113, 17571, 17760, 
    17781, 17809, 17956, 18018, 18068, 17732, 17981, 17431, 17511, 
    17928, 17591, 17892, 17995, 18318, 17498, 17521, 17554, 17610, 
    17610, 17683, 17815, 17829, 17918, 18025, 18059, 18263, 18298, 
    18074, 18089, 18149, 18135, 18178, 18200, 17739, 17808, 17851, 
    17886, 17906, 18089, 18267, 17994, 17969, 18120, 18165, 17822, 
    17854, 18061, 18094, 18120, 18009, 18159, 18016, 18205, 18327, 
    17934, 17918, 18277, 17998, 18009, 18044, 18047, 18200, 18233, 
    18320, 18107, 18207, 18214, 18285, 18285, 18165, 18277, 18197
    ), label = "Randomization Date", class = "Date", format.sas = "DATE")), class = c("tbl_df", 
"tbl", "data.frame"), row.names = c(NA, -312L))

Fit:

fit <- survfit(Surv(time_to_trt, status) ~ 1, data = data)

«Плохая» попытка с использованием ggsurvplot:

a <- ggsurvplot(fit, 
           data = data, 
           conf.int = FALSE,
           risk.table = TRUE,
           cumevents = TRUE,
           fun = "event",
           ggtheme = theme_minimal(),
           risk.table.y.text = FALSE,
           risk.table.y.text.col = TRUE,
           legend = "none",
           legend.title = "")

Вывод: enter image description here

1 Ответ

3 голосов
/ 22 апреля 2020

Вот одно из возможных решений. Это что-то вроде монстра Франкенштейна, но я думаю, что это хорошо показывает, что происходит за кулисами в функции ggsurvplot. Пожалуйста, дайте мне знать, если это не то, что вы имели в виду. Я мог бы порекомендовать вам использовать принципы из этого ответа, чтобы создать свой график напрямую, используя только функции ggplot, но я оставлю это на ваше усмотрение

ggsurvplot возвращает 3 ggplot объектов как часть список: скаттер и 2 таблицы. Вы можете извлечь эти объекты, отредактировать их напрямую и вернуть обратно в список. Сначала я генерирую график с указанным вами xlim и с break.x.by = 30

a <- ggsurvplot(fit, 
                data = data, 
                conf.int = FALSE,
                risk.table = TRUE,
                cumevents = TRUE,
                fun = "event",
                ggtheme = theme_minimal(),
                risk.table.y.text = FALSE,
                risk.table.y.text.col = TRUE,
                legend = "none",
                legend.title = "",
                xlim = c(0, 360),
                break.x.by = 30)
a

enter image description here

Затем я вручную редактирую каждый из ggplot объектов для установки разрывов по оси X. Две таблицы генерируются с geom_text. Если вы хотите удалить метки, которые не соответствуют тикам, для меня проще всего было стереть метки, а затем вручную воссоздать их с правильными значениями x.

# extract ggplot object from ggsurvplot
p <- a$plot 
p <- p + scale_x_continuous(breaks = c(0, 30, 360))

# extract table object from ggsurvplot
tab <- a$table
tab$layers = NULL # clear labels
tab <- tab + 
  geom_text(aes(x = time, y = rev(strata), label = llabels), data = tab$data[tab$data$time %in% c(0, 30, 360),]) +
  scale_x_continuous(breaks = c(0, 30, 360))

# extract cumevents object from ggsurvplot
tab2 <- a$cumevents
tab2$layers = NULL # clear labels
tab2 <- tab2 + 
  geom_text(aes(x = time, y = rev(strata), label = cum.n.event), data = tab$data[tab$data$time %in% c(0, 30, 360),]) +
  scale_x_continuous(breaks = c(0, 30, 360))


# Add plots back
a$plot <- p
a$table <- tab
a$cumevents <- tab2

a

enter image description here

...