crossValidationNeRIFeatureSelection: NeRI-based selection of a linear, logistic, or Cox proportional hazards regression model from a set of candidate variables

Description

This function performs a cross-validation analysis of a feature selection algorithm based on net residual improvement (NeRI) to return a predictive model. It is composed of a NeRI-based feature selection followed by an update procedure, ending with a bootstrapping backwards feature elimination. The user can control how many train and blind test sets will be evaluated.

Usage

crossValidationNeRIFeatureSelection(size = 10,
	                                    fraction = 1.0,
	                                    pvalue = 0.05,
	                                    loops = 100,
	                                    covariates = "1",
	                                    Outcome,
	                                    timeOutcome = "Time",
	                                    variableList,
	                                    data,
	                                    maxTrainModelSize = 10,
	                                    type = c("LM", "LOGIT", "COX"),
	                                    testType = c("Binomial",
	                                                 "Wilcox",
	                                                 "tStudent",
	                                                 "Ftest"),
	                                    loop.threshold = 10,
	                                    startOffset = 0,
	                                    elimination.bootstrap.steps = 25,
	                                    trainFraction = 0.67,
	                                    trainRepetition = 9,
	                                    elimination.pValue = 0.05,
	                                    setIntersect = 1,
	                                    interaction = c(1,1),
	                                    update.pvalue = c(0.05,0.05),
	                                    unirank = NULL,
	                                    print=TRUE,
	                                    plots=TRUE)

Arguments

size

The number of candidate variables to be tested (the first size variables from variableList)

fraction

The fraction of data (sampled with replacement) to be used as train

pvalue

The maximum p-value, associated to the NeRI, allowed for a term in the model

loops

The number of bootstrap loops

covariates

A string of the type "1 + var1 + var2" that defines which variables will always be included in the models (as covariates)

Outcome

The name of the column in data that stores the variable to be predicted by the model

timeOutcome

The name of the column in data that stores the time to event (needed only for a Cox proportional hazards regression model fitting)

variableList

A data frame with two columns. The first one must have the names of the candidate variables and the other one the description of such variables

data

A data frame where all variables are stored in different columns

maxTrainModelSize

Maximum number of terms that can be included in the model

type

Fit type: Logistic ("LOGIT"), linear ("LM"), or Cox proportional hazards ("COX")

testType

Type of non-parametric test to be evaluated by the improvedResiduals function: Binomial test ("Binomial"), Wilcoxon rank-sum test ("Wilcox"), Student's t-test ("tStudent"), or F-test ("Ftest")

loop.threshold

After loop.threshold cycles, only variables that have already been selected in previous cycles will be candidates to be selected in posterior cycles

startOffset

Only terms whose position in the model is larger than the startOffset are candidates to be removed

elimination.bootstrap.steps

The number of bootstrap loops for the backwards elimination procedure

trainFraction

The fraction of data (sampled with replacement) to be used as train for the cross-validation procedure

setIntersect

The intersect of the model (To force a zero intersect, set this value to 0)

trainRepetition

The number of cross-validation folds (it should be at least equal to $1/trainFraction$ for a complete cross-validation)

elimination.pValue

The maximum p-value, associated to the NeRI, allowed for a term in the model by the backward elimination procedure

interaction

A vector of size two. The terms are used by the search and update procedures, respectively. Set to either 1 for first order models, or to 2 for second order models

update.pvalue

The maximum p-value, associated to the NeRI, allowed for a term in the model by the update procedure

unirank

A list with the results yielded by the uniRankVar function, required only if the rank needs to be updated during the cross-validation procedure

Logical. If TRUE, information will be displayed

plots

Logical. If TRUE, plots are displayed

Value

formula.listA list containing objects of class formula with the formulas used to fit the models found at each cycle
Models.testPredictionA data frame with the blind test set predictions made at each fold of the cross validation, where the models used to generate such predictions (formula.list) were generated via a feature selection process which included only the train set. It also includes a column with the Outcome of each prediction, and a column with the number of the fold at which the prediction was made.
FullModel.testPredictionA data frame similar to Models.testPrediction, but where the model used to generate the predictions was the full model, generated via a feature selection process which included all data.
backNeRIEliminationA list containing the values returned by bootstrapVarNeRIElimination using all data and the model from updateNeRISelection
varNeRISelectionA list containing the values returned by NeRIBasedFRESA.Model using all data
updateNeRISelectionA list containing the values returned by updateNeRIModel using all data and the model from varNeRISelection
testRMSEThe global blind test root-mean-square error (RMSE) of the cross-validation procedure
testPearsonThe global blind test Pearson r product-moment correlation coefficient of the cross-validation procedure
testSpearmanThe global blind test Spearman $\rho$ rank correlation coefficient of the cross-validation procedure
fulltestRMSEThe global blind test RMSE of the full model
fullTestPearsonThe global blind test Pearson r product-moment correlation coefficient of the full model
fullTestSpearmanThe global blind test Spearman $\rho$ rank correlation coefficient of the full model
trainRMSEThe train RMSE at each fold of the cross-validation procedure
trainPearsonThe train Pearson r product-moment correlation coefficient at each fold of the cross-validation procedure
trainSpearmanThe train Spearman $\rho$ rank correlation coefficient at each fold of the cross-validation procedure
fullTrainRMSEThe train RMSE of the full model at each fold of the cross-validation procedure
fullTrainPearsonThe train Pearson r product-moment correlation coefficient of the full model at each fold of the cross-validation procedure
fullTrainSpearmanThe train Spearman $\rho$ rank correlation coefficient of the full model at each fold of the cross-validation procedure
testRMSEAtFoldThe blind test RMSE at each fold of the cross-validation procedure
fullTestRMSEAtFoldThe blind test RMSE of the full model at each fold of the cross-validation procedure
fullenetAn object of class cv.glmnet containing the results of an elastic net cross-validation fit
enet.testPredictionsA data frame similar to Models.testPrediction, but where the predictions were made by the elastic net model
enetVariablesA list with the elastic net full model and the models found at each cross-validation fold
byFoldTestMSA vector with the Mean Square error for each blind fold
byFoldTestSpearmanA vector with the Spearman correlation between prediction and outcome for each blind fold
byFoldTestPearsonA vector with the Pearson correlation between prediction and outcome for each blind fold
byFoldCstatA vector with the C-index (Somers' Dxy rank correlation :rcorr.cens) between prediction and outcome for each blind fold
CVBlindPearsonA vector with the Pearson correlation between the outcome and prediction for each repeated experiment
CVBlindSpearmanA vector with the Spearm correlation between the outcome and prediction for each repeated experiment
CVBlindRMSA vector with the RMS between the outcome and prediction for each repeated experiment

Details

This function produces a set of data and plots that can be used to inspect the degree of over-fitting or shrinkage of a model. It uses bootstrapped data, cross-validation data, and, if possible, retrain data.

Examples

Run this code

# Start the graphics device driver to save all plots in a pdf format
	pdf(file = "Example.pdf")
	# Get the stage C prostate cancer data from the rpart package
	library(rpart)
	data(stagec)
	# Split the stages into several columns
	dataCancer <- cbind(stagec[,c(1:3,5:6)],
	                    gleason4 = 1*(stagec[,7] == 4),
	                    gleason5 = 1*(stagec[,7] == 5),
	                    gleason6 = 1*(stagec[,7] == 6),
	                    gleason7 = 1*(stagec[,7] == 7),
	                    gleason8 = 1*(stagec[,7] == 8),
	                    gleason910 = 1*(stagec[,7] >= 9),
	                    eet = 1*(stagec[,4] == 2),
	                    diploid = 1*(stagec[,8] == "diploid"),
	                    tetraploid = 1*(stagec[,8] == "tetraploid"),
	                    notAneuploid = 1-1*(stagec[,8] == "aneuploid"))
	# Remove the incomplete cases
	dataCancer <- dataCancer[complete.cases(dataCancer),]
	# Load a pre-stablished data frame with the names and descriptions of all variables
	data(cancerVarNames)
	# Rank the variables:
	# - Analyzing the raw data
	# - According to the NeRI
	rankedDataCancer <- univariateRankVariables(variableList = cancerVarNames,
	                                            formula = "Surv(pgtime, pgstat) ~ 1",
	                                            Outcome = "pgstat",
	                                            data = dataCancer,
	                                            categorizationType = "Raw",
	                                            type = "COX",
	                                            rankingTest = "NeRI",
	                                            description = "Description")
	# Get a Cox proportional hazards model using:
	# - The top 7 ranked variables
	# - 10 bootstrap loops in the feature selection procedure
	# - The Wilcoxon rank-sum test as the feature inclusion criterion
	# - 5 bootstrap loops in the backward elimination procedure
	# - A 5-fold cross-validation in the feature selection, 
	#           update, and backward elimination procedures
	# - First order interactions in the update procedure
	cancerModel <- crossValidationNeRIFeatureSelection(size = 7,
	                                                   loops = 10,
	                                                   Outcome = "pgstat",
	                                                   timeOutcome = "pgtime",
	                                                   variableList = rankedDataCancer,
	                                                   data = dataCancer,
	                                                   type = "COX",
	                                                   testType = "Wilcox",
	                                                   elimination.bootstrap.steps = 5,
	                                                   trainRepetition = 5,
	                                                   interaction = c(1,2))
	# Shut down the graphics device driver
	dev.off()

Run the code above in your browser using DataLab