# Significance testing of cross-validated classification accuracy: shuffling vs. binomial test

I have a dataset with 2 classes and a certain way to build a binary classifier. I want to measure its performance and to test if it is significantly above the chance level. I measure its performance with repeated cross-validation (see below). My question is: how to test the significance?

Here is my cross-validation procedure. I use 100-fold stratified Monte Carlo cross-validation (I am not exactly sure that's the correct term though; some people seem to call it boostrap, or out-of-bootstrap, or leave-group-out cross-validation): on each fold I randomly select $K=4$ test cases, 2 from each class, train the classifier on the remaining data, and classify these 4 cases. Let's say I get $a_i$ correct classifications. This is repeated $N=100$ times, and so I get an overall number of correct classifications $A = \Sigma a_i$. I report mean classification accuracy $A/400$ and the standard deviation of individual accuracies $a_i/4$.

• Note 1: For the reasons that I think are not very important here I cannot increase $K$ and cannot use the usual k-fold cross-validation, so this Monte Carlo approach is the only possible one for me. The variance of my estimator is quite large, but I have nothing else to do.
• Note 2: @FrankHarrel would say that classification accuracy is a "discontinuous improper scoring rule". I know, but in this particular case I am fine with it. I am not optimizing any model here, my classifier is already given.

Now, naively I would think that a random classifier would predict each case with probability 50%, so the number of correct classifications under null hypothesis of chance level classifier would be $\mathrm{Binom}(400, 0.5) \approx \mathcal{N}(200,100)$, so I can simply test if my $A$ is in the upper $\alpha$% (say 1%) percentile of this binomial/normal distribution.

However, I decided to do a shuffling test. I shuffle my labels, then use the whole above procedure with 100 folds to get a mean shuffled accuracy $B_j$, and repeat this shuffling $M=100$ times. The purpose is to sample accuracies under null hypothesis. After I obtain 100 values $B_j$ I look at their distribution. The mean is very close to 200, which is good. However, the variance is much larger than 100, it's around 1500. I don't understand how it is possible.

After I looked closer, I noticed that inside each shuffle the variance of correct classifications over 100 folds is around 1, as expected: $4*0.5*(1-0.5)=1$. But inside some shuffles the mean number of correct classifications is quite a bit below 2, and inside some other shuffles it is quite a bit over 2. This additional variation causes the variance of $B_j$ to be so high. In contrast, if inside each shuffle I use a truly random classifier instead of the classifier built on a training set according to my method, then I get $B_j$ nicely following $\mathcal{N}(200,100)$.

How is it possible? Should I use binomial or empirical distribution for statistical testing? The difference is very large. If I should use the empirical one, is there a way to somehow approximate it without actually performing the shuffles (which takes ages)?

Update

There is actually a third way to test for significance that comes to mind: I can look at the confidence interval of my $A$ and see if it excludes 200. No need to do any binomial testing, no need for any shuffling. However, I am confused about what variance I should build the confidence interval upon: is it the standard deviation of $a_i$ or the standard error of the mean of $a_i$? The difference is of course huge. If the total number of samples is large, then my individual Monte Carlo folds can be assumed to be independent, and I guess I can take standard error of the mean of $a_i$. But if the total number of samples is not so large, they are not independent anymore. On the other hand, standard deviation is too large to be useful for testing. Seems I would need to divide the variance of $a_i$ by the square root of effective number of folds (which would be less than 100), but I have no idea how to do estimate it.

• Mar 5 '14 at 18:01
• What seems to be a very relevant paper: Permutation Tests for Classification. Mar 5 '14 at 22:55
• I don't see how the confidence interval approach is that distinct from the binomial testing: if you believe your random variable to be binomially distributed you'd also construct a binomial c.i. - but that would throw you back to the same problem, woudn't it? Mar 6 '14 at 20:50
• @cbeleites: To test significance with the binomial confidence interval is of course equivalent to binomial testing. In my update I was referring to the confidence interval that one would try to construct from iterations of CV (even though I am not sure how to do it correctly in case of Monte Carlo CV that I described)... Anyway, please see the reply to this question that I just posted, I provided a very simple simulation there that I hope illustrates my observations. Mar 6 '14 at 23:52

Radmacher and colleagues (J. Comput. Biol. 9:505-511) describe a process for computing the significance of an error rate. We permute the class labels (a few thousand times) and repeat the entire cross-validation procedure to assess the probability of producing a cross-validated error rate as small as the observed one. One concept related to the current question is that the summaries from folds of cross-validation are not independent, and this presents challenges in obtaining confidence intervals and p-values. This is discussed by Jiang and colleagues (Stat Appl Genet Mol Biol. 2008;7(1)). Also note that a 'significant' error-rate measure is a bare-minimum requirement for a prediction rule and does not say much about its usefulness. However, given that small samples are often used and there are many problems with prediction rules, it is still helpful as a sanity check.

• Thank you! It took me some time to find the papers you cite, so I would like to leave the links here: A Paradigm for Class Prediction Using Gene Expression Profiles and Calculating Confidence Intervals for Prediction Error in Microarray Classification Using Resampling. I will read the papers and get back to you after that. Mar 5 '14 at 22:54
• I read your Radmacher et al. 2002 paper. I appreciate your Monte Carlo procedure, but I am also confused because it seems that you can achieve the same conclusions with a MUCH simpler binomial testing. You wrote: "The question is this: is a cross-validated error rate of 3 out of 14 small enough to ensure the class labels for this data set are not artifactual?" Well, 95% binomial CI given 11 correct classifications out of 14 is [0.49, 0.95] which does not exclude 0.5, so the direct answer is no. My question here is precisely this: is permutation test still preferred, and why? Mar 5 '14 at 23:16
• This relates to the dependence I was referring to. The binomial test assumes independence, which you have if you also have an independent validation set. But using cross-validation, the results from one fold are dependent because the training set shared observations in common. Mar 6 '14 at 0:27
• @julieth: whether the dependence of the training sets across the folds (and iterations) is a problem or wanted depends on what exactly you want to measure. If you want to estimate performance for unknown cases for the model built on the whole data set, this is desired behaviour. If you want to estimate performance of models built on $n$ cases of the underlying population, it is a problem. Mar 6 '14 at 13:44
• @cbeleites: For the current case, it was about why the usual binomial test cannot be used to assess significance and about putting an interval around the binomial proportion. Do you agree that in this case the dependence makes it so that the classical formula cannot be used? Thanks. Mar 6 '14 at 20:37

This is not really an answer to my question, but I would like to provide an explicit and very simple simulation demonstrating the symptoms I described, and I don't want to clutter the question too much.

Let's consider the absolutely simplest case possible. I will take 400 one-dimensional samples: 200 are equal to -1 and 200 are equal to +1. All samples located at -1 belong to class A, and all samples located at +1 belong to class B. The classifier will measure the mean (centroid) of each class and assign every test sample to the class whose centroid is closer. Nothing can be simpler than that. Here is an illustration ("A" means 200 points from class A, "B" means 200 points from class B):

---------A-----------B--------->

I can do exactly the same Monte Carlo cross-validation as I described above: I randomly select $K=4$ test cases, 2 from each class, train the classifier on the remaining training set, and classify these four; this is repeated 100 times. Obviously the number of correct classification is 400, i.e. 100% accuracy. Especially for @cbeleites I can also run a usual 100-fold (stratified) CV, the accuracy is also 100%. Note that it does not make sense to iterate this CV, because nothing will change.

And now we do the shuffles. I randomly shuffle the labels, repeat exactly the same procedure and get the number of correct classifications $B_j$. Then shuffle the labels again, and repeat it 100 times. Results: mean of $B_j$ is around 200 (very close, like 198-202 on different runs), so chance level. But variance of $B_j$ is in the range of 400-900 (on different runs). This is true for both Monte Carlo CV and standard 100-fold CV. The variance is always MUCH larger than the expected binomial variance that should be equal to 400*0.5*(1-0.5) = 100.

Now either I am overlooking a completely stupid mistake (which is absolutely possible!), or we have a big problem with all the reasoning from Confidence interval for cross-validated classification accuracy, because the binomial intervals don't make any sense. For example, if I spoil my ideal class separation by relabeling 80 points from class A to B and vice versa, then my actual number of decoded samples becomes 240. The stability over iterations of CV is perfect. The binomial confidence interval binofit(240,400) is [0.55, 0.65] which excludes 0.5 so we would conclude that the decoding is significant. But the variance of shuffled correct decoders is still on average around 500-600, so standard deviation is around let's say 22, so 95% interval for the null hypothesis of random decoding is around 200$\pm$45, which includes 240, which means not significant.

As far as I can see this problem has nothing to do with different CV folds not being independent, it's entirely different problem that has to do with finite sample size. The larger the sample size, the smaller the variance of $B_j$ (now I am back to Monte Carlo cross-validation, where I can still classify 400 cases even if the sample size is much larger). But I have to go to sample sizes above 10000 to get variance close to 100. Such sample sizes are way beyond realistic.

Update:

In comments above @julieth quoted a paper by Jiang et al. Calculating Confidence Intervals for Prediction Error in Microarray Classification Using Resampling: "... the test [I think they mean "training" -- amoeba] set on which prediction of the $i$-th case has $n-2$ specimens in common with the training set on which prediction of the $j$-th is based, hence, the number of prediction errors is not binomial". In other words, they claim that the reason for non-binomiality is that training sets are not mutually exclusive. Turns out, Nadeau and Bengio have a mammoth 49-page long paper about it called Inference for the Generalization Error where they discuss exactly this issue in great detail.

I did not believe it at first, so I used the simulation above to check this claim. If I increase the total number of samples to 4000, I can use the Monte Carlo CV procedure with 100 folds (each time classifying 4 test cases) to get 400 predictions. On the shuffled data (I increased the number of shuffles to 1000) the mean number of correct classification is 199 and the variance is 346: still a whole lot more than 100 even though the sample size is now as large as 4000.

But now I can also do the following: split my 4000 samples in 100 stratified parts of 40, and in each part use 36 samples to predict 4. I will also get 400 predictions, but this time all training sets are mutually exclusive. The outcome (also after 1000 shuffles): mean 199, variance 98.

Wow! Nadeau, Bengio, and @julieth seem to be right. And binomial assumption seems to be dead wrong. I wonder how many papers there are out there using binomial confidence intervals and tests...

• Sorry, I don't have time right now to think about this (nor at the weekend) but I think you need to post the actual code. It is also possible that for that question code review is the more suitable stackexchange. Mar 7 '14 at 0:34
• @cbeleites: So you think I have a mistake? I will double check and can also ask a colleague of mine to reproduce the simulation in Python. However, it does seem plausible to me. Binomial variance assumes that on each shuffling probability $p$ of correct decoding is 50%. But some shuffles will result in a bit better class separability, some in worse, so $p$ will fluctuate around $p=0.5$ and this will inflate the overall variance. Mar 7 '14 at 8:16
• Consider 5 observations. If we hold out "1", 2-5 make up the training set. If we hold out "2" 1, 3-5 make up the training set. So the training sets have a majority of observations in common and therefore predictions will not be independent. There are many ways to go wrong in cross-val, and many strong opinions out there about what is right and wrong. I commend your detailed study of a question you sought an answer for. Mar 8 '14 at 0:01
• "I wonder how many papers there are out there using binomial confidence intervals and tests" My applied impression is: far less than papers having no thought whatsoever on confidence intervals and reporting 4 digits of sensitivity on a basis of < 20 independent cases... Mar 11 '14 at 16:53
• @julieth: One thing I am still confused about though. Taking your example with 5 observations: if instead of classifying the held-out sample we choose to classify a random sample, then the variance will be strictly binomial. Even though trainings sets still have a majority of observations in common. So I am once again not sure that non-independence of the training sets is the correct (or full) explanation of over-dispersion of CV results under shuffling. Any thoughts about that? Mar 14 '14 at 15:43