Как построить кривую выживаемости, вызванную выживанием (пакетная выживаемость R)?

Question

Я пытаюсь подогнать и построить модель Вейбулла для данных о выживании. Данные имеют только одну ковариацию, когорту, которая работает с 2006 по 2010 год. Итак, какие-либо идеи о том, что добавить к двум строкам кода, которые следует, чтобы построить кривую выживаемости когорты 2010 года?

library(survival) 
s <- Surv(subSetCdm$dur,subSetCdm$event) 
sWei <- survreg(s ~ cohort,dist='weibull',data=subSetCdm) 

Выполнение этого же с моделью Cox PH довольно просто, со следующими строками. Проблема в том, что survfit() не принимает объекты типа survreg.

sCox <- coxph(s ~ cohort,data=subSetCdm) 
cohort <- factor(c(2010),levels=2006:2010) 
sfCox <- survfit(sCox,newdata=data.frame(cohort)) 
plot(sfCox,col='green') 

Использование данных легкого (из пакета выживания), вот что я пытаюсь выполнить.

#create a Surv object 
s <- with(lung,Surv(time,status)) 

#plot kaplan-meier estimate, per sex 
fKM <- survfit(s ~ sex,data=lung) 
plot(fKM) 

#plot Cox PH survival curves, per sex 
sCox <- coxph(s ~ as.factor(sex),data=lung) 
lines(survfit(sCox,newdata=data.frame(sex=1)),col='green') 
lines(survfit(sCox,newdata=data.frame(sex=2)),col='green') 

#plot weibull survival curves, per sex, DOES NOT RUN 
sWei <- survreg(s ~ as.factor(sex),dist='weibull',data=lung) 
lines(survfit(sWei,newdata=data.frame(sex=1)),col='red') 
lines(survfit(sWei,newdata=data.frame(sex=2)),col='red') 

Answer 1

Надеется, что это помогает, и я не сделал какую-то вводит в заблуждении ошибки:

скопирован из выше:

#create a Surv object 
    s <- with(lung,Surv(time,status)) 

    #plot kaplan-meier estimate, per sex 
    fKM <- survfit(s ~ sex,data=lung) 
    plot(fKM) 

    #plot Cox PH survival curves, per sex 
    sCox <- coxph(s ~ as.factor(sex),data=lung) 
    lines(survfit(sCox,newdata=data.frame(sex=1)),col='green') 
    lines(survfit(sCox,newdata=data.frame(sex=2)),col='green') 

для Вейбулла, использование предсказать, повторно комментарий от Vincent:

#plot weibull survival curves, per sex, 
    sWei <- survreg(s ~ as.factor(sex),dist='weibull',data=lung) 

    lines(predict(sWei, newdata=list(sex=1),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red") 
    lines(predict(sWei, newdata=list(sex=2),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red") 

Фокус здесь заключался в том, чтобы изменить порядок квантилей для построения графика и прогноза. Вероятно, лучший способ сделать это, но он работает здесь. Удачи!

Answer 2

Альтернативный вариант - использовать пакет flexsurv. Это предлагает некоторые дополнительные функции по сравнению с пакетом survival - в том числе, что функция параметрической регрессии flexsurvreg() имеет прекрасный метод построения графика, который делает то, что вы просите.

Использование легких, как указано выше;

#create a Surv object 
s <- with(lung,Surv(time,status)) 

require(flexsurv) 
sWei <- flexsurvreg(s ~ as.factor(sex),dist='weibull',data=lung) 
sLno <- flexsurvreg(s ~ as.factor(sex),dist='lnorm',data=lung) 

plot(sWei) 
lines(sLno, col="blue") 

Вы можете построить по кумулятивной шкале опасности или опасности, используя type аргумент, и добавить доверительные интервалы с ci аргумента.

Answer 3

Это просто примечание уточнения Tim Riffe's answer, который использует следующий код:

lines(predict(sWei, newdata=list(sex=1),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red") 
lines(predict(sWei, newdata=list(sex=2),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red") 

Причина для двух последовательностей зеркальных, seq(.01,.99,by=.01) и seq(.99,.01,by=-.01), происходит потому, что предсказать() метод, дают квантили для распределение событий f (t), то есть значения обратного CDF f (t), а кривая выживаемости составляет 1 (CDF по f) по сравнению с t. Другими словами, если вы построите p по сравнению с прогнозом (p), вы получите CDF, и если вы построите 1-p против прогноза (p), вы получите кривую выживаемости, которая равна 1-CDF. Следующий код более прозрачен и обобщен на произвольные векторы значений р:

pct <- seq(.01,.99,by=.01) 
lines(predict(sWei, newdata=list(sex=1),type="quantile",p=pct),1-pct,col="red") 
lines(predict(sWei, newdata=list(sex=2),type="quantile",p=pct),1-pct,col="red")