Last data update: 2014.03.03

R: Predict method for iCoxBoost fits
predict.iCoxBoostR Documentation

Predict method for iCoxBoost fits

Description

Obtains predictions at specified boosting steps from a iCoxBoost object fitted by iCoxBoost.

Usage

## S3 method for class 'iCoxBoost'
predict(object,newdata=NULL,
        subset=NULL,at.step=NULL,times=NULL,
        type=c("lp","logplik","risk","CIF"),...)

Arguments

object

fitted CoxBoost object from a CoxBoost call.

newdata

data frame with new covariate values (for n.new observations). If just prediction for the training data is wanted, it can be omitted. If the predictive partial log-likelihood is wanted (type=logplik), this frame also has to contain the response information.

subset

an optional vector specifying a subset of observations to be used for evaluation.

at.step

scalar or vector of boosting step(s) at which prediction is wanted. If type="risk" is used, only one step is admissible. If no step is given, the final boosting step is used.

times

vector with T time points where prediction is wanted. Only needed for type="risk"

type

type of prediction to be returned: "lp" gives the linear predictor, "logplik" the partial log-likelihood, "risk" the predicted probability of not yet having had the event at the time points given in times, and "CIF" the predicted cumulative incidence function, i.e., the predicted probability of having had the event of interest.

...

miscellaneous arguments, none of which is used at the moment.

Value

For type="lp" and type="logplik" a vector of length n.new (at.step being a scalar) or a n.new * length(at.step) matrix (at.step being a vector) with predictions is returned. For type="risk" or type="CIF" a n.new * T matrix with predicted probabilities at the specific time points is returned.

Author(s)

Harald Binder binderh@uni-mainz.de

Examples

n <- 200; p <- 100
beta <- c(rep(1,2),rep(0,p-2))
x <- matrix(rnorm(n*p),n,p)
actual.data <- as.data.frame(x)
real.time <- -(log(runif(n)))/(10*exp(drop(x %*% beta)))
cens.time <- rexp(n,rate=1/10)
actual.data$status <- ifelse(real.time <= cens.time,1,0)
actual.data$time <- ifelse(real.time <= cens.time,real.time,cens.time)

#   define training and test set

train.index <- 1:100
test.index <- 101:200

#   Fit a Cox proportional hazards model by iCoxBoost

cbfit <- iCoxBoost(Surv(time,status) ~ .,data=actual.data[train.index,],
				   stepno=300,cv=FALSE)

#   mean partial log-likelihood for test set in every boosting step

step.logplik <- predict(cbfit,newdata=actual.data[test.index,],
                        at.step=0:300,type="logplik")

plot(step.logplik)

#   names of covariates with non-zero coefficients at boosting step
#   with maximal test set partial log-likelihood

print(coef(cbfit,at.step=which.max(step.logplik)-1))

Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(CoxBoost)
Loading required package: survival
Loading required package: Matrix
Loading required package: prodlim
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/CoxBoost/predict.iCoxBoost.Rd_%03d_medium.png", width=480, height=480)
> ### Name: predict.iCoxBoost
> ### Title: Predict method for iCoxBoost fits
> ### Aliases: predict.iCoxBoost
> ### Keywords: models regression survial
> 
> ### ** Examples
> 
> n <- 200; p <- 100
> beta <- c(rep(1,2),rep(0,p-2))
> x <- matrix(rnorm(n*p),n,p)
> actual.data <- as.data.frame(x)
> real.time <- -(log(runif(n)))/(10*exp(drop(x %*% beta)))
> cens.time <- rexp(n,rate=1/10)
> actual.data$status <- ifelse(real.time <= cens.time,1,0)
> actual.data$time <- ifelse(real.time <= cens.time,real.time,cens.time)
> 
> #   define training and test set
> 
> train.index <- 1:100
> test.index <- 101:200
> 
> #   Fit a Cox proportional hazards model by iCoxBoost
> 
> ## No test: 
> cbfit <- iCoxBoost(Surv(time,status) ~ .,data=actual.data[train.index,],
+ 				   stepno=300,cv=FALSE)
> ## End(No test)
> 
> #   mean partial log-likelihood for test set in every boosting step
> 
> ## No test: 
> step.logplik <- predict(cbfit,newdata=actual.data[test.index,],
+                         at.step=0:300,type="logplik")
> 
> plot(step.logplik)
> ## End(No test)
> 
> #   names of covariates with non-zero coefficients at boosting step
> #   with maximal test set partial log-likelihood
> 
> ## No test: 
> print(coef(cbfit,at.step=which.max(step.logplik)-1))
         V1          V2          V3          V4          V5          V6 
 1.09211021  0.84349939  0.00000000  0.00000000  0.00000000  0.00000000 
         V7          V8          V9         V10         V11         V12 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V13         V14         V15         V16         V17         V18 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 -0.01223960 
        V19         V20         V21         V22         V23         V24 
 0.00000000  0.00000000  0.03937006  0.00000000  0.01192784  0.00000000 
        V25         V26         V27         V28         V29         V30 
 0.05123715 -0.04547263  0.00000000 -0.09812249 -0.05796674  0.00000000 
        V31         V32         V33         V34         V35         V36 
 0.08551553  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V37         V38         V39         V40         V41         V42 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V43         V44         V45         V46         V47         V48 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V49         V50         V51         V52         V53         V54 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V55         V56         V57         V58         V59         V60 
 0.00000000  0.00000000 -0.04847894  0.00000000  0.00000000  0.00000000 
        V61         V62         V63         V64         V65         V66 
 0.00000000  0.00000000  0.00000000  0.00000000  0.04079100  0.00000000 
        V67         V68         V69         V70         V71         V72 
 0.00000000  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V73         V74         V75         V76         V77         V78 
 0.01317528  0.00000000  0.00000000  0.00000000  0.00000000  0.00000000 
        V79         V80         V81         V82         V83         V84 
 0.00000000  0.00000000  0.06976195  0.00000000  0.00000000  0.00000000 
        V85         V86         V87         V88         V89         V90 
 0.09730206  0.00000000  0.00000000  0.02507007  0.00000000  0.00000000 
        V91         V92         V93         V94         V95         V96 
 0.00000000  0.00000000  0.00000000  0.00000000  0.01194670  0.00000000 
        V97         V98         V99        V100 
 0.00000000  0.00000000  0.00000000  0.00000000 
> ## End(No test)
> 
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>