Test parametrici e non parametrici per dati appaiati

Il confronto tra dati appaiati si applica quando la stessa variabile viene misurata nello stesso caso [1] in due occasioni diverse. In campo medico la situazione tipica è quella di una misura che viene effettuata sullo stesso soggetto prima e dopo uno specifico trattamento.

Il Volume Espiratorio Massimo nel primo Secondo (VEMS), impiegato nella diagnostica della capacità respiratoria, è il volume di aria espirata nel corso del primo secondo di una espirazione massima forzata e indica il grado di pervietà delle grandi vie aeree. Viene anche denominato FEV1 dall'acronimo inglese di Forced Expiratory Volume in the 1st second.

I seguenti dati, ricavati da Campbell [2], riportano i valori di VEMS (FEV1) misurati in 5 soggetti asmatici prima (al tempo t0) e dopo (al tempo t1) l'assunzione di un broncodilatatore, e la loro differenza (t0 - t1).

t0	t1	Differenza
1.5	1.7	-0.2
1.7	1.9	-0.2
2.1	2.2	-0.1
1.6	1.9	-0.3
2.4	2.4	0

Per proseguire ora è necessario:

→ effettuare il download del file FEV1.csv 

→ salvare il file nella cartella C:\Rdati\

Per il file di dati vedere istruzioni e link riportati alla pagina Dati, ma potete anche semplicemente copiare i dati riportati qui sotto aggiungendo un ↵ Invio al termine dell'ultima riga e salvarli in C:\Rdati\ in un file di testo denominato FEV1.csv (assicuratevi che il file sia effettivamente salvato con l'estensione .csv). 

t0;t1;Differenza

1.5;1.7;-0.2

1.7;1.9;-0.2

2.1;2.2;-0.1

1.6;1.9;-0.3

2.4;2.4;0

Copiate e incollate nella Console di R questo script e premete ↵ Invio:

# CONFRONTO TRA DUE CAMPIONI PER DATI APPAIATI

#

mydata <- read.table("c:/Rdati/FEV1.csv", header=TRUE, sep=";") # importa i dati

attach(mydata) # consente di impiegare direttamente i nomi delle variabili

mydata # mostra i dati

#

library(nortest) # carica il pacchetto

lillie.test(Differenza) # test di normalità di Lilliefors (Kolmogorov-Smironv) sulla Differenza (t0 - t1)

#

t.test(t0, t1, paired=TRUE) # test t di Student per dati appaiati (test parametrico)

#

wilcox.test(t0, t1, paired=TRUE, exact=FALSE) # test di Wilcoxon per dati appaiati (test non parametrico)

#

detach(mydata) # termina l'impiego diretto dei nomi delle variabili

#

Dopo avere importato i dati del file nell'oggetto mydata e dopo avere eseguito la funzione attach() che consente nel codice successivo di impiegare direttamente i nomi delle variabili, con mydata sono mostrati i dati importati:

> mydata # mostra i dati 

     t0   t1 Differenza

1   1.5  1.7       -0.2

2   1.7  1.9       -0.2

3   2.1  2.2       -0.1

4   1.6  1.9       -0.3

5   2.4  2.4        0.0

La scelta tra il test parametrico, il test t di Student, e il test non parametrico, il test di Wilcoxon (o Wilcoxon Signed Rank Test) viene effettuata come al solito sulla base dei risultati dei test di normalità (gaussianità).

Qui viene impiegato come test di normalità (gaussianità) il test di Lilliefors (Kolmogorov-Smirnov) che viene applicato alla variabile Differenza (t0 - t1) presa con il segno:

> lillie.test(mydata$Differenza) # test di normalità di Lilliefors (Kolmogorov-Smironv) sulla Differenza (t0 - t1)

        Lilliefors (Kolmogorov-Smirnov) normality test

data:  Differenza

D = 0.23714, p-value = 0.4686

Il test conferma una distribuzione gaussiana dei dati (p = 0.4686). Si procede quindi a impiegare le conclusioni tratte dal test t di Student

> t.test(t0, t1, paired=TRUE) # test t di Student per dati appaiati

        Paired t-test

data:  t0 and t1

t = -3.1379, df = 4, p-value = 0.03492

alternative hypothesis: true difference in means is not equal to 0

95 percent confidence interval:

 -0.30157148 -0.01842852 

sample estimates:

mean of the differences 

                  -0.16 

che con un valore di p = 0.03492 indica significatività (anche se abbastanza risicata, di poco inferiore a 0.05) della differenza tra VEMS misurato prima e VEMS misurato dopo la somministrazione del farmaco.

Soggiacente al test t di Student è l'ipotesi che la media delle differenze t0 – t1 sia zero (ipotesi H₀ o ipotesi nulla). Ma il test eseguito con R verifica anche l'ipotesi alternativa, cioè che la media delle differenze non sia uguale a 0 (zero). 

L'intervallo di confidenza al 95% calcolato per questa ipotesi va da -0.30157148 a -0.01842852 e indica che la media delle differenze, uguale a -0.16, è significativamente diversa da 0 (non sarebbe significativamente diversa da 0 se l'intervallo di confidenza della media delle differenze includesse lo 0). Le due soluzioni si confermano l'una con l'altra.

Lo script prevede comunque di eseguire anche il test di Wilcoxon, questo è il risultato:

> wilcox.test(t0, t1, paired=TRUE, exact=FALSE) # test di Wilcoxon  per dati appaiati (test non parametrico)

        Wilcoxon signed rank test with continuity correction

data:  t0 and t1

V = 0, p-value = 0.09751

alternative hypothesis: true location shift is not equal to 0

Dato il ridotto numero di casi l'ipotesi alternativa non viene sviluppata. 

Da notare che assumendo come soglia di significatività il classico p = 0.05 il test di Wilcoxon con un valore di p = 0.09751 indica una differenza non significativa. Una cosa che non deve meravigliare, in quanto è noto che le statistiche non parametriche aumentano il valore di p ovvero rendono la differenza meno significativa rispetto a quella ottenuta con le statistiche parametriche (che hanno evidenziato una significatività risicata).

Il fatto che approcci statistici differenti possono portare a conclusioni opposte, e il fatto che una significatività statistica (risicata) non implica necessariamente una rilevanza pratica, meritano sempre attenzione nella fase di valutazione critica dei risultati di un'analisi statistica.

Una interessante rappresentazione grafica dei dati appaiati qui riportati può essere effettuata mediante un grafico prima-dopo (slopegraph) al quale si rimanda.

----------

[1] Un caso è per definizione un soggetto/oggetto univocamente identificato, come illustrato anche nel post Importazione dei dati da un file .csv.

[2] Campbell MJ, Machin D. Medical Statistics. A Commonsense Approach. John Wiley & Sons, New York, 1993, ISBN 0-471-93764-9, p. 142.

Test parametrici e non parametrici per due campioni indipendenti

Attenzione: nel caso di confronti multipli, cioè se dovete effettuare il confronto tra più campioni indipendenti come ad esempio, nel caso più semplice di tre campioni, confrontare il primo campione con il secondo campione, il primo campione con il terzo e il secondo campione con il terzo, dovete impiegare uno dei test riportati nel post Test parametrici e non parametrici per più campioni indipendenti. 

Qui ci occupiamo invece del confronto tra due campioni indipendenti [1] che viene effettuato verificando se le medie (la tradizionale misura di posizione, parametrica) o le mediane (la misura di posizione non parametrica alternativa alla media) dei due campioni differiscono significativamente tra loro. 

Come nel caso delle statistiche elementari anche in questo caso è necessario effettuare una analisi esplorativa dei dati e i test di normalità (gaussianità) [2] per decidere quale sia il test appropriato da impiegare:

→ il tradizionale test parametrico per il confronto tra medie, rappresentata dal test t di Student, se i dati sono distribuiti in modo gaussiano;

→ uno dei test non parametrici per il confronto tra mediane, che devono essere impiegati quando i dati non sono distribuiti in modo gaussiano, come il test di Wilcoxon per campioni indipendenti (in genere meglio noto come test U di Mann-Whitney) e il test di Kruskal-Wallis. 

Lo script è diviso in cinque parti al fine di sviluppare in modo analitico il percorso logico da seguire, che prevede:

→ verifica della distribuzione dei dati mediante i test di normalità (gaussianità);

→ se la distribuzione è gaussiana in entrambi i campioni messi a confronto, esecuzione del test F di Fisher per verificare se le varianze dei due campioni sono omogenee;

→ se la distribuzione è gaussiana e con il test di Fisher le varianze sono omogenee (non differiscono significativamente) esecuzione dei test t di Student per varianze omogenee;

→ se la distribuzione è gaussiana e con il test di Fisher le varianze non sono omogenee (differiscono significativamente) esecuzione dei test t di Student per varianze non omogenee;

→ se la distribuzione non è gaussiana, esecuzione dei test non parametrici (test di Wilcoxon per campioni indipendenti e/o test di Kruskall-Wallis).

Vediamo un esempio con i dati rilevati in atleti di sesso femminile e di sesso maschile e la domanda: altezza, peso, percentuale di grasso corporeo, emoglobina e gli altri dati rilevati, differiscono significativamente nei due sessi, o sono simili tra loro? I dati sono contenuti nella tabella ais del pacchetto DAAG, accertatevi di avere installato il pacchetto o in alternativa procedete come indicato in [3].

Copiate questa prima parte dello script, incollatela nella Console di R e premete ↵ Invio:

# CONFRONTO TRA DUE CAMPIONI INDIPENDENTI (1/5)

#

library(DAAG) # carica il pacchetto DAAG incluso il set di dati ais

str(ais) # mostra la struttura di ais

attach(ais) # consente di impiegare direttamente i nomi delle variabili

#

# test di normalità

#

library(nortest) # carica il pacchetto

lillie.test(wt[sex == "f"]) # test di Lilliefors (Kolmogorov-Smirnov) su peso corporeo donne

lillie.test(wt[sex == "m"]) # idem su peso corporeo uomini

lillie.test(pcBfat[sex == "f"]) # test di Lilliefors (Kolmogorov-Smirnov) su percentuale di grasso corporeo donne

lillie.test(pcBfat[sex == "m"]) # idem su percentuale di grasso corporeo uomini

#

Dopo avere caricato i dati con il pacchetto DAAG, avere mostrato la struttura di ais e avere consentito con attach(ais) di impiegare direttamente nel codice che segue i nomi delle variabili del set di dati ais, viene caricato il pacchetto nortest che consente di effettuare i test di normalità.

Scegliamo uno dei più classici e collaudati, il test di Lilliefors (Kolmogorov-Smirnov) e lo applichiamo al peso corporeo (variabile wt) e alla percentuale di grasso corporeo (variabile pcBfat) di uomini (m) e donne (f).

Nel caso del peso corporeo wt il test di Lilliefors conferma, sia per donne sia per uomini, che i dati non si discostano significativamente da una distribuzione gaussiana (p = 0.6532 e p = 0.6384 rispettivamente).

> lillie.test(wt[sex == "f"]) # test di Lilliefors (Kolmogorov-Smirnov) su peso corporeo donne

        Lilliefors (Kolmogorov-Smirnov) normality test

data:  wt[sex == "f"]

D = 0.05469, p-value = 0.6532

> lillie.test(wt[sex == "m"]) # idem su peso corporeo uomini

        Lilliefors (Kolmogorov-Smirnov) normality test

data:  wt[sex == "m"]

D = 0.054681, p-value = 0.6384

Nel caso della percentuale di grasso corporeo pcBfat invece ci troviamo di fronte a dati che se nelle donne non si discostano significativamente da una distribuzione gaussiana (p = 0.3291), negli uomini invece non sono distribuiti in modo gaussiano (p = 2.254e-08).

> lillie.test(pcBfat[sex == "f"]) # test di Lilliefors (Kolmogorov-Smirnov) su percentuale di grasso corporeo donne

        Lilliefors (Kolmogorov-Smirnov) normality test

data:  pcBfat[sex == "f"]

D = 0.066987, p-value = 0.3291

> lillie.test(pcBfat[sex == "m"]) # idem su percentuale di grasso corporeo uomini

        Lilliefors (Kolmogorov-Smirnov) normality test

data:  pcBfat[sex == "m"]

D = 0.17702, p-value = 2.254e-08

Per poter confrontare due campioni impiegando un test parametrico, è necessario che i dati siano distribuiti in modo gaussiano in entrambi, pertanto:

→ per il confronto del peso corporeo wt in uomini e donne possiamo impiegare un test parametrico, il test t di Student;

→ per il confronto della percentuale di grasso corporeo pcBfat in uomini e donne dobbiamo impiegare un test non parametrico (in realtà ne vediamo due, il test di Wilcoxon per campioni indipendenti e il test di Kruskal-Wallis).

Al peso corporeo, distribuito in modo gaussiano, applichiamo innanzitutto mediante la funzione var.test() il test F di Fisher per il controllo dell'omogeneità tra le varianze dei due campioni da mettere a confronto, che sono ricavati aggregando i valori della variabile peso corporeo (wt) in base al valore della variabile sesso (sex) mediante l'argomento wt~sex:

# test parametrico (confronto tra medie) per il peso corporeo, verifica preliminare (2/5)

#

var.test(wt~sex) # controllo dell'omogeneità delle varianze con il test F di Fisher

#

Il test F di Fisher indica che le varianze osservate in uomini e donne non sono significativamente diverse (la probabilità di osservare per caso una differenza tra le medie pari a quella effettivamente osservata è p = 0.2029) e anche l'ipotesi alternativa conferma questo:

> var.test(wt~sex) # controllo dell'omogeneità delle varianze con il test F di Fisher

        F test to compare two variances

data:  wt by sex

F = 0.77423, num df = 99, denom df = 101, p-value = 0.2029

alternative hypothesis: true ratio of variances is not equal to 1

95 percent confidence interval:

 0.5221795 1.1488587

sample estimates:

ratio of variances 

         0.7742343 

Qualora le varianze dei due campioni fossero state significativamente diverse, sarebbe stato necessario procedere con il test t di Student per varianze non omogenee, impostando nella funzione t.test() l'argomento var.equal=FALSE.

Nel nostro caso invece procediamo con il test t di Student per varianze omogenee impostando nella funzione t.test() l'argomento var.equal=TRUE e aggregando di nuovo i valori della variabile peso corporeo (wt) in base ai valori della variabile sesso (sex) mediante l'argomento wt~sex.

# test parametrico (confronto tra medie) per il peso corporeo (3/5)

#

t.test(wt~sex, var.equal=TRUE) # test t di Student per varianze omogenee

#

Questi sono i risultati:

> t.test(wt~sex, var.equal=TRUE) # test t di Student per varianze omogenee

        Two Sample t-test

data:  wt by sex

t = -9.2272, df = 200, p-value < 2.2e-16

alternative hypothesis: true difference in means is not equal to 0

95 percent confidence interval:

 -18.42589 -11.93717

sample estimates:

mean in group f mean in group m 

       67.34200        82.52353 

La media del peso corporeo nelle donne è 67.342, negli uomini è 82.52353 e le due medie sono significativamente diverse, in quanto la probabilità di osservare per caso una tale differenza è inferiore a 0.00000000000000022 (p < 2.2e-16) [1].

Soggiacente al test t di Student è l'ipotesi che le due medie siano uguali ovvero che la differenza tra le due medie sia zero (ipotesi H₀ o ipotesi nulla). Ma il test eseguito con R verifica anche l'ipotesi alternativa, cioè che la differenza tra le due medie non sia uguale a 0 (zero). L'intervallo di confidenza al 95% della differenza tra medie va da -18.42589 a -11.93717 e non include lo 0 indicando che la differenza tra le due medie, pari a -15.18153 (67.342 - 82.52353) è significativamente diversa da 0 (la differenza non sarebbe significativamente diversa se l'intervallo di confidenza della differenza tra le due medie includesse lo 0). Le due soluzioni si confermano l'una con l'altra e ci consentono di affermare che le donne pesano in media significativamente meno degli uomini.

Alla percentuale di grasso corporeo, che non è distribuita in modo gaussiano, applichiamo i test non parametrici per il confronto tra mediane con le funzioni wilcox.test() e kruskal.test(), in ciascuna delle quali i due campioni da mettere a confronto sono ottenuti aggregando i valori della variabile percentuale di grasso corporeo (pcBfat) in base ai valori della variabile sesso (sex) mediante l'argomento pcBfat~sex:

# test non parametrici (confronto tra mediane) per la percentuale di grasso corporeo (4/5)

#

wilcox.test(pcBfat~sex) # test di Wilcoxon per campioni indipendenti

#

kruskal.test(pcBfat~sex) # test di Kruskal-Wallis

#

median(pcBfat[sex == "f"]) # mediana della percentuale di grasso corporeo per sesso = f

median(pcBfat[sex == "m"]) # mediana della percentuale di grasso corporeo per sesso = m

#

Diversamente dalla funzione t.test(), che riporta al termine dei calcoli il valore della media dei due campioni, le funzioni wilcox.test() e kruskal.test() non forniscono la mediana, per cui allo script sono state aggiunte due righe di codice per calcolarla sia nelle donne sia negli uomini.

La mediana della percentuale di grasso corporeo nelle donne è 17.94, negli uomini è 8.625 ed entrambi i test confermano che le due mediane sono significativamente diverse, la probabilità di osservare per caso una tale differenza è per entrambi i test inferiore a 0.00000000000000022 (p < 2.2e-16) [4].

> wilcox.test(pcBfat~sex) # test di Wilcoxon per campioni indipendenti

        Wilcoxon rank sum test with continuity correction

data:  pcBfat by sex

W = 9417.5, p-value < 2.2e-16

alternative hypothesis: true location shift is not equal to 0

> kruskal.test(pcBfat~sex) # test di Kruskal-Wallis

        Kruskal-Wallis rank sum test

data:  pcBfat by sex

Kruskal-Wallis chi-squared = 108.03, df = 1, p-value < 2.2e-16

> median(pcBfat[sex == "f"]) # mediana della percentuale di grasso corporeo per sesso = f

[1] 17.94

> median(pcBfat[sex == "m"]) # mediana della percentuale di grasso corporeo per sesso = m

[1] 8.625

Quindi possiamo affermare che nelle donne la mediana della percentuale di grasso corporeo è significativamente maggiore di quella rilevata uomini.

Un interessante integrazione ai risultati numerici ci viene ancora una volta dalla grafica, e in particolare dai grafici a scatola con i baffi:

# analisi grafica non parametrica delle distribuzioni (5/5)

#

windows() # apre una nuova finestra

par(mfrow=c(1,2)) # suddivide la finestra in due quadranti, uno per grafico

#

boxplot(wt~sex, data=ais, xlab="Sesso", ylab="Peso corporeo in kg", notch=TRUE, col="green") # traccia i boxplot per sesso del peso corporeo

#

boxplot(pcBfat~sex, data=ais, xlab="Sesso", ylab="Grasso corporeo in %", notch=TRUE, col="green") # traccia i boxplot per sesso della percentuale di grasso corporeo

#

Questo è il grafico

nel quale i boxplot per sesso del peso corporeo e della percentuale di grasso corporeo sono tracciati con una incisura (notch=TRUE) che rappresenta i limiti di confidenza al 95% della mediana. Se i limiti di confidenza delle mediane confrontate non si sovrappongono, le mediane sono significativamente diverse.

Qui il fatto interessante è che andiamo ben oltre l'impiego della grafica come complemento dei risultati dell'analisi numerica: abbiamo un esempio di come sia possibile realizzare un test statistico (non parametrico) mediante una rappresentazione grafica dei dati. La conclusione è che sia per il peso corporeo sia per la percentuale di grasso corporeo le incisure dei boxplot di donne e uomini non si sovrappongono, e pertanto le corrispondenti mediane risultano significativamente diverse in entrambi i casi.

----------

[1] Se i due campioni non sono indipendenti impiegare uno dei test riportati nel post Test parametrici e non parametrici per dati appaiati.

[2] Vedere il post Analisi esplorativa dei dati e i post a questo collegati.

[3] Vedere il post Il set di dati ais nel quale trovate anche come caricare i dati della tabella senza impiegare il pacchetto DAAG. 

[4] Viene riportato questo valore quando R non è più in grado di approssimare numericamente il risultato.