# NOT RUN {
##########################
# Example with cancer data
library(survival)
head(cancer)
# Data preparation and convertion to 30 intervals
cancerPrep <- cancer
cancerPrep$status <- cancerPrep$status-1
intLim <- quantile(cancerPrep$time, prob=seq(0, 1, length.out=30))
intLim [length(intLim)] <- intLim [length(intLim)] + 1
# Cut discrete time in smaller number of intervals
cancerPrep <- contToDisc(dataSet=cancerPrep, timeColumn="time", intervalLimits=intLim)
# Generate training and test data
set.seed(753)
TrainIndices <- sample (x=1:dim(cancerPrep) [1], size=dim(cancerPrep) [1]*0.75)
TrainCancer <- cancerPrep [TrainIndices, ]
TestCancer <- cancerPrep [-TrainIndices, ]
TrainCancer$timeDisc <- as.numeric(as.character(TrainCancer$timeDisc))
TestCancer$timeDisc <- as.numeric(as.character(TestCancer$timeDisc))
# Convert to long format
LongTrain <- dataLong(dataSet=TrainCancer, timeColumn="timeDisc", censColumn="status")
LongTest <- dataLong(dataSet=TestCancer, timeColumn="timeDisc", censColumn="status")
# Convert factors
LongTrain$timeInt <- as.numeric(as.character(LongTrain$timeInt))
LongTest$timeInt <- as.numeric(as.character(LongTest$timeInt))
LongTrain$sex <- factor(LongTrain$sex)
LongTest$sex <- factor(LongTest$sex)
# Estimate, for example, a generalized, additive model in discrete survival analysis
library(mgcv)
gamFit <- gam (formula=y ~ s(timeInt) + s(age) + sex + ph.ecog, data=LongTrain, family=binomial())
summary(gamFit)
# 1. Specification of predicted discrete hazards
# Estimate survival function of each person in the test data
testPredHaz <- predict(gamFit, newdata=LongTest, type="response")
# 1.1. Calculate integrated prediction error
evalIntPredErr(hazPreds=testPredHaz, survPreds=NULL,
newTimeInput=TestCancer$timeDisc, newEventInput=TestCancer$status,
trainTimeInput=TrainCancer$timeDisc, trainEventInput=TrainCancer$status,
testDataLong=LongTest)
# 2. Alternative specification
oneMinusPredHaz <- 1-testPredHaz
predSurv <- aggregate(formula=oneMinusPredHaz ~ obj, data=LongTest, FUN=cumprod, na.action=NULL)
# 2.1. Calculate integrated prediction error
evalIntPredErr(survPreds=predSurv[[2]],
newTimeInput=TestCancer$timeDisc, newEventInput=TestCancer$status,
trainTimeInput=TrainCancer$timeDisc, trainEventInput=TrainCancer$status)
# }
Run the code above in your browser using DataLab