# PRACTICAL EXERCISE 9
# ====================

# Read cirrhosis data:
cirrhosis=read.table("http://www.uio.no/studier/emner/matnat/math/STK4080/h12/cirrhosis.txt",header=T)

# We convert time to years and introduce a binary covariate for ascites
cirrhosis$time=cirrhosis$time/365.25
cirrhosis$newasc=as.numeric(cirrhosis$asc>=1)

# We use the survival library, so this has to be loaded.

# Question a
# ----------
# We fit additive regression models for one covariate at a time.

# Treatment
fit.treat=aareg(Surv(time,status)~factor(treat),data=cirrhosis)
par(mfrow=c(1,2))
print(fit.treat)
plot(fit.treat)

# Sex
fit.sex=aareg(Surv(time,status)~factor(sex), data=cirrhosis)
par(mfrow=c(1,2))
print(fit.sex)
plot(fit.sex)


# Ascites
fit.newasc=aareg(Surv(time,status)~factor(newasc), data=cirrhosis)
par(mfrow=c(1,2))
print(fit.newasc)
plot(fit.newasc)


# Age (centered at the mean age of 59.4 years)
cirrhosis$cage=cirrhosis$age-mean(cirrhosis$age)
fit.cage=aareg(Surv(time,status)~cage, data=cirrhosis)
par(mfrow=c(1,2))
print(fit.cage)
plot(fit.cage)


# Prothrombin index (centered at the mean value of 69.3)
cirrhosis$cprot=cirrhosis$prot-mean(cirrhosis$prot)
fit.cprot=aareg(Surv(time,status)~cprot, data=cirrhosis)
par(mfrow=c(1,2))
print(fit.cprot)
plot(fit.cprot)


# Question b
# -----------
# We then fit a model with all covariates:
fit.all=aareg(Surv(time,status)~factor(treat)+factor(sex)+factor(newasc)+cage+cprot, data=cirrhosis)
par(mfrow=c(2,3))
print(fit.all)
plot(fit.all)


# Question c
# -----------
# We fit a model with all covariates and interaction between treatment and acites
fit.int=aareg(Surv(time,status)~factor(treat)+factor(sex)+factor(newasc)+factor(treat):factor(newasc)+cage+cprot, data=cirrhosis)
par(mfrow=c(2,4))
print(fit.int)
plot(fit.int)