Survival analysis in plant pathology

Nesi, C.N.; Shimakura, S.E.; Ribeiro Jr, P.J.; May De Mio, L.L.
IDESIA 31(2)

Data

R code for analysis

require(survival) source("Rfun.r") dados ← read.csv('data.csv', header=T, dec=",") dados$ano ← as.factor(dados$ano) dados$arvore ← as.factor(dados$arvore)

(s1 ← survfit(Surv(tempo, status)~ cult, data=subset(dados, ano=="1"))) (s2 ← survfit(Surv(tempo, status)~ cult, data=subset(dados, ano=="2"))) (s ← survfit(Surv(tempo, status)~ cult, data=dados))

survdiff(Surv(tempo, status)~ cult + strata(ano), data=subset(ensac, ano==1)) survdiff(Surv(tempo, status)~ cult + strata(ano), data=subset(ensac, ano==2)) survdiff(Surv(tempo, status)~ cult + strata(ano), data=ensac)

fit ← coxph(Surv(tempo, status) ~ strata(ano)+ cult + frailty(arvore), data=dados) summary(fit)

(residuo.sch ← cox.zph(fit))

### Plots

par(mfrow=c(2,2),mar=c(5,4,1,1), oma=c(1,1,1,1)) # kaplan meier par(mar=c(4,4,1,1)) plot(s, conf.int=F, lty=c(2:4), ylim=c(0.5, 1),

   xlab="Survival time (days)",
   font.lab=2, font=2, cex=1.5, ylab="Probability of survival",
   main="", lwd=c(2, 3, 4))

legend("bottomright", bty="n", cex=2, "a") legend(0.5, 0.7, c('Cultivar A','Cultivar B',' Cultivar C'),

     cex=1.2, bty="n",x.intersp=0.2, y.intersp=1,
     xjust=0, lty=2:4, lwd=c(2, 3, 4))

# frailty par(mar=c(4,4,1,1)) plot.frail(ensac$arvore, fit, xlab="Peach tree",

          font=2, font.lab=2, ylim=c(0.3,2.3),
         main="")

legend("bottomright", bty="n", cex=2, "b") legend("topleft", bty="n", cex=1.2, "Frailty (p = 0,17)")

# residuals

par(mar=c(4,4,1,1)) plot(residuo.sch, var="cultB",font=2, font.lab=2, main="",

   ylim=c(-21,21))

legend("bottomright", bty="n", cex=2, "c") legend("topleft", bty="n", cex=1.2, "r = 0,14 (p = 0,0925)")

par(mar=c(4,4,1,1)) plot(residuo.sch, var="cultC",font=2, font.lab=2, main="",

    ylim=c(-21,21))

legend("bottomright", bty="n", cex=2, "d") legend("topleft", bty="n", cex=1.2, "r = 0,12 (p = 0,1442)")