0
$\begingroup$

My data have 119 cases and we did ROC for x (continuous variable) to predict postoperative y (categorical variable) available here, we got a comment from a reviewer asking:

Please provide statistical evidence that the AUC was not overfitted to the model. With N=119, C-stat = 0.81 seems optimistic. Optimism-adjusted?

I found some answers: basically I split my cohort into train set (n=107) and test set (n=12), did 2 ROC curves and compared AUC for both using this command:

testobj <- roc.test(rocobj0, rocobj1);testobj 

and this was the output

DeLong's test for two ROC curves

data: rocobj0 and rocobj1
D = 0.44822, df = 12.488, p-value = 0.6617
alternative hypothesis: true difference in AUC is not equal to 0
sample estimates:
AUC of roc1 AUC of roc2
80.26961 72.22222

Regarding the optimism adjustment, I did it using the following code copied from a blog post on cainarchaeology.

library (pROC);library (kimisc) 
  auc.adjust <- function(data, fit, B){
  fit.model <- fit
  data$pred.prob <- fitted(fit.model)
  auc.app <- roc(data\[,1\], data$pred.prob, data=data)$auc # require 'pROC'
  auc.boot <- vector (mode = "numeric", length = B)
  auc.orig <- vector (mode = "numeric", length = B)
  o <- vector (mode = "numeric", length = B)
  for(i in 1:B){    
    boot.sample <- sample.rows(data, nrow(data), replace=TRUE) # require 'kimisc'
    fit.boot <- glm(formula(fit.model), data = boot.sample, family = "binomial")
    boot.sample$pred.prob <- fitted(fit.boot)
    auc.boot\[i\] <- roc(boot.sample\[,1\], boot.sample$pred.prob, data=boot.sample)$auc
    data$pred.prob.back <- predict.glm(fit.boot, newdata=data, type="response")
    auc.orig\[i\] <- roc(data\[,1\], data$pred.prob.back, data=data)$auc
    o\[i\] <- auc.boot\[i\] - auc.orig\[i\]
  }
  auc.adj <- auc.app - (sum(o)/B)
  boxplot(auc.boot, auc.orig, names=c("auc.boot", "auc.orig"))
  title(main=paste("Optimism-adjusted AUC", "\nn of bootstrap resamples:", B), sub=paste("auc.app (blue line)=", round(auc.app, digits=4),"\nadj.auc (red line)=", round(auc.adj, digits=4)), cex.sub=0.8)
  abline(h=auc.app, col="blue", lty=2)
  abline(h=auc.adj, col="red", lty=3)
}

model <- glm(data$y ~  data$x, data = data, family = "binomial")
auc.adjust(data, model, B=200)

But it gave me this weird figure below in which AUC is > 1:

**Weird AUC >1 after bootstrapping and AUC =1 originally!!**

  1. I should conclude that AUC was not overfitted based on the splitting of the data that I did and P value of 0.66 from DeLong's test for two ROC curves, Right?
  2. How to fix the optimism adjusted AUC figure to get it <=1?
$\endgroup$
1
$\begingroup$

A single test/train breakdown, your first suggestion, is not going to work well for evaluating overfitting in your situation with only 119 cases. Separating out one test and one training set doesn't usually work well unless there are thousands of cases, as otherwise the results depend too heavily on the vagaries both of the sample from the population and the particular choice of test set out of your sample.

The second approach, based on bootstrapping, is well accepted and works well with limited data sets like yours. For each bootstrap sample you develop a model and compare the AUC on the bootstrap sample against the AUC from applying that model to the original data set. The difference is a measure of the optimism of the model that was developed on the bootstrap sample. The idea is that bootstrap sampling from your data set is like taking the original data set from the full population, so this optimism (averaged over multiple bootstrap samples) estimates the optimism in your original model.

It's not clear why your code produced such strange looking results, and coding questions per se are not on topic here. (I do note that in one place you call your data set data and in other places df.) You can avoid coding problems by using well documented statistical functions. The rms package in R, for example, provides an lrm() function to do logistic regression and a validate() function that evaluates optimism (of several measures) directly on lrm output and provides optimism-adjusted values. It might not show AUC directly, but you can get that from the $D_{xy}$ value it reports: AUC = $0.5 + (D_{xy}/2)$. The rms package is from Frank Harrell, acknowledged in the blog post you cite as the originator of this approach to evaluating and correcting for optimism.

$\endgroup$
3
  • $\begingroup$ Thanks for your invaluable input. I edited previously provided code here but it was OK in my lab computer. I followed your advice and here is the used code and output $\endgroup$ Aug 15 '19 at 19:50
  • $\begingroup$ @MohamedRahouma coding-specific questions and troubleshooting of code aren't on topic on this site. From the statistics perspective of this site, you should go with the bootstrapping-based estimate of optimism rather than the train/test approach. You can get your optimism estimate with established tools, in the rms package for example, and avoid the risks of such coding errors. $\endgroup$
    – EdM
    Aug 15 '19 at 20:15
  • $\begingroup$ Thanks a lot for your invaluable help. $\endgroup$ Aug 15 '19 at 21:10
1
$\begingroup$

Cross-validation is another approach to avoid over-optimism. I ran a logistic model with your data, using 4-fold cross-validation (i.e. fit the model on 75% of the cases, and then measure the error on the remaining 25% -- then do this 3 more times).

I get a cross-validated ROC area of .788, so your .81 looks pretty accurate. I think EdM's suggestion of using Frank Harrell's work is going to work better with the reviewers, though.

$\endgroup$
2
  • 1
    $\begingroup$ Thanks for adding in about cross validation as an extension of the single test/train breakdown approach. The rms validate() function does provide a cross-validation option to the default bootstrap. The cross validation should be extensive, though, to pass muster. Harrell notes 100 repeats of 10-fold cross-validation as an alternative to bootstrapping. $\endgroup$
    – EdM
    Aug 15 '19 at 20:30
  • $\begingroup$ Thanks a lot I used validate ( ) and it gives the same AUC value obtained in pRROC. $\endgroup$ Aug 15 '19 at 21:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.