[R-br] Extraindo o "cumulative hazard" de objeto survfit
Cesar Rabak
cesar.rabak em gmail.com
Sex Dez 10 20:38:38 -03 2021
Pedro,
A função plot.surfit não coloca os erros padrão da função cumulativa.
Ademais, pegue os resultados do summary e veja se se o cálculo do intervalo
do CI 95% "bate" com a conta 1,96 × d.p.
HTH
--
Cesar Rabak
On Fri, Dec 10, 2021 at 4:02 PM Pedro Emmanuel Alvarenga Americano do
Brasil <emmanuel.brasil em gmail.com> wrote:
> Na verdade eu estou fuçando ainda como plot.survfit faz pq não esta
> escrito em qualquer lugar que eu tenha achado como essa função faz e a
> summary.survfit retorna o erro padrão. Mas como a duvida era como extrair o
> cumhaz... achei que poderia postar. Se voce souber como a função plot faz
> posta aí que eu acerto.
>
> Pedro Brasil
>
>
> Em sex., 10 de dez. de 2021 às 13:32, Cesar Rabak <cesar.rabak em gmail.com>
> escreveu:
>
>> Esse cálculo dos limites superiores e inferiores do risco acumulado me
>> parece estranho..
>>
>> A formulação é baseada em alguma referência?
>>
>>
>> On Fri, Dec 10, 2021 at 1:23 PM Pedro Emmanuel Alvarenga Americano do
>> Brasil por (R-br) <r-br em listas.c3sl.ufpr.br> wrote:
>>
>>> > fit <- survfit(Surv(time, status) ~ x, data = aml)
>>> > y <- summary(fit, times = c(14,28,35))
>>> > y
>>> Call: survfit(formula = Surv(time, status) ~ x, data = aml)
>>>
>>> x=Maintained
>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>> 14 8 2 0.818 0.116 0.619 1.000
>>> 28 6 2 0.614 0.153 0.377 0.999
>>> 35 3 2 0.368 0.163 0.155 0.875
>>>
>>> x=Nonmaintained
>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>> 14 7 5 0.583 0.142 0.3616 0.941
>>> 28 4 2 0.389 0.147 0.1854 0.816
>>> 35 2 2 0.194 0.122 0.0569 0.664
>>>
>>> > cumhaz.summ <- function(y, digits = 3){
>>> + cond <- y$strata == levels(y$strata)[1]
>>> + y$cumhaz.upper <- pmax(y$cumhaz + 1.96 * y$std.chaz,0)
>>> + y$cumhaz.lower <- pmax(y$cumhaz - 1.96 * y$std.chaz,0)
>>> + output <- list()
>>> + for(i in seq_along(levels(y$strata))){
>>> + cond <- y$strata == levels(y$strata)[i]
>>> + output[[i]] <- round(data.frame(Time = y$time[cond],
>>> + 'N risk' = y$n.risk[cond],
>>> + 'N events' = y$n.event[cond],
>>> + 'Cum Hazard' = y$cumhaz[cond],
>>> + lower = y$cumhaz.lower[cond],
>>> + upper =
>>> y$cumhaz.upper[cond]),digits)
>>> + }
>>> + names(output) <- levels(y$strata)
>>> + output
>>> + }
>>> > cumhaz.summ(y)
>>> $`x=Maintained`
>>> Time N.risk N.events Cum.Hazard lower upper
>>> 1 14 8 2 0.191 0.000 0.456
>>> 2 28 6 2 0.459 0.002 0.915
>>> 3 35 3 2 0.909 0.133 1.685
>>>
>>> $`x=Nonmaintained`
>>> Time N.risk N.events Cum.Hazard lower upper
>>> 1 14 7 5 0.492 0.056 0.928
>>> 2 28 4 2 0.858 0.187 1.530
>>> 3 35 2 2 1.442 0.385 2.499
>>>
>>> >
>>>
>>> Em qui., 9 de dez. de 2021 às 10:43, Cid Póvoas por (R-br) <
>>> r-br em listas.c3sl.ufpr.br> escreveu:
>>>
>>>> library(survival)
>>>> library(broom)
>>>> library(tidyverse)
>>>>
>>>> fit <- survfit(Surv(time, status) ~ x, data = aml)
>>>> df<-tidy(fit)
>>>> df$cumhaz <- fit$cumhaz
>>>> df %>% filter(time%in%c(9,8,28,33,48))
>>>>
>>>> fit1 <- summary(fit, times = c(14,28,35))
>>>> fit1
>>>>
>>>> tab <- data.frame(time = fit1$time,
>>>> n.risk = fit1$n.risk,
>>>> n.event = fit1$n.event,
>>>> survival = fit1$surv,
>>>> std.err = fit1$std.err,
>>>> `lower 95% CI` = fit1$lower,
>>>> `upper 95% CI` = fit1$upper,
>>>> cumhaz = fit1$cumhaz,
>>>> strata = fit1$strata)
>>>>
>>>> tab
>>>>
>>>> *Cid Edson Mendonça Póvoas*
>>>>
>>>> *AnovAgro <http://www.anovagro.com/>*
>>>> *Engenheiro Agrônomo - **Data Scientist*
>>>> *CREA : 051984991-4*
>>>> *Técnico em Segurança do Trabalho *
>>>> *Nº: **0012669/BA*
>>>> *Tel: +55 73 99151-9565*
>>>> *Lattes : *http://lattes.cnpq.br/2303498368142537
>>>> *LinkedIn :* http://br.linkedin.com/in/cidedson/
>>>> *Whatsapp :* https://wa.me/5573991519565
>>>>
>>>>
>>>> Em qui., 9 de dez. de 2021 às 09:26, Pedro Emmanuel Alvarenga Americano
>>>> do Brasil por (R-br) <r-br em listas.c3sl.ufpr.br> escreveu:
>>>>
>>>>> Ei Cesar,
>>>>>
>>>>> Ei sei que o cumhaz esta la. Eu so queria uma saida parecida com a do
>>>>> survival. Parece que vou ter que montar uma um esquema aqui pegar esses
>>>>> valores e montar uma tabela de saida.
>>>>>
>>>>> Valeu.
>>>>>
>>>>> Pedro Emmanuel Brasil
>>>>> (:)=
>>>>>
>>>>> Em qua, 8 de dez de 2021 20:42, Cesar Rabak <cesar.rabak em gmail.com>
>>>>> escreveu:
>>>>>
>>>>>> Isto não é suficiente?
>>>>>>
>>>>>> > fit$cumhaz
>>>>>> [1] 0.09090909 0.19090909 0.31590909 0.45876623 0.45876623 0.65876623
>>>>>> [7] 0.90876623 0.90876623 1.40876623 1.40876623 0.16666667 0.36666667
>>>>>> [13] 0.49166667 0.49166667 0.65833333 0.85833333 1.10833333 1.44166667
>>>>>> [19] 1.94166667 2.94166667
>>>>>> >
>>>>>> HTH
>>>>>> --
>>>>>> Cesar Rabak
>>>>>>
>>>>>> On Wed, Dec 8, 2021 at 5:12 PM Pedro Emmanuel Alvarenga Americano do
>>>>>> Brasil por (R-br) <r-br em listas.c3sl.ufpr.br> wrote:
>>>>>>
>>>>>>> Saudações amigos do R,
>>>>>>>
>>>>>>> Estou as voltas de estimar taxas de eventos e estou batendo cabeça.
>>>>>>> Antes de escrever uma função eu mesmo para fazer uma tabela com alguns
>>>>>>> valores gostaria de uma luz dos amigos de R.
>>>>>>>
>>>>>>> O banco aml está no pacote survival então bastaria carregar o pacote
>>>>>>> pra reproduzir o exemplo. Eu gostaria de ver em formato de tabela alguns
>>>>>>> momentos específicos da tabela de sobrevivência. Só que o summary.survfit
>>>>>>> só retorna a sobrevivência e não a taxa cumulativa.
>>>>>>>
>>>>>>> library("survival")
>>>>>>> fit <- survfit(Surv(time, status) ~ x, data = aml)
>>>>>>> > fit
>>>>>>> Call: survfit(formula = Surv(time, status) ~ x, data = aml)
>>>>>>>
>>>>>>> n events median 0.95LCL 0.95UCL
>>>>>>> x=Maintained 11 7 31 18 NA
>>>>>>> x=Nonmaintained 12 11 23 8 NA
>>>>>>> # Summary com todos os momentos de eventos
>>>>>>> > summary(fit)
>>>>>>> Call: survfit(formula = Surv(time, status) ~ x, data = aml)
>>>>>>>
>>>>>>> x=Maintained
>>>>>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>>>>>> 9 11 1 0.909 0.0867 0.7541 1.000
>>>>>>> 13 10 1 0.818 0.1163 0.6192 1.000
>>>>>>> 18 8 1 0.716 0.1397 0.4884 1.000
>>>>>>> 23 7 1 0.614 0.1526 0.3769 0.999
>>>>>>> 31 5 1 0.491 0.1642 0.2549 0.946
>>>>>>> 34 4 1 0.368 0.1627 0.1549 0.875
>>>>>>> 48 2 1 0.184 0.1535 0.0359 0.944
>>>>>>>
>>>>>>> x=Nonmaintained
>>>>>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>>>>>> 5 12 2 0.8333 0.1076 0.6470 1.000
>>>>>>> 8 10 2 0.6667 0.1361 0.4468 0.995
>>>>>>> 12 8 1 0.5833 0.1423 0.3616 0.941
>>>>>>> 23 6 1 0.4861 0.1481 0.2675 0.883
>>>>>>> 27 5 1 0.3889 0.1470 0.1854 0.816
>>>>>>> 30 4 1 0.2917 0.1387 0.1148 0.741
>>>>>>> 33 3 1 0.1944 0.1219 0.0569 0.664
>>>>>>> 43 2 1 0.0972 0.0919 0.0153 0.620
>>>>>>> 45 1 1 0.0000 NaN NA NA
>>>>>>>
>>>>>>> # Summary com os momentos desejados
>>>>>>> > summary(fit, times = c(14,28,35))
>>>>>>> Call: survfit(formula = Surv(time, status) ~ x, data = aml)
>>>>>>>
>>>>>>> x=Maintained
>>>>>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>>>>>> 14 8 2 0.818 0.116 0.619 1.000
>>>>>>> 28 6 2 0.614 0.153 0.377 0.999
>>>>>>> 35 3 2 0.368 0.163 0.155 0.875
>>>>>>>
>>>>>>> x=Nonmaintained
>>>>>>> time n.risk n.event survival std.err lower 95% CI upper 95% CI
>>>>>>> 14 7 5 0.583 0.142 0.3616 0.941
>>>>>>> 28 4 2 0.389 0.147 0.1854 0.816
>>>>>>> 35 2 2 0.194 0.122 0.0569 0.664
>>>>>>>
>>>>>>> > plot(fit)
>>>>>>> > plot(fit, cumhaz = T)
>>>>>>> >
>>>>>>> Reparem que há uma opção para o gráfico cumhaz = T, o que significa
>>>>>>> que o cumhaz está depositado no objeto, inclusive dentro do summary também.
>>>>>>> Tipo summary(fit)$cumhaz. Só que não há uma opção summary(fit, cumhaz = T)
>>>>>>> que retorne o cumhaz ao invés da sobrevivência. Alguém tem algum bizu pra
>>>>>>> fazer isso organizado, extrair a mesma tabela só que com o cumhaz, como summary(fit,
>>>>>>> times = c(14,28,35)) sem muito trabalho?
>>>>>>>
>>>>>>> Abraço forte,
>>>>>>> Pedro Brasil
>>>>>>> _______________________________________________
>>>>>>> R-br mailing list
>>>>>>> R-br em listas.c3sl.ufpr.br
>>>>>>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>>>>>>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e
>>>>>>> forneça código mínimo reproduzível.
>>>>>>>
>>>>>> _______________________________________________
>>>>> R-br mailing list
>>>>> R-br em listas.c3sl.ufpr.br
>>>>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>>>>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
>>>>> código mínimo reproduzível.
>>>>>
>>>> _______________________________________________
>>>> R-br mailing list
>>>> R-br em listas.c3sl.ufpr.br
>>>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>>>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
>>>> código mínimo reproduzível.
>>>>
>>> _______________________________________________
>>> R-br mailing list
>>> R-br em listas.c3sl.ufpr.br
>>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
>>> código mínimo reproduzível.
>>>
>>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20211210/70b84523/attachment.htm>
Mais detalhes sobre a lista de discussão R-br