24
$\begingroup$

I have a test set of 100 cases and two classifiers.

I generated predictions and computed ROC AUC, sensitivity and specificity for both classifiers.

Question 1: How can I compute p-value to check if one is significantly better than the other with respect to all scores (ROC AUC, sensitivity, specificity)?


Now, for the same test set of 100 cases, I have different and independent feature assignments for each case. This is because my features are fixed but subjective and provided by multiple (5) subjects.

So, I evaluated my two classifiers again for 5 "versions" of my test set and obtained 5 ROC AUCs, 5 sensitivities and 5 specificities for both classifiers. Then, I computed the mean of each performance measure for 5 subjects (mean ROC AUC, mean sensitivity and mean specificity) for both classifiers.

Question 2: How can I compute p-value to check if one is significantly better than the other with respect to mean scores (mean ROC AUC, mean sensitivity, mean specificity)?


Answers with some example python (preferably) or MatLab code are more than welcome.

$\endgroup$
6
  • $\begingroup$ Do direct comparison of precision, accuracy, AuC to get the best classifier among the two. P-value doesn't make sense here. The p-value is used in the context of evaluating if the model is doing better than random/50-50 assignments (as a null/alternate hypothesis test) $\endgroup$
    – wololo
    Commented Jul 20, 2018 at 5:48
  • 2
    $\begingroup$ First, I don’t agree that comparison of two performance measures using p-value doesn’t make sense here. I see that one classifier has AUC 0.80 and the other 0.85. My null hypothesis would be that there is no difference in performance of both classifieds. I want to know if the difference is statistically significant. $\endgroup$
    – kostek
    Commented Jul 20, 2018 at 14:03
  • 2
    $\begingroup$ Second, I don’t make 5 versions of my model. I have two models trained on a separate training set and now I evaluate them on 5 different “versions” of my test set. I have a mean performance for both classifiers (e.g. 0.81 AUC and 0.84 AUC) and want to check if the difference is statistically significant. $\endgroup$
    – kostek
    Commented Jul 20, 2018 at 14:07
  • 1
    $\begingroup$ I wouldn't say that what I'm doing is close to cross validation. In my case, the values of features depend on the subject that is providing them. I know that AUC can be used to compare models, but I want to know if, in my setting, the result of my comparison is statistically significant. I am sure that it can be done and that it makes a lot of sense to do it. My question is how to do it. $\endgroup$
    – kostek
    Commented Jul 21, 2018 at 4:27
  • 3
    $\begingroup$ I'm not sure what @Nishad is getting at, you can and should use a hypothesis test to determine if your models are significantly different from one another. The standard deviations of your metrics do exist, and get smaller as sample size increases (all other things being equal). An AUC difference between 0.8 and 0.9 may not be significant if you only have 10 samples, but might be very significant if you have 10M samples. I fail to see any relationship to cross-validation, as well. Would down-vote the comments if I could. $\endgroup$ Commented Jul 23, 2018 at 18:14

3 Answers 3

18
+100
$\begingroup$

Wojtek J. Krzanowski and David J. Hand ROC Curves for Continuous Data (2009) is a great reference for all things related to ROC curves. It collects together a number of results in what is a frustratingly broad literature base, which often uses different terminology to discuss the same topic.

Additionally, this book offers commentary and comparisons of alternative methods which have been derived to estimate the same quantities, and points out that some methods make assumptions which may be untenable in particular contexts. This is one such context; other answers report the Hanley & McNeil method, which assumes the binormal model for distributions of scores, which may be inappropriate in cases where the distribution of class scores are not (close to) normal. The assumption of normally distributed scores seems especially inappropriate in modern contexts, typical common models such as tend to produce scores with a "bathtub" distribution for classification tasks (that is, distributions with high densities in the extremes near 0 and 1).

Question 1 - AUC

Section 6.3 discusses comparisons of ROC AUC for two ROC curves (pp 113-114). In particular, my understanding is that these two models are correlated, so the information about how to compute $r$ is critically important here; otherwise, your test statistic will be biased because it doesn't account for the contribution of correlation.

For the case of uncorrelated ROC curves not based on any parametric distributional assumptions, statistics for tets and confidence intervals comparing AUCs can be straightforwardly based on estimates $\widehat{\text{AUC}}_1$ and $\widehat{\text{AUC}}_2$ of the AUC values, and estimates of their standard deviations $S_1$ and $S_2$, as given in section 3.5.1:

$$ Z = \frac{\widehat{\text{AUC}}_1 - \widehat{\text{AUC}}_2}{\sqrt{S_1^2 + S_2^2}} $$

To extend such tests to the case in which the same data is used for both classifiers, we need to take account of the correlation between the AUC estimates: $$ z=\frac{\widehat{\text{AUC}}_1 - \widehat{\text{AUC}}_2}{\sqrt{S_1^2 + S_2^2 - rS_1S_2}} $$

where $r$ is the estimate of this correlation. Hanley and McNeil (1983) made such an extension, basing their analysis on the binormal case, but only gave a table showing how to calculate the estimated correlation coefficient $r$ from the correlation $r_P$ of the two classifiers within class P, and the correlation of $r_n$ of the two classifiers within class N, saying that the mathematical derivation was available upon request. Various other authors (e.g. Zou, 2001) have developed tests based on the binormal model, assuming that an appropriate transformation can be found which will simultaneously transform the score distributions of classes P and N to normal.

DeLong et al (1988) took advantage of the identity between AUC and the Mann-Whitney test statistic, together with results from the theory of generalized $U$-statistics due to Sen (1960), to derive an estiamte of the correlation between the AUCs that does not rely on the binormal assumption. In fact, DeLong et al (1988) presented the following results for comparisons between $k\ge 2$ classifiers.

In Section 3.5.1, we showed that the area under the empirical ROC curve was equal to the Mann-Whitney $U$-statistic, and was given by

$$ \widehat{AUC}=\frac{1}{n_N n_P} \sum_{i=1}^{n_N} \sum_{j=1}^{n_P} \left[ I(s_{P_j} > s_{N_i}) + \frac{1}{2}I(s_{P_j} = s_{N_i}) \right] $$ where $s_{P_i}, i = 1, \dots,n_P$ are the score for the class $P$ objects and $s_{N_j}, j = 1, \dots,n_N$ are the scores for the class $N$ objects in the sample. Suppose that we have $k$ classifiers, yielding scores $s^r_{N_j}, j=1\dots n_N$ and $s_{P_i}^r, j = 1, \dots,n_P$ [I corrected an indexing error in this part - Sycorax], and $\widehat{AUC}_r, r = 1, \dots, k$. Define

$$ V^r_{10}=\frac{1}{n_N}\sum_{j=1}^{n_N} \left[ I(s_{P_i}^r > s_{N_j}^r) + \frac{1}{2}I(s_{P_i}^r = s_{N_j}^r) \right] , i=1,\dots,n_P $$ and $$ V^r_{01} = \frac{1}{n_P}\sum_{i=1}^{n_P} \left[ I(s_{P_i}^r > s_{N_j}^r) + \frac{1}{2}I(s_{P_i}^r = s_{N_j}^r) \right] , j=1,\dots,n_N $$

next, define the $k \times k$ matrix $\mathbf{W}_{10}$ with $(r,s)$th element $$ w_{10}^{r,s} = \frac{1}{n_P - 1}\sum_{i=1}^{n_P} \left[ V_{10}^r(s_{P_i}) - \widehat{AUC}_r \right] \left[ V_{10}^s(s_{P_i}) - \widehat{AUC}_s \right] $$ and the $k \times k$ matrix $\mathbf{W}_{01}$ with $(r,s)$th element $$ w_{01}^{r,s} = \frac{1}{n_N - 1}\sum_{i=1}^{n_N} \left[ V_{01}^r(s_{N_i}) - \widehat{AUC}_r \right] \left[ V_{01}^s(s_{N_i}) - \widehat{AUC}_s \right] $$ Then the estiamted covariance matrix for the vector $(\widehat{AUC}_1, \dots, \widehat{AUC}_k)$ of the estimated areas under the curves is $$ \mathbf{W} = \frac{1}{n_P}\mathbf{W}_{10} + \frac{1}{n_N}\mathbf{W}_{01} $$ with elements $w^{r,s}$. This is a generalization of the result for the estimated variance of a single estiamted AUC, also given in section 3.5.1. In the case of two classifiers, the estiamted correlation $r$ between the estimated AUCs is thus given by $\frac{w^{1,2}}{\sqrt{w^{1,1}w^{2,2}}}$ which can be used in $z$ above.

Since another answers gives the Hanley and McNeil expressions for estimators of AUC variance, here I'll reproduce the DeLong estimator from p. 68:

The alternative approach due to DeLong et al (1988) and exemplified by Pepe (2003) gives perhaps a simpler estimate, and one that introduces the extra useful concept of a placement value. The placement value of a score $s$ with reference to a specified population is that population's survivor function at $s$. This the placement value for $s$ in population N is $1 - F(s)$ and for $s$ in population P it is $1 - G(s)$. Empirical estimates of placement values are given by the obvious proportions. Thus the placement value of observation $s_{N_i}$ in population P denoted $s^P_{N_i}$, is the proportion of sample values from P that exceed $s_{N_i}$, and $\text{var}(s_{P_i}^N)$ is the variance of the placement values of each observation from N with respect to population P...

The DeLong et al (1988) estimate of variance of $\widehat{AUC}$ is given in terms of these variances: $$ s^2(\widehat{AUC}) = \frac{1}{n_P} \text{var}\left(s_{P_i}^N\right) + \frac{1}{n_N}\text{var}\left(s_{N_i}^P\right) $$

Note that $F$ is the cumulative distribution function of the scores in population N and $G$ is the cumulative distribution function of the scores in population P. A standard way to estimate $F$ and $G$ is to use the . The book also provides some alternative methods to the ecdf estimates, such as kernel density estimation, but that is outside the scope of this answer.

The statistics $Z$ and $z$ may be assumed to be standard normal deviates, and statistical tests of the null hypothesis proceed in the usual way. (See also: )

This is a simplified, high-level outline of how hypothesis testing works:

  • Testing, in your words, "whether one classifier is significantly better than the other" can be rephrased as testing the null hypothesis that the two models have statistically equal AUCs against the alternative hypothesis that the statistics are unequal.

  • This is a two-tailed test.

  • We reject the null hypothesis if the test statistic is in the critical region of the reference distribution, which is a standard normal distribution in this case.

  • The size of the critical region depends on the level $\alpha$ of the test. For a significance level of 95%, the test statistic falls in the critical region if $z > 1.96$ or $z < -1.96$. (These are the $\alpha/2$ and $1 - \alpha/2$ quantiles of the standard normal distribution.) Otherwise, you fail to reject the null hypothesis and the two models are statistically tied.

Question 1 - Sensitivity and Specificity

The general strategy for comparing sensitivity and specificity is to observe that both of these statistics amount to performing statistical inference on proportions, and this is a standard, well-studied problem. Specifically, sensitivity is the proportion of population P that has a score greater than some threshold $t$, and likewise for specificity wrt population N: $$ \begin{align} \text{sensitivity} = tp &= \mathbb{P}(s_P > t) \\ 1 - \text{specificity} = fp &= \mathbb{P}(s_N > t) \end{align} $$

The main sticking point is developing the appropriate test given that the two sample proportions will be correlated (as you've applied two models to the same test data). This is addressed on p. 111.

Turning to particular tests, several summary statistics reduce to proportions for each curve, so that standard methods for comparing proportions can be used. For example, the value of $tp$ for fixed $fp$ is a proportion, as is the misclassification rate for fixed threshold $t$. We can thus compare curves, using these measures, by means of standard tests to compare proportions. For example, in the unpaired case, we can use the test statistic $(tp_1 - tp_2) / s_{12}$, where $tp_i$ is the true positive rate for curve $i$ as the point in question, and $s_{12}^2$ is the sum of the variances of $tp_1$ and $tp_2$...

For the paired case, however, one can derive an adjustment that allows for the covariance between $tp_1$ and $tp_2$, but an alternative is to use McNemar's test for correlated proportions (Marascuilo and McSweeney, 1977).

The is appropriate when you have $N$ subjects, and each subject is tested twice, once for each of two dichotomous outcomes. Given the definitions of sensitivity and specificity, it should be obvious that this is exactly the test that we seek, since you've applied two models to the same test data and computed sensitivity and specificity at some threshold.

The McNemar test uses a different statistic, but a similar null and alternative hypothesis. For example, considering sensitivity, the null hypothesis is that the proportion $tp_1 = tp_2$, and the alternative is $tp_1 \neq tp_2$. Re-arranging the proportions to instead be raw counts, we can write a contingency table $$ \begin{array}{c|c|c|} & \text{Model 1 Positive at } t & \text{Model 1 Negative at } t \\ \hline \text{Model 2 Positive at } t & a & b \\ \hline \text{Model 2 Negative at } t & c & d \\ \hline \end{array} $$ where cell counts are given by counting the true positives and false negatives according to each model

$$ \begin{align} a &= \sum_{i=1}^{n_P} I(s_{P_i}^1 > t) \cdot I(s_{P_i}^2 > t) \\ b &= \sum_{i=1}^{n_P} I(s_{P_i}^1 \le t) \cdot I(s_{P_i}^2 > t) \\ c &= \sum_{i=1}^{n_P} I(s_{P_i}^1 > t) \cdot I(s_{P_i}^2 \le t) \\ d &= \sum_{i=1}^{n_P} I(s_{P_i}^1 \le t) \cdot I(s_{P_i}^2 \le t) \\ \end{align} $$

and we have the test statistic $$ M = \frac{(b-c)^2}{b + c} $$ which is distributed as $\chi^2_1$ a chi-squared distribution with 1 degree of freedom. With a level $\alpha=95\%$, the null hypothesis is rejected for $M > 3.841459$.

For the specificity, you can use the same procedure, except that you replace the $s^r_{P_i}$ with the $s^r_{N_j}$.

Question 2

It seems that it is sufficient to merge the results by averaging the prediction values for each respondent, so that for each model you have 1 vector of 100 averaged predicted values. Then compute the ROC AUC, sensitivty and specificity statistics as usual, as if the original models didn't exist. This reflects a modeling strategy that treats each of the 5 respondents' models as one of a "committee" of models, sort of like an ensemble.

$\endgroup$
6
  • $\begingroup$ Thanks for your answer and provided references. What about p-values for sensitivity and specificity? $\endgroup$
    – kostek
    Commented Jul 23, 2018 at 21:26
  • $\begingroup$ For Q1, does it mean that there is no difference between computing p-value for sensitivity and specificity and that they both always have the same p-value and I simply make a contingency table and run McNemar test on it? $\endgroup$
    – kostek
    Commented Jul 25, 2018 at 3:58
  • $\begingroup$ No, you’d do one test for each. $\endgroup$
    – Sycorax
    Commented Jul 25, 2018 at 4:19
  • $\begingroup$ That is a very detailed answer, thank you. About McNemar-test; what are exactly $a,b,c,d$? What proportions are these? $\endgroup$
    – Drey
    Commented Jul 25, 2018 at 5:30
  • $\begingroup$ @Drey They're not proportions; they're counts. I make this explicit in a revision. $\endgroup$
    – Sycorax
    Commented Jul 25, 2018 at 14:24
5
$\begingroup$

Let me keep the answer short, because this guide does explain a lot more and better.

Basically, you have your number of True Postives ($nTP$) and number of True Negatives ($nTN$). Also you have your AUC, A. The standard error of this A is:

$\texttt{SE}_A = \sqrt{\frac{A(1-A) + (nTP-1)(Q_1 - A^2)+(nTN-1)(Q_2 - A^2)}{nTP \cdot nTN}}$

with $Q_1 = A / (2 - A)$ and $Q_2 = 2A^2 / (1 + A)$.

To compare two AUCs you need to compute the SE of them both using:

$\texttt{SE}_{A_1 - A_2} = \sqrt{(SE_{A_1})^2 + (SE_{A_2})^2 - 2r\cdot (SE_{A_1})(SE_{A_2})}$

where $r$ is a quantity that represents the correlation induced between the two areas by the study of the same set of cases. If your cases are different, then $r=0$; otherwise you need to look it up (Table 1, page 3 in freely available article).

Given that you compute the $z$-Score by

$z = (A_1 - A_2) / SE_{A_1-A_2}$

From there you can compute p-value using probability density of a standard normal distribution. Or simply use this calculator.

This hopefully answers Question 1. - at least the part comparing AUCs. Sens/Spec is already covered by the ROC/AUC in some way. Otherwise, the answer I think lies in the Question 2.

As for Question 2, Central Limit Theorem tells us that your summary statistic would follow a normal distribution. Hence, I would think a simple t-test would suffice (5 measures of one classifier against 5 measures of the second classifier where measures could be AUC, sens, spec)

Edit: corrected formula for $\texttt{SE}$ ($\ldots - 2r \ldots$)

$\endgroup$
3
  • $\begingroup$ Thanks for provided links. For Question 1, If I set A to be sensitivity or specificity, would the equations for SE and z-Score hold? $\endgroup$
    – kostek
    Commented Jul 23, 2018 at 2:42
  • $\begingroup$ No, because sens only handles TPs and spec handles TNs. It is possible to compute confidence intervals for sens/spec with Binomial proportion CI, but be vigilant (small sample size?). Your $\hat{p}$ would be sens or spec. If CIs overlap in your comparison, then the difference would be not statistically significant under the alpha-level. $\endgroup$
    – Drey
    Commented Jul 23, 2018 at 4:27
  • $\begingroup$ Instead of nTP and nTN that you mentioned for the formula, one should use actual #[positive] and #[negative] (independently of what model predicts). $\endgroup$ Commented May 18, 2020 at 11:22
0
$\begingroup$

For Question 1, @Sycorax provided a comprehensive answer.

For Question 2, to the best of my knowledge, averaging predictions from subjects is incorrect. I decided to use bootstrapping to compute p-values and compare models.

In this case, the procedure is as follows:

For N iterations:
  sample 5 subjects with replacement
  sample 100 test cases with replacement
  compute mean performance of sampled subjects on sampled cases for model M1
  compute mean performance of sampled subjects on sampled cases for model M2
  take the difference of mean performance between M1 and M2
p-value equals to the proportion of differences smaller or equal than 0

This procedure performs one-tailed test and assumes that M1 mean performance > M2 mean performance.

A Python implementation of bootstrapping for computing p-values comparing multiple readers can be found in this GitHub repo: https://github.com/mateuszbuda/ml-stat-util

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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