pdf("g02.pdf");options(width=64) #setwd("C:\\Users\\kolassa\\Class542") setwd("~/Taught1/960-542/Data")
# Block 1
#Effect of censoring on log rank statistic
library(PHInfiniteEstimates)#For checkcensor
checkcensor()
#*******************************************************/
# Larynx cancer data from Klelin and Moeschberger. See */
# class web page for link.  Cancers are graded on an   */
# integer scale, with increasing numbers indicating    */
# more sever.  Variables are stage, time to death, year*/
# of diagnosis, and status.                            */
#*******************************************************/
larynx<-read.table("larynx.txt",header=TRUE)
#***************************************************/
# Kaplan--Meier curve for the combined larynx data */
#***************************************************/
# The next command loads the appropriate package in R.
library(survival)#For Surv
#The inner function Surv takes a time variable and a
#censoring variable, and sets them up as input for a
#survival analysis.  survfit does the Kaplan-Meier.
sf.larynx<-survfit(Surv(time,delta)~1,
  data=larynx, conf.type="none")
cat("Summary of survival fit\n")
summary(sf.larynx)
cat("Print survival fit\n")
print(sf.larynx)
plot(sf.larynx,main="Larynx Tumor Survival",
   conf.int="none", xlab="Time (months)",ylab="Survival")
#***************************************************/
# Weaning data from Klein and Moeschberger.       */
# See class web page for link.  Variables are      */
#duration of breast feeding, weeks                 */
#Completed breast feeding? (1=yes, 0=no)           */
#Race of mother (1=white, 2=black, 3=other)        */
#Mother in poverty (1=yes, 0=no)                   */
#Mother smoked at birth of child (1=yes, 0=no)     */
#Mother used alcohol at birth of child (1=yes 0=no)*/
#Age of mother at birth of child                   */
#Year of birth                                     */
#Education level of mother (years of school)       */
#Prenatal care after 3rd month (1=yes, 0=no)       */
#***************************************************/
bfeed<-read.table("bfeed.txt",header=TRUE)
#****************************************************/
# Plot weaning times for smoking and non-smoking    */
# groups separately.                                */
#****************************************************/
plot(survfit(Surv(duration,delta)~strata(smoke),data=bfeed),
   main="Time Until Weaning, by Smoking Status",
   ylab="Probability of Not Weaning",xlab="Time(days)",
   lty=1:2)
legend("topright",lty=1:2, legend=c("Nonsmoker","Smoker"))
# Block 2
#***********************************************/
# Survival functions with Confidence Intervals */
#***********************************************/
# Default confidence interval in survit is "log".
# Do something more primitive here.
plot(survfit(Surv(duration,delta)~1,data=bfeed,
   conf.type="plain"),
   main="Time Until Weaning, All Mothers",
   sub="Confidence intervals, raw scale",
   ylab="Probability of Not Weaning",xlab="Time(days)")
cat('\n Larynx survival with conf. limits, raw scale\n')
#***************************************************/
# Kaplan--Meier curve for the combined larynx data */
#***************************************************/
sfp.larynx<-survfit(Surv(time,delta)~1,data=larynx,
   conf.type="plain")
plot(sfp.larynx,main="Larynx Tumor Survival",
   xlab="Time (months)",ylab="Survival",
   sub="No transformation")
# Application to high and low risk weeding data.;
smdrk<-bfeed[(bfeed$smoke==1)&(bfeed$alcohol==1),]
weanout<-survfit(Surv(duration,delta)~1,data=smdrk,
   conf.type="plain")
# Choose a small subset of mothers.  R truncates the
# confidence interval at zero, but it would be negative
# without truncation.
plot(weanout, sub="Plain Confidence Interval",
   xlab="Weaning time (days)",ylab="Survival",
   main="Weaning times for Mothers who Smoke and Drink")