library(survival)  
args(coxph)
Rossi <- read.table("http://www.stat.unm.edu/~ghuerta/stat574/Rossi.txt", header=T)
Rossi[1:5, 1:10]
mod.allison <- coxph(Surv(week, arrest) ~ fin + age + race + wexp + mar + paro + prio,
    data=Rossi)
summary(mod.allison)

plot(survfit(mod.allison), ylim=c(.7, 1), xlab='Weeks', 
    ylab='Proportion Not Rearrested')

 # diagnostics
mod.allison.2 <- coxph(Surv(week, arrest) ~ fin + age + prio,
    data=Rossi)
summary(mod.allison.2)

    # proportional hazards
    
cox.zph(mod.allison.2)
par(mfrow=c(2,2))
plot(cox.zph(mod.allison.2))

   # influence
    
dfbeta <- residuals(mod.allison.2, type='dfbeta')
par(mfrow=c(2,2))
for (j in 1:3) {
    plot(dfbeta[,j], ylab=names(coef(mod.allison.2))[j])
    abline(h=0, lty=2)
    }

    # nonlinearity
        
res <- residuals(mod.allison.2, type='martingale')
X <- as.matrix(Rossi[,c("age", "prio")]) # matrix of covariates
par(mfrow=c(1,2))
for (j in 1:2) {  # residual plots
    plot(X[,j], res, xlab=c("age", "prio")[j], ylab="residuals")
    abline(h=0, lty=2)
    lines(lowess(X[,j], res, iter=0))
    }

b <- coef(mod.allison.2)[c(2,3)]  # regression coefficients
for (j in 1:2) {  # partial-residual plots
    plot(X[,j], b[j]*X[,j] + res, xlab=c("age", "prio")[j], 
        ylab="component+residual")
    abline(lm(b[j]*X[,j] + res ~ X[,j]), lty=2)
    lines(lowess(X[,j], b[j]*X[,j] + res, iter=0))
    }