#  R code to fit the Cox model to data 
#  from 50 rat litters. This code is posted as  
#        rats.ssc   

# Variables
# Litter number: 
#     even litter numbers are males, 
#     odd litter numbers are female
#  Exposure: 0=no, 1=yes
#  Time: Follow-up time
#  Status: 0-tumor, 1-censored (due to animal's death) 

# Enter the data into a data frame.

rats <- read.table("c:/documents and settings/kkoehler.IASTATE/My Documents/courses/st565/sas/therneau/rats.dat",
      header=F, col.names=c("litter", "rx", 
        "time", "status")) 

library(survival)

fitm <- coxph(Surv(time, status) ~ rx +  
                cluster(litter), data = rats,
                robust=T, method="efron")
summary(fitm)


# Fit a frailty model

fitf <- coxph(Surv(time, status) ~ rx +  
                frailty(litter, distribution="gamma"), 
                data = rats, method="efron")
summary(fitf)

# Fit a frailty model

fitg <- coxph(Surv(time, status) ~ rx +  
                frailty(litter, distribution="gaussian"), 
                data = rats, method="efron")
summary(fitg)