giovedì 30 luglio 2009

Test del chi-quadro per il trend (test di Armitage)

Test di verifica di un trend lineare in proporzioni e frequenze.

Supponiamo di disporre di una tabella di contingenza 2xK, ad esempio relativa al numero di cavie morte dopo somministrazione di dosi crescenti di un certo farmaco sperimentale, oppure relativa alll'incidenza di patologie in gruppi di individui appartenenti a classi di età sempre maggiori. Con il test di Armitage (Cochran-Armitage test for trend; in alcuni software come SPSS è chiamato anche Mantel-Haenszel linear-by-linear association chi-squared test) possiamo verificare se è presente o meno un trend lineare delle proporzioni. L'ipotesi nulla è la presenza di un trend lineare, contro l'ipotesi alternativa di allontanamento da esso. Vediamo come procedere analizzando la seguente tabella:



Per prima cosa dobbiamo verificare se le proporzioni sono tra loro simili, con il test del Chi-quadro per l'indipendenza delle variabili, già presentato in questo post.

morti <- c(19, 29, 24, 20)
tot <- c(516, 589, 293, 167)
prop.test(morti, tot)

4-sample test for equality of proportions without continuity correction

data: con out of tot
X-squared = 19.5233, df = 3, p-value = 0.0002131
alternative hypothesis: two.sided
sample estimates:
    prop 1     prop 2     prop 3     prop 4
0.03682171 0.04923599 0.08191126 0.11976048


Il Chi-quadro ottenuto è pari a 19.5233, e p-value < 0.05. Quindi rifiutiamo l'ipotesi nulla di uguaglianza delle proporzioni. Le 4 proporzioni sono statisticamente differenti.

Procediamo con un test di correlazione per variabili qualitative, con il test del chi-quadro di Pearson illustrato in questo post. Procediamo in questo modo:

tabella = matrix(c(19, 29, 24, 20, 497, 580, 269, 147), nrow=2, byrow=TRUE)
chisq.test(tabella)

        Pearson's Chi-squared test

data:  tabella
X-squared = 20.1478, df = 3, p-value = 0.0001582


Il risultato di questo test ci dice che esiste una correlazione tra le variabili. Quindi è lecito supporre un andamento lineare di queste, e perciò esaminiamo la presenza di un trend lineare con il test di Armitage:

morti <- c(19, 29, 24, 20)
tot <- c(516, 589, 293, 167)
prop.trend.test(morti, tot)

        Chi-squared Test for Trend in Proportions

data:  morti out of tot ,
using scores: 1 2 3 4
X-squared = 18.2109, df = 1, p-value = 1.977e-05


Abbiamo ottenuto un chi-quadro elevato, e p-value < 0.05. Rifiuto l'ipotesi nulla, e quindi posso affermare che esiste un trend lineare significativo delle proporzioni.
Abbiamo quindi individuato che la componente lineare spiega l'andamento delle proporzioni. Ma è questa l'unica spiegazione all'andamento delle proporzioni, o esiste un'altra regola sottesa? Analizziamo la deviazione dalla linearità, calcolando il chi-quadro residuo, che è pari al chi-quadro totale (cioè il chi-quadro di Pearson), meno il chi-quadro lineare (cioè il chi-quadro di Armitage):

corr <- chisq.test(tabella)
trend <- prop.trend.test(morti, tot)

chi_quadro_residuo <- corr$statistic - trend$statistic

> chi_quadro_residuo

X-squared
1.936932


Il chi-quadro residuo ha gradi di libertà pari a k-2, quindi 2 in questo caso. Calcoliamo quindi il chi-quadro tabulato:

> qchisq(0.950, 2)
[1] 5.991465


Poichè il chi-quadro residuo è inferiore al valore tabulato, accetto l'ipotesi nulla: la componente residua non è significativa, quindi la componente lineare è la sola a spiegare l'andamento delle proporzioni (se avessi rifiutato, significherebbe che la componente lineare non è l'unica variabile esplicativa, ma ne esistono anche altre da ricercare).

9 commenti:

  1. Ciao,
    penso che i dati portati in tabella colonna dose2 siano errati ( 29+580 !=589).
    Davide

    RispondiElimina
  2. @ Davide
    Hai ragione, grazie per la segnalazione. Provvedo subito ad aggiustare il tutto :)

    RispondiElimina
    Risposte
    1. Il database sembra parzialmente derivare da quello riportato qui:
      http://www.sixsigmain.it/ebook/Capu21-17.html

      L'errore in questo caso sarebbe il 580 che doveva essere 560...

      Elimina
  3. Correggendo i dati, il risultato di prop.test e chisq.test è lo stesso.
    Da una mail di R-help sembra non ci sia una differenza tra i due test.
    Come mai li hai riportati in questa sequenza logica? Ho l'impressione che mi sfugga qualcosa...

    RispondiElimina
    Risposte
    1. Forse intendeva controllare che le proporzioni all'interno dei gruppi fossero diverse?

      In questo caso l'analisi avrebbe dovuto essere:

      rowSums(tabella)
      [1] 92 1473

      totali<-rowSums(tabella)

      prop.test(tabella1[1,],rep(totali[1],4))

      4-sample test for equality of proportions without continuity
      correction

      data: tabella[1, ] out of rep(totali[1], 4)
      X-squared = 3.5942, df = 3, p-value = 0.3087
      alternative hypothesis: two.sided
      sample estimates:
      prop 1 prop 2 prop 3 prop 4
      0.2065217 0.3152174 0.2608696 0.2173913

      > prop.test(tabella1[2,],rep(totali[2],4))

      4-sample test for equality of proportions without continuity correction

      data: tabella1[2, ] out of rep(totali[2], 4)
      X-squared = 406.05, df = 3, p-value < 0.00000000000000022
      alternative hypothesis: two.sided
      sample estimates:
      prop 1 prop 2 prop 3 prop 4
      0.33740665 0.38017651 0.18262050 0.09979633

      Significatività per le proporzioni di vivi tra le dosi crescenti...

      Elimina
  4. Ciao! Probabilmente è una domanda molto stupida ma sono ancora alle prime armi con R..
    Devo testare la presenza del trend su un data set composto da 6 gruppi sperimentali di numerosità diversa, per ognuno è riportato il conteggio dei tumori per tipologia e sede: esiste un modo, caricando l'intero set di dati, di non dover creare per ogni tumore la singola tabella di contingenza ed eseguire il test singolarmente??

    RispondiElimina
    Risposte
    1. Direi di sì.

      Devi far ricorso però alla funzione
      independence_test() del pacchetto coin.

      Es.
      library(coin)
      independence_test(TUMOR.type-site ~ GROUPS, teststat = "quad",
      scores = list(GROUPS = c(1, 2, 3, 4, 5, 6))))
      Con il parametro scores indichi l'ordinamento dei gruppi per cui deve essere verificato il trend lineare.

      pdeninis

      Elimina
  5. prima si legge che "L'ipotesi nulla è la presenza di un trend lineare, contro l'ipotesi alternativa di allontanamento da esso" mentre poi si legge che "Abbiamo ottenuto un chi-quadro elevato, e p-value < 0.05. Rifiuto l'ipotesi nulla, e quindi posso affermare che esiste un trend lineare significativo delle proporzioni". C'è qualcosa che non va...

    RispondiElimina
  6. in effetti, l'ipotesi nulla dovrebbe essere assenza di un trend. Probabilmente è un refuso.

    RispondiElimina