venerdì 7 agosto 2009

Test di analisi di omogeneità della varianza

Una delle applicazioni più classiche di questi test, è l'analisi della varianza con il metodo ANOVA. Per poter procedere con una ANOVA omoschedastica, è necessario fare l'assunzione che le varianze dei k gruppi in analisi siano omogenee (altrimenti occorre procedere con una ANOVA eteroschedastica). I test per verificare l'omogeneità delle varianze sono numerosi; alcuni vengono applicati su campioni che seguono una distribuzione normale, altri vengono utilizzati per l'analisi non parametrica (campioni non gaussiani). Nel post sull'ANOVA parametrica omoschedastica, abbiamo già visto 2 test di verifica dell'omogeneità delle varianze. Vediamoli nuovamente in questo post assieme ad altri test (sono tutti test parametrici, che pertanto possono essere applicati se i dati seguono una distribuzione normale; la verifica di normalità è discussa in questo post):

  • Test di Bartlett

  • Test di Fligner-Killeen

  • Brown e Forsyth test

  • Test di Hartley

  • Metodo di Cochran

  • Test di Levene parametrico


Un test per la verifica dell'omogeneità delle varianze, non parametrico (quando i gruppi non seguono una distribuzione normale), è il test di Levene non parametrico, che sfrutta il test di Kruskal-Wallis (ne parlo in fondo al post).

Innanzitutto vediamo come presentare i dati. In generale si deve creare un vettore contenente tutti i dati, e una seconda variabile di tipo factor, per mantenere la suddivisione in gruppi dei dati.
Consideriamo ad esempio i seguenti 4 gruppi, costituiti ciascuno da 5 elementi, e creiamo i vettori dati e gruppi:


a <- c(191, 211, 204, 209, 207)
b <- c(137, 148, 127, 137, 153)
c <- c(172, 165, 186, 178, 189)
d <- c(199, 206, 233, 183, 194)

dati <- c(a,b,c,d)
gruppi = factor(rep(letters[1:4], each = 5))


Se tutti i gruppi hanno la stessa numerosità campionaria, è sufficiente modificare il parametro each, dichiarando la lunghezza dei gruppi. Se i gruppi sono invece di numerosità campionarie differenti, la procedura è leggermente diversa: ne parlerò più giù, a proposito del test di Levene.

Può essere utile creare un dataframe (alcuni test accettano anche questa classe di dati):


data <- data.frame(dati, gruppi)
data

dati gruppi
1 191 a
2 211 a
3 204 a
4 209 a
5 207 a
6 137 b
7 148 b
8 127 b
9 137 b
10 153 b
11 172 c
12 165 c
13 186 c
14 178 c
15 189 c
16 199 d
17 206 d
18 233 d
19 183 d
20 194 d





Test di Bartlett

Questo è sicuramente uno dei più utilizzati, specie in biostatica. E' disponibile nel pacchetto standard (stats) di R, e quindi può essere richiamato velocemente:


# potremmo anche digitare: bartlett.test(dati ~ gruppi, data=data)
> bartlett.test(dati, gruppi)

Bartlett test of homogeneity of variances

data: dati and gruppi
Bartlett's K-squared = 3.3243, df = 3, p-value = 0.3443


In questo caso, otteniamo p-value > 0.05. Quindi accettiamo (o per meglio dire, non rifiutiamo) l'ipotesi nulla: le varianze sono omogenee.




Test di Fligner-Killeen

Anche questo test è disponibile nel pacchetto stats. Riceve i dati con le stesse modalità della funzione per il test di Bartlett:


fligner.test(dati ~ gruppi, data=data)

Fligner-Killeen test of homogeneity of variances

data: dati by gruppi
Fligner-Killeen:med chi-squared = 1.3866, df = 3, p-value = 0.7087


Anche qui abbiamo ottenuto p-value > 0.05, quindi le varianze risultano essere omogenee.




Test di Brown e Forsyth

Per poter utilizzare questo test occorre scaricare, installare e quindi aprire la libreria HH. Quindi si procede come seguen:


library(HH)
hov(dati ~ gruppi, data=data)

hov: Brown-Forsyth

data: dati
F = 0.6509, df:gruppi = 3, df:Residuals = 16, p-value = 0.5939
alternative hypothesis: variances are not identical


I dati possono essere presentati solamente sottoforma di una formula (come quelle da fittare per una regressione). Il p-value risulta maggiore di 0.05, quindi le varianze sono omogenee.




Test di Hartley

Non ho trovato questo test implementato in nessun pacchetto. Si procede in questo modo: si calcolano le varianze di tutti i gruppi e quindi si calcola il rapporto (varianza massima) / (varianza minima):


varianze <- c(var(a), var(b), var(c), var(d))
varianze
[1] 62.8 104.8 97.5 351.5

max(varianze) / min(varianze)
[1] 5.597134


Il valore calcolato va confrontato con le tavole tabulate della distribuzione F-max (che non è la ben più nota distribuzione F di Fisher-Snedocor), essendo noti il numero di gruppi k e i gradi di libertà di ciascun gruppo con stessa numerosità campionaria (n-1).
Potete trovare una tavola tabulata a questo link.

Andiamo a cercare in orizzontale il numero k di varianze a confronto (4) e in verticale di gradi di libertà (4), e osserviamo che il valore tabulato è 20.6 per alpha = 0.05.
Poichè il valore calcolato è inferiore al valore tabulato, accettiamo l'ipotesi nulla: le varianze sono omogenee.




Metodo di Cochran

Questo test è disponibile nel pacchetto outliers. Dopo averlo scaricato e installato, possiamo applicare la funzione cochran.test:


library(outliers)
cochran.test(dati ~ gruppi)

Cochran test for outlying variance

data: dati ~ gruppi
C = 0.5701, df = 5, k = 4, p-value = 0.1117
alternative hypothesis: Group d has outlying variance
sample estimates:
a b c d
62.8 104.8 97.5 351.5


Il p-value è maggiore di 0.05, quindi anche secondo questo test le varianze sono omogenee. Se non vogliamo scaricare il pacchetto, possiamo facilmente calcolare il valore della statistica test, da confrontare con le tavole tabulate della distribuzione per le varianze di Cochran.
Se sospettiamo che è presente una varianza più grande delle altre, si calcola il rapporto (varianza massima) / (somma di tutte le varianze). In alternativa si calcola il rapporto (varianza minima) / (somma di tutte le varianze). In questo caso è presente una varianza più grande delle altre, quindi calcoliamo il rapporto massimo:


varianze <- c(var(a), var(b), var(c), var(d))
varianze
[1] 62.8 104.8 97.5 351.5

max(varianze) / sum(varianze)
[1] 0.5700616


Questo valore va confrontato con le tavole di Cochran, disponibili in questo pdf.
Il valore tabulato per 4 gruppi, e 4 gradi di libertà, con alpha = 0.05, è 0.6287, che è maggiore del valore calcolato. Quindi accettiamo l'ipotesi nulla: le varianze sono omogenee.




Test di Levene, per campioni normali

Esistono diverse rielaborazioni di questo test, che si basa sull'analisi della varianza degli scarti dei valori campionari dalla media, o dalla mediana o dalla media-trimmed (gli ultimi due sono anche chiamati test di Levene modificati secondo Brown-Forsyth). A differenze degli altri test visti finora, il test di Levene può essere applicato anche quando i k gruppi hanno diverse numerosità campionarie.
Esaminiamo ad esempio i seguenti 4 gruppi (sono gli stessi di prima, cui è stato tolto l'ultimo elemento dai gruppi a e c). Per prima cosa creiamo il vettore contenente tutti i dati, e il fattore gruppi, nel modo seguente:


a <- c(191, 211, 204, 209)
b <- c(137, 148, 127, 137, 153)
c <- c(172, 165, 186, 178)
d <- c(199, 206, 233, 183, 194)

dati <- c(a,b,c,d)
gruppi = factor(rep(letters[1:4], times=c(length(a),length(b),length(c),length(d))))


A questo punto dobbiamo scaricare la libreria lawstat e applicare la funzione levene.test, scegliendo quale dei 3 test applicare. Personalmente utilizzo il metodo classico basato sugli scarti dalla media, ma possiamo cambiare questa scelta, modificando il parametro location (che può assumere i valori mean, median e trim.mean):


library(lawstat)
levene.test(dati, gruppi, location="mean")

classical Levene's test based on the absolute deviations from the mean ( none not applied
because the location is not set to median )

data: dati
Test Statistic = 0.8531, p-value = 0.4879


Il p-value è ancora una volta maggiore di 0.05, quindi le varianze risultano omogenee.
Se volessimo calcolare manualmente il valore della statistica test (0.8531), dovremmo procedere come segue. Dapprima si calcola la media dei 4 gruppi (oppure la mediana, o la media trimmed). Quindi si calcolano la differenze in valore assoluto (funzione abs) tra le misure e il valore calcolato; infine si avvia una ANOVA a una via, e il risultato della statistica F è il valore che cercavamo:


a <- c(191, 211, 204, 209)
b <- c(137, 148, 127, 137, 153)
c <- c(172, 165, 186, 178)
d <- c(199, 206, 233, 183, 194)

dati <- c(a,b,c,d)
gruppi = factor(rep(letters[1:4], times=c(length(a),length(b),length(c),length(d))))

a.mean <- mean(a)
b.mean <- mean(b)
c.mean <- mean(c)
d.mean <- mean(d)

a.scarto <- abs(a - a.mean)
b.scarto <- abs(b - b.mean)
c.scarto <- abs(c - c.mean)
d.scarto <- abs(d - d.mean)

dati.scarto <- c(a.scarto, b.scarto, c.scarto, d.scarto)

anova(lm(dati.scarto ~ gruppi))

Analysis of Variance Table

Response: dati.scarto
Df Sum Sq Mean Sq F value Pr(>F)
gruppi 3 139.71 46.57 0.8531 0.4879
Residuals 14 764.26 54.59


La statistica F ha valore 0.8531, esattamente quello calcolato con il test di Levene; anche il p-value è lo stesso.
In alternativa, possiamo calcolare le differenze tra i valori e la mediana (median()), o tra i valori e la media trimmed.




Test di Levene non parametrico per campioni non normali

Fino abbiamo visto una lunga serie di test che può essere applicata quando tutti i gruppi seguono una distribuzione normale. Se però siamo interessati a verificare l'omogeneità delle varianze in gruppi che non seguono una disitribuzione non normale, dobbiamo applicare un test non parametrico di verifica della omogeneità delle varianze: è il test di Levene non parametrico.
Fondamentalmente la teoria è molto semplice: calcolati gli scarti dalla media dei valori rilevati, anzichè applicare una ANOVA parametrica ai dati ottenuti (vedi sopra il test di Levene parametrico), si applica l'analisi della varianza non-parametrica, ossia il test di Kruskal-Wallis.

Osserviamo il seguente esempio: le varianze di 4 gruppi, con numerosità campionarie diverse, vengono confrontate:


a <- c(2, 8, 8, 4)
b <- c(8, 7, 14)
c <- c(33, 59, 48, 56)
d <- c(60, 101, 67)

dati <- c(a,b,c,d)
gruppi = factor(rep(letters[1:4], times=c(length(a),length(b),length(c),length(d))))

levene.test(dati, gruppi, location="mean", kruskal.test=T)

rank-based (Kruskal-Wallis) classical Levene's test based on the absolute deviations from the
mean ( none not applied because the location is not set to median )

data: dati
Test Statistic = 6.5705, p-value = 0.08692


Essendo il p-value > 0.05, accettiamo l'ipotesi di omogeneità delle varianze.

Alla luce delle considerazioni fatte sul test di Levene parametrico, osservate il codice seguente, che porta ad ottenere lo stesso risultato della statistica test:


a.mean <- mean(a)
b.mean <- mean(b)
c.mean <- mean(c)
d.mean <- mean(d)

a.scarto <- abs(a - a.mean)
b.scarto <- abs(b - b.mean)
c.scarto <- abs(c - c.mean)
d.scarto <- abs(d - d.mean)

dati = list(g1=a.scarto, g2=b.scarto, g3=c.scarto, g4=d.scarto)

kruskal.test(dati)

Kruskal-Wallis rank sum test

data: dati
Kruskal-Wallis chi-squared = 6.5705, df = 3, p-value = 0.08692

Nessun commento:

Posta un commento