#loading the survival library
library(survival)

# Fig. 1.1, p. 6
# Note: The first two numbers are the coordinates of the upper left corner of the legend.

windows(6,6)
print( plot(hmohiv$age[hmohiv$censor==1], hmohiv$time[hmohiv$censor==1], xlab="age", 
       ylab="time" , pch=2, ylim=c(0, 65), xlim=c(15,55), cex=.7, cex.lab=.7, cex.axis=.7) )
points(hmohiv$age[hmohiv$censor==0], hmohiv$time[hmohiv$censor==0], pch=18, cex=.7)
legend(40, 60, c("Censor=1", "Censor=0") , pch=c(2,18), cex=.7 ) 

# Fig. 1.2, p. 6

hmohiv$age1 <- 1000/hmohiv$age

windows(6,6)
print( plot(hmohiv$age1[hmohiv$censor==1], hmohiv$time[hmohiv$censor==1], xlab="1000/Age", 
       ylab="time" , pch=2, ylim=c(0, 65), xlim=c(15,55), cex=.7, cex.lab=.7, cex.axis=.7) )
points(hmohiv$age1[hmohiv$censor==0], hmohiv$time[hmohiv$censor==0], pch=18, cex=.7)
legend(15, 60, c("censor=1", "censor=0") , pch=c(2,18), cex=.7 ) 

# table 1.2, p. 14
# Note: In order to get the complete output we must first create the object "test" which is the 
# fitted survival regression and then get the summary of the object by using the summary command.

test <- survreg(Surv(time, censor) ~ age, hmohiv, dist="exponential") 
print(summary(test))


# The list of all the distribution available in R.

print(names(survreg.distributions))

# Fig. 1.3, p. 16

hmohiv$p <- predict(test, type="response"  )

# sort by age 

hmohiv$age.1<-hmohiv$age[order(hmohiv$age)]
hmohiv$time.1<-hmohiv$time[order(hmohiv$age)]
hmohiv$censor.1<-hmohiv$censor[order(hmohiv$age)]
hmohiv$p.1<-hmohiv$p[order(hmohiv$age)]

windows(6,6)
print( plot(hmohiv$age.1[hmohiv$censor.1==1], hmohiv$time.1[hmohiv$censor.1==1], xlab="Age", 
       ylab="time" , pch=2, ylim=c(0, 65), xlim=c(15,55), cex=.7, cex.lab=.7, cex.axis=.7) )
points(hmohiv$age.1[hmohiv$censor.1==0], hmohiv$time.1[hmohiv$censor.1==0], pch=18, cex=.7)
lines(hmohiv$age.1, hmohiv$p.1)
legend(40, 60, c("Censor=1", "Censor=0") , pch=c(2,18) , cex=.7)