Learn R Programming

ModelGood (version 1.0.7)

RocBrierScore: Comparing prediction models with Receiver operating characteristics and Brier's scores

Description

Evaluating the performance of risk prediction models with binary status response variable (case/control or similar). Roc curves are based on a single marker, which is the probability predicted case probability by a model. The area under the curve and the Brier score is used to summarize the performance.

Bootstrap-crossvalidation techniques are implemented to estimate the generalization performance of the model(s), i.e. the performance which can be expected in new subjects.

Usage

Roc(object,...)
## S3 method for class 'list':
Roc(object,
                     formula,
                     data,
                     splitMethod="noSplitMethod",
                     noinf.method=c("simulate"),
                     simulate="reeval",
                     B,
                     M,
                     breaks,
                     crRatio=1,
                     RocAverageMethod="vertical",
                     RocAverageGrid=switch(RocAverageMethod,
		       "vertical"=seq(0,1,.01),
		       "horizontal"=seq(1,0,-.01)),
                     model.args=NULL,
                     model.parms=NULL,
                     keepModels=FALSE,
                     keepSampleIndex=FALSE,
                     keepCrossValRes=FALSE,
                     keepNoInfSimu,
                     slaveseed,
                     na.accept=0,
                     verbose=TRUE,
                     ...)
Brier(object,...)
## S3 method for class 'list':
Brier(object,
                     formula,
                     data,
                     splitMethod="noSplitMethod",
                     noinf.method=c("simulate"),
                     simulate="reeval",
                     crRatio=1,
                     B,
                     M,
                     model.args=NULL,
                     model.parms=NULL,
                     keepModels=FALSE,
                     keepSampleIndex=FALSE,
                     keepCrossValRes=FALSE,
                     na.accept=0,
                     verbose=TRUE,
                     ...)
## S3 method for class 'glm':
Brier(object,formula,data,...)
## S3 method for class 'lrm':
Brier(object,formula,data,...)
## S3 method for class 'rpart':
Brier(object,formula,data,...)
## S3 method for class 'randomForest':
Brier(object,formula,data,...)
## S3 method for class 'glm':
Roc(object,formula,data,...)
## S3 method for class 'lrm':
Roc(object,formula,data,...)
## S3 method for class 'rpart':
Roc(object,formula,data,...)
## S3 method for class 'randomForest':
Roc(object,formula,data,...)

Arguments

object
A named list of prediction models. Each entry is either an R-object for which a predictStatusProb method exists (see details) or a call that can be evaluated to such an R-object. For c
formula
A formula. If missing, use the formula of the (first) model in object. The left hand side is used to find the status response variable in data.
data
A data frame to validate the prediction models If missing, use the data of the (first) model in object.
splitMethod
Method for estimating the generalization error.

none:Assess the models in the same data where they are fitted. Yields the apparent or re-substitution performance. Often overestimates the generalization performance.

noinf.method
Method for computing the no-information Roc curve.
simulate
If equal to "reeval" then the models are re-evaluated in the current permuted data for computing the no-information Roc curve.
B
Number of bootstrap samples. The default depends on argument splitMethod. When splitMethod in c("Bootcv","Boot632","Boot632plus" the default is 100. For splitMethod="cvK" B is the number of cross
M
The size of the bootstrap samples for cross-validation without replacement.
breaks
Break points for computing the Roc curve. Defaults to seq(0,1,.01) for the Roc.list method and to sort(unique(breaks)) for the default method.
crRatio
Experimental and not yet tested. Cost/benefit ratio. Default value is to 1, meaning that misclassified cases are as bad as misclassified controls.
RocAverageMethod
Method for averaging the bootstrapped Roc curves. Either "vertical" or "horizontal". See reference 1 below for details.
RocAverageGrid
Grid points for the averaging of boostrapped Roc curves.
model.args
List of extra arguments that can be passed to the predictStatusProb methods. The list must have an entry for each entry in object.
model.parms
List with exactly one entry for each entry in object. Each entry names parts of the value of the fitted models that should be extracted and added to the output (see value).
keepModels
If FALSE keep only the names of the elements of object. If "Call" then keep the call of the elements of object. Else, add the object as it is to the output.
keepSampleIndex
Logical. If FALSE remove the cross-validation index (which tells who was in the learn and who in the validation set) from the output list which otherwise is included in the method part of the output list.
keepCrossValRes
Logical. If TRUE add all B crossvalidation results to the output (see value). Defaults to TRUE.
keepNoInfSimu
Logical. If TRUE add the B results in permuted data (for no-information performance) to the output (see value). Defaults to FALSE.
slaveseed
Vector of seeds, as long as B, to be given to the slaves in parallel computing.
na.accept
Works only for "Bootcv" estimate of performance. The maximal number of failures during training the models to the bootstrap samples. Usually a small number relative to B.
verbose
if TRUE the procedure is reporting details of the progress, e.g. it prints the current step in cross-validation procedures.
...
Difficult to explain

Value

  • Object of class Roc or class Brier for which print.Roc, summary.Roc, and plot.Roc (only Roc) methods are available. The object includes the following components:
  • RocThe Roc curve(s) estimated according to the splitMethod. A matrix where each column represents the estimated prediction error of a fit at the time points in time.
  • AppRocThe Roc curve(s) estimated when the model(s) are evaluated in the same data where they were fitted. Only if splitMethod is one of "NoInf", "Bootcv", "Boot632" or "Boot632plus", since otherwise repla is "apparent" and then this is stored in Roc as explained above.
  • BootcvRocThe prediction error when the model(s) are trained in the bootstrap sample and evaluated in the data that are not in the bootstrap sample. Only if splitMethod is one of "Boot632" or "Boot632plus". When splitMethod="Bootcv" then the BootcvRoc is stored in the component PredRoc.
  • NoInfRocThe prediction error when the model(s) are evaluated in the permuted data. Only if splitMethod is one of "Bootcv", "Boot632", or "Boot632plus". For splitMethod="NoInf" the NoInfRoc is stored in the component PredRoc.
  • weightThe weight used to linear combine the AppRoc and the BootcvRoc Only if splitMethod is one of "Boot632", or "Boot632plus".
  • overfitEstimated overfit of the model(s). See references. Only if splitMethod is one of "Boot632", or "Boot632plus".
  • callThe call that produced the object
  • modelsSee keepModels
  • methodThe method used for estimation of the overfitting bias.

Details

Missing data in the response or in the input matrix cause a failure. For each prediction model there must be a predictStatusProb method: for example, to assess a prediction model which evaluates to a myclass object one defines a function called predictStatusProb.myclass with arguments object,newdata,cutpoints,train.data,..., like this

myFit=myModel(Y~X,data=dat) class(myFit)="myclass"

predictStatusProb.myclass <- function(object,newdata,cutpoints,train.data,...){ predict(object, data=newdata,method="probabilities") out }

Such a function takes the object which was fitted with train.data and derives a matrix with predicted event status probabilities for each subject in newdata (rows) and each cutpoint (column) of the response variable that defines an event status. Currently implemented are predictStatusProb methods for the following R-functions: [object Object],[object Object],[object Object],[object Object]

References

Fawcett, T. (2006). An introduction to ROC analysis. Pattern Recognition Letters, 27, 861-874. Gerds, Cai & Schumacher (2008). The Performance of Risk Prediction Models. Biometrical Journal, Vol 50, 4, 457-479.

Efron, Tibshirani (1997) Journal of the American Statistical Association 92, 548--560 Improvement On Cross-Validation: The .632+ Bootstrap Method.

Wehberg, S and Schumacher, M (2004) A comparison of nonparametric error rate estimation methods in classification problems Biometrical Journal, Vol 46, 35--47

Examples

Run this code
## Generate some data with binary response Y
## depending on X1 and X2 and X1*X2 
N <- 400
X1 <- rnorm(N)
X2 <- rbinom(N,1,.4)
expit <- function(x) exp(x)/(1+exp(x))
lp <- expit(1 + X1 + X2 - X1*X2)
Y <- factor(rbinom(N,1,lp))
dat <- data.frame(Y=Y,X1=X1,X2=X2)

## fit a logistic model
lm1 <- glm(Y~X1,data=dat,family="binomial")
lm2 <- glm(Y~X1+X2,data=dat,family="binomial")
r1=Roc(list(lm1,lm2),verbose=0,crRatio=1)
summary(r1)
Brier(list(lm1,lm2),verbose=0,crRatio=1)


# crossing curves
set.seed(18)
N=500
Y=rbinom(N,1,.5)
X1=rnorm(N)
X1[Y==1]=rnorm(sum(Y==1),mean=rbinom(sum(Y==1),1,.5))
X2=rnorm(N)
X2[Y==0]=rnorm(sum(Y==0),mean=rbinom(sum(Y==0),1,.5))
dat <- data.frame(Y=Y,X1=X1,X2=X2)
lm1 <- glm(Y~X1,data=dat,family="binomial")
lm2 <- glm(Y~X2,data=dat,family="binomial")
plot(Roc(list(lm1,lm2),data=dat,verbose=0,crRatio=1))
library(randomForest)
dat$Y=factor(dat$Y)
rf <- randomForest(Y~X2,data=dat)
rocCV=Roc(list(RandomForest=rf,LogisticRegression=lm2),
  data=dat,
  verbose=TRUE,
  splitMethod="bootcv",
  B=10,
  crRatio=1)
plot(rocCV,diag=TRUE)

Run the code above in your browser using DataLab