Does anybody know how to arrange a cross validation so I can compare two models (negative binomial with quasi-Poisson)? I know some theory beyond cross validation, but don't know what kind of cross validation makes sense (leave one out, simple cross validation, k-fold) when comparing two glms with about 1000 observations. And how to do this using R?
Thanks for your help!
Here's an example of cross-validating poisson and quasi-poisson. You're on your own for negative-binomial.
#Setup
rm(list = ls(all = TRUE))
set.seed(1)
#Construct Data
counts <- c(18,17,15,20,10,20,25,13,12)
outcome <- as.numeric(gl(3,1,9))
treatment <- as.numeric(gl(3,3))
X <- data.frame(treatment, outcome)
#Create 10 Cross-Validation folds
library(caret)
tmp <- createResample(counts,times = 10)
myCtrl <- trainControl(method = "cv", index = tmp, timingSamps = 10)
#Run models
Pmodel <- train(X,counts,method='glm',trControl=myCtrl,family = poisson(link = "log"))
QPmodel <- train(X,counts,method='glm',trControl=myCtrl,family = quasipoisson(link = "log"))
#Assess Cross-Validated Error
resamps <- resamples(
list(Poisson = Pmodel,QuasiPoisson = QPmodel))
summary(resamps)
The example data is taken from the help page for glm. The tmp object stores the indexes of the cross-validation folds. You can write your own code to loop through the list of folds, fit a negative binomial model to the indexes in each fold, and then testing it on the indexes not in the fold.
Finally, you aggregate the error from all models and look at it's mean, medain, 95% CI, etc.
