Consiglierei il & articolo di McNeil del 1982 “ Il significato e luso dellarea sotto una caratteristica operativa del ricevitore (ROC ) curva .

Esempio

Hanno la seguente tabella dello stato della malattia e del risultato del test (corrispondente, ad esempio, al rischio stimato da un modello logistico). Il primo numero a destra è il numero di pazienti con stato di malattia vero “normale” e il secondo numero è il numero di pazienti con stato di malattia vero “anormale”:

(1) decisamente normale: 33/3
(2) probabilmente normale: 6/2
(3) discutibile: 6/2
(4) probabilmente anormale: 11/11
(5) Decisamente anormale: 2/33

Quindi ci sono in totale 58 pazienti “normali” e “51” anormali. Vediamo che quando il predittore è 1, “decisamente normale”, il paziente è solitamente normale (vero per 33 dei 36 pazienti), e quando è 5, “decisamente anormale” i pazienti sono solitamente anormali (vero per 33 dei 36 pazienti) 35 pazienti), quindi il predittore ha senso. Ma come dovremmo giudicare un paziente con un punteggio di 2, 3 o 4? Il valore limite impostato per giudicare un paziente come anormale o normale determina la sensibilità e la specificità del test risultante.

Sensibilità e specificità

Possiamo calcolare stimato sensibilità e specificità per differenti cutoff. (Dora in poi scriverò solo “sensibilità” e “specificità”, lasciando che la natura stimata dei valori sia implicita.)

Se scegliamo il nostro limite in modo da classificare tutto i pazienti come anormali, indipendentemente da ciò che dicono i risultati del test (cioè, scegliamo il cutoff 1+), otterremo una sensibilità di 51/51 = 1. La specificità sarà 0/58 = 0. Non suona così bene.

OK, quindi scegliamo un taglio meno rigoroso. Classifichiamo i pazienti come anormali solo se hanno un risultato del test di 2 o superiore. Quindi perdiamo 3 pazienti anormali e abbiamo una sensibilità di 48/51 = 0,94. Ma abbiamo una specificità molto maggiore, di 33/58 = 0,57.

Ora possiamo continuare, scegliendo vari valori limite (3, 4, 5,> 5). (Nellultimo caso, non classificheremo nessun paziente come anormale, anche se ha il punteggio di test più alto possibile di 5.)

La curva ROC

Se lo facciamo per tutti i possibili valori limite e tracciamo il grafico della sensibilità contro 1 meno la specificità, otteniamo la curva ROC. Possiamo utilizzare il seguente codice R:

 # Data norm = rep(1:5, times=c(33,6,6,11,2)) abnorm = rep(1:5, times=c(3,2,2,11,33)) testres = c(abnorm,norm) truestat = c(rep(1,length(abnorm)), rep(0,length(norm))) # Summary table (Table I in the paper) ( tab=as.matrix(table(truestat, testres)) )  

Loutput è:

  testres truestat 1 2 3 4 5 0 33 6 6 11 2 1 3 2 2 11 33  

Possiamo calcolare varie statistiche:

 ( tot=colSums(tab) ) # Number of patients w/ each test result ( truepos=unname(rev(cumsum(rev(tab[2,])))) ) # Number of true positives ( falsepos=unname(rev(cumsum(rev(tab[1,])))) ) # Number of false positives ( totpos=sum(tab[2,]) ) # The total number of positives (one number) ( totneg=sum(tab[1,]) ) # The total number of negatives (one number) (sens=truepos/totpos) # Sensitivity (fraction true positives) (omspec=falsepos/totneg) # 1 − specificity (false positives) sens=c(sens,0); omspec=c(omspec,0) # Numbers when we classify all as normal  

E usando questo, possiamo tracciare la curva ROC (stimata):

 plot(omspec, sens, type="b", xlim=c(0,1), ylim=c(0,1), lwd=2, xlab="1 − specificity", ylab="Sensitivity") # perhaps with xaxs="i" grid() abline(0,1, col="red", lty=2)  

Calcolo manuale della AUC

Possiamo calcolare molto facilmente larea sotto la curva ROC, usando la formula per larea di un trapezio:

 height = (sens[-1]+sens[-length(sens)])/2 width = -diff(omspec) # = diff(rev(omspec)) sum(height*width)  

Il risultato è 0,8931711.

Una misura di concordanza

LAUC può anche essere vista come una misura di concordanza.Se prendiamo tutte le coppie possibili di pazienti in cui uno è normale e laltro è anormale, possiamo calcolare con quale frequenza è quello anormale che ha il risultato del test più alto (più anormale) hanno lo stesso valore, lo consideriamo “mezza vittoria”):

 o = outer(abnorm, norm, "-") mean((o>0) + .5*(o==0))  

La risposta è ancora 0,8931711, larea sotto la curva ROC. Questo sarà sempre il caso.

Una visione grafica della concordanza

Come sottolineato da Harrell nella sua risposta, anche questa ha uninterpretazione grafica. Tracciamo il punteggio del test (stima del rischio) sullasse y e il vero stato della malattia sullasse x (qui con un po di tremolio, per mostrare punti sovrapposti):

 plot(jitter(truestat,.2), jitter(testres,.8), las=1, xlab="True disease status", ylab="Test score")  

Tracciamo ora una linea tra ogni punto a sinistra (un paziente “normale”) e ogni punto a destra (un paziente “anormale”). La proporzione di linee con una pendenza positiva (cioè la proporzione di coppie concordanti ) è lindice di concordanza (le linee piatte contano come “concordanza del 50%”).

È un po difficile visualizzare le linee effettive per questo esempio, a causa del numero di legami (punteggio di rischio uguale), ma con un po di jitter e trasparenza possiamo ottenere una trama ragionevole:

 d = cbind(x_norm=0, x_abnorm=1, expand.grid(y_norm=norm, y_abnorm=abnorm)) library(ggplot2) ggplot(d, aes(x=x_norm, xend=x_abnorm, y=y_norm, yend=y_abnorm)) + geom_segment(colour="#ff000006", position=position_jitter(width=0, height=.1)) + xlab("True disease status") + ylab("Test\nscore") + theme_light() + theme(axis.title.y=element_text(angle=0))  

Vediamo che la maggior parte delle linee inclina verso lalto, quindi lindice di concordanza sarà alto. Vediamo anche il contributo allindice di ogni tipo di coppia di osservazioni. La maggior parte proviene da pazienti normali con un punteggio di rischio di 1 in coppia con pazienti anormali con un punteggio di rischio di 5 (1-5 paia), ma parecchio proviene anche da 1–4 paia e 4-5 paia. Ed è molto facile calcolare lindice di concordanza effettivo in base alla definizione della pendenza:

 d = transform(d, slope=(y_norm-y_abnorm)/(x_norm-x_abnorm)) mean((d$slope > 0) + .5*(d$slope==0))  

La risposta è ancora 0,8931711, cioè lAUC.

Il test di Wilcoxon – Mann – Whitney

Esiste una stretta connessione tra la misura di concordanza e il Wilcoxon – Mann – Whitney test. In realtà, questultimo verifica se la probabilità di concordanza (cioè che sia il paziente anormale in una coppia casuale normale-anormale che avrà il risultato del test più “anormale”) è esattamente 0,5. E la sua statistica di test è solo una semplice trasformazione della probabilità di concordanza stimata:

 > ( wi = wilcox.test(abnorm,norm) ) Wilcoxon rank sum test with continuity correction data: abnorm and norm W = 2642, p-value = 1.944e-13 alternative hypothesis: true location shift is not equal to 0  

La statistica del test (W = 2642) conta il numero di coppie concordanti. Se lo dividiamo per il numero di coppie possibili, otteniamo un numero familiare:

 w = wi$statistic w/(length(abnorm)*length(norm))  

Sì, è 0,8931711, larea sotto la curva ROC.

Metodi più semplici per calcolare lAUC (in R)

Ma rendiamoci la vita più facile. Ci sono vari pacchetti che calcolano lAUC per noi automaticamente.

Il pacchetto Epi

Il pacchetto Epi crea una bella curva ROC con vari statistiche (inclusa lAUC) incorporate:

 library(Epi) ROC(testres, truestat) # also try adding plot="sp"  

Il pacchetto pROC

Mi piace anche il pacchetto pROC, dato che può levigare la stima ROC (e calcolare una stima AUC basata sul ROC livellato):

(La linea rossa è il ROC originale e la linea nera è il ROC smussato. Notare anche il rapporto di aspetto 1: 1 predefinito. Ha senso usarlo, poiché sia la sensibilità che la specificità hanno un valore 0-1 intervallo.)

LAUC stimato dal ROC livellato è 0,9107, simile, ma leggermente superiore, allAUC del ROC non livellato (se si osserva un nella figura, puoi facilmente vedere perché è più grande). (Anche se abbiamo davvero troppo pochi possibili valori distinti dei risultati del test per calcolare un AUC uniforme).

Il pacchetto rms

Pacchetto rms di Harrell può calcolare varie statistiche di concordanza correlate utilizzando la funzione rcorr.cens(). C Index nelloutput è lAUC:

 > library(rms) > rcorr.cens(testres,truestat)[1] C Index 0.8931711  

Il pacchetto caTools

Infine, abbiamo il pacchetto caTools e la sua funzione colAUC(). Ha alcuni vantaggi rispetto ad altri pacchetti (principalmente velocità e capacità di lavorare con dati multidimensionali – vedi ?colAUC) che a volte possono essere utili.Ma ovviamente fornisce la stessa risposta che abbiamo calcolato più e più volte:

 library(caTools) colAUC(testres, truestat, plotROC=TRUE) [,1] 0 vs. 1 0.8931711  

Parole finali

Molte persone sembrano pensare che lAUC ci dica quanto “buono” un test è. E alcune persone pensano che lAUC sia la probabilità che il test classifichi correttamente un paziente. non è . Come puoi vedere dallesempio e dai calcoli sopra, lAUC ci dice qualcosa su una famiglia di test, un test per ogni possibile cutoff.

E lAUC è calcolata in base a tagli che non si userebbero mai in pratica. Perché dovremmo preoccuparci della sensibilità e della specificità dei valori limite “privi di senso”? Tuttavia, questo è ciò su cui si basa (parzialmente) lAUC. (Ovviamente, se lAUC è molto vicino a 1, quasi tutti i test possibili avranno un grande potere discriminatorio e saremmo tutti molto felici.)

– linterpretazione anomala della coppia dellAUC è piacevole (e può essere estesa, ad esempio, ai modelli di sopravvivenza, dove vediamo se è la persona con il rischio (relativo) più alto che muore prima). Ma non lo si userebbe mai in pratica. È un raro caso in cui sa di avere una persona sana e una malata, non sa quale persona è malata e deve decidere quale di loro trattare. (In ogni caso, la decisione è facile; tratta quello con il rischio stimato più elevato.)

Quindi penso che studiare la curva ROC effettiva sarà più utile che semplicemente guardare la misura sommaria dellAUC. E se utilizzi il ROC insieme alle (stime dei) costi di falsi positivi e falsi negativi, insieme alle tariffe base di ciò che stai studiando, puoi ottenere da qualche parte.

Notare inoltre che lAUC misura solo la discriminazione , non la calibrazione. Cioè, misura se è possibile discriminare tra due persone (una malata e una sana), in base al punteggio di rischio. Per questo, considera solo i valori di rischio relativo (o classifica, se vuoi, vedi linterpretazione del test di Wilcoxon – Mann – Whitney), non quelli assoluti, che dovresti essere interessato. Ad esempio, se dividi per 2 ciascuna stima del rischio dal tuo modello logistico, otterrai esattamente la stessa AUC (e ROC).

Quando si valuta un modello di rischio, anche la calibrazione è molto importante. Per esaminarlo, esaminerai tutti i pazienti con un punteggio di rischio di circa 0,7, e vedrai se circa il 70% di questi era effettivamente malato. Fatelo per ogni possibile punteggio di rischio (possibilmente utilizzando una sorta di attenuazione / regressione locale). Traccia i risultati e otterrai una misura grafica della calibrazione .

Se hai un modello con entrambe una buona calibrazione e una buona discriminazione, allora iniziare ad avere un buon modello. 🙂

Commenti

• Grazie, @Karl Ove Hufthammer, questa è la risposta più completa che abbia mai ricevuto. Apprezzo particolarmente la tua sezione ” Final Words “. Lavoro eccellente! Grazie ancora!
• Grazie mille per questa risposta dettagliata. Sto lavorando con un set di dati in cui Epi :: ROC () v2.2.6 è convinto che lAUC sia 1,62 (no, non è uno studio mentalista), ma secondo il ROC, credo molto di più nello 0,56 che il codice sopra riportato risulta in.
• Penso che ci sia un piccolo errore in sens=c(sens,0); omspec=c(omspec,0), non dovrebbe ‘ essere sens=c(0, sens); omspec=c(0, omspec)? Viene stampato correttamente con il 0 iniziale ma non nel modo in cui è attualmente nella risposta.
• No, la definizione corrente è, AFAICS, corretta, @steveb, e si traduce in una trama corretta. Penso che ciò che forse confonde sia che la curva ROC è disegnata da destra a sinistra (cioè dallangolo in alto a destra allangolo in basso a sinistra), non da da sinistra a giusto , come la maggior parte dei grafici. Questo è solo il risultato di come ho definito le variabili; si potrebbe anche averlo tracciato da sinistra a destra (invertendo sia sens e omspec prima di tracciare).

Dai unocchiata a questa domanda: Capire la curva ROC

Ecco come costruire una curva ROC (da quella domanda):

Disegnare una curva ROC

dato un set di dati elaborato dal tuo classificatore di classificazione

• esempi di test di classificazione in base al punteggio decrescente
• inizia da $ (0, 0) $
• per ogni esempio $ x $ (nel ordine decrescente)
  • se $ x $ è positivo, sposta $ 1 / \ text {pos} $ su
  • se $ x $ è negativo, sposta $ 1 / \ text {neg} $ a destra

dove $ \ text {pos} $ e $ \ text {neg} $ sono rispettivamente le frazioni di esempi positivi e negativi.

Puoi usare questa idea per calcolare manualmente AUC ROC usando il seguente algoritmo:

auc = 0.0 height = 0.0 for each training example x_i, y_i if y_i = 1.0: height = height + tpr else auc = auc + height * fpr return auc 

Questa bella immagine animata con gif dovrebbe illustrare questo processo più chiaro

Commenti

• Grazie @Alexey Grigorev, questa è unottima grafica e probabilmente si rivelerà utile in futuro! +1
• Potrei spiegare qualcosa sulle ” frazioni di esempi positivi e negativi “, intendi il valore unitario minimo di due assi?
• @Allan Ruin: pos qui indica il numero di dati positivi. Supponiamo che tu abbia 20 punti dati, in cui 11 punti sono 1. Quindi, quando si disegna il grafico, abbiamo un rettangolo 11×9 (altezza x larghezza). Alexey Grigorev ha scalato ma lascia che sia ‘, se lo desideri. Ora, sposta 1 sul grafico ad ogni passaggio.

Il post di Karl ha molto di informazioni eccellenti, ma negli ultimi 20 anni non ho ancora visto un esempio di curva ROC che abbia cambiato il pensiero di qualcuno in una buona direzione. Lunico valore di una curva ROC a mio modesto parere è che la sua area è uguale a una probabilità di concordanza molto utile. La stessa curva ROC induce il lettore a usare cutoff, che è una cattiva pratica statistica.

Per quanto riguarda il calcolo manuale dellindice $ c $, tracciare un grafico con $ Y = 0,1 $ sul $ asse x $ e il predittore continuo o probabilità prevista che $ Y = 1 $ sullasse $ y $. Se colleghi ogni punto con $ Y = 0 $ con ogni punto con $ Y = 1 $, la proporzione delle linee che hanno una pendenza positiva è la probabilità di concordanza.

Qualsiasi misura che abbia un denominatore di $ n $ in questa impostazione sono regole di punteggio di accuratezza improprie e dovrebbero essere evitate. Ciò include proporzione classificata correttamente, sensibilità e specificità.

Per la funzione R Hmisc package rcorr.cens, stampare il lintero risultato per visualizzare più informazioni, in particolare un errore standard.

Commenti

• Grazie, @Frank Harell, apprezzo la tua prospettiva. Uso semplicemente la statistica c come probabilità di concordanza, poiché non ‘ mi piacciono i valori limite. Grazie ancora!

Ecco unalternativa al modo naturale di calcolare lAUC utilizzando semplicemente il trapezoidale regola per ottenere larea sotto la curva ROC.

LAUC è uguale alla probabilità che unosservazione positiva campionata casualmente abbia una probabilità prevista (di essere positiva) maggiore di unosservazione negativa campionata casualmente. Puoi usarlo per calcolare lAUC abbastanza facilmente in qualsiasi linguaggio di programmazione esaminando tutte le combinazioni a coppie di osservazioni positive e negative. Potresti anche campionare in modo casuale le osservazioni se la dimensione del campione fosse troppo grande. Se vuoi calcolare lAUC usando carta e penna, questo potrebbe non essere lapproccio migliore a meno che tu non abbia un campione molto piccolo / molto tempo. Ad esempio in R:

n <- 100L x1 <- rnorm(n, 2.0, 0.5) x2 <- rnorm(n, -1.0, 2) y <- rbinom(n, 1L, plogis(-0.4 + 0.5 * x1 + 0.1 * x2)) mod <- glm(y ~ x1 + x2, "binomial") probs <- predict(mod, type = "response") combinations <- expand.grid(positiveProbs = probs[y == 1L], negativeProbs = probs[y == 0L]) mean(combinations$positiveProbs > combinations$negativeProbs) [1] 0.628723 

Possiamo verificare utilizzando il pROC pacchetto:

library(pROC) auc(y, probs) Area under the curve: 0.6287 

Utilizzo del campionamento casuale:

mean(sample(probs[y == 1L], 100000L, TRUE) > sample(probs[y == 0L], 100000L, TRUE)) [1] 0.62896 

1. Hai un vero valore per le osservazioni.
2. Calcola la probabilità a posteriori e quindi classifica le osservazioni in base a questa probabilità.
3. Ipotizzando una probabilità di interruzione di $ P $ e un numero di osservazioni $ N $:
   $$ \ frac {\ text {Sum of true ranks} -0.5PN (PN + 1)} { PN (N-PN)} $$

Commenti

• @ user73455 … 1) Sì, ho il valore vero per le osservazioni. 2) Probabilità a posteriori è sinonimo di probabilità predette per ciascuna delle osservazioni? 3) Capito; tuttavia, cosè la ” Somma dei veri ranghi ” e come si calcola questo valore? Forse un esempio potrebbe aiutarti a spiegare questa risposta in modo più approfondito? Grazie!