Come eseguire unANCOVA in R

Lo strumento di base per questo è lm; tieni presente che aov è un wrapper per lm.

In particolare, se hai qualche variabile di raggruppamento (fattore), $ g $ e una covariata continua $ x $ , il modello y ~ x + g si adatterebbe a un modello ANCOVA con effetti principali, mentre y ~ x * g si adatterebbe a un modello che include linterazione con la covariata. aov utilizzerà le stesse formule.

Presta particolare attenzione a Note nella guida su aov.

Per quanto riguarda + vs *, russellpierce lo copre praticamente, ma “ti consiglio di consultare ?lm e ?formula e soprattutto la sezione 11.1 del manuale Unintroduzione a R che viene fornito con R (oppure puoi trovarlo online se non hai capito come trovarlo sul tuo computer; più facilmente, questo implica trovare il menu a discesa “Guida” in R o RStudio).

Commenti

• supponiamo di avere due fattori di gruppo $ g_1, g_2 $ e due covariate $ x_1, x_2 $, con il mio modello $$ y_ {ij} = \ mu + \ alpha_i + \ eta_j + x_ {ij1} \ gamma_1 + x_ {ij2} \ gamma_2 + \ epsilon_ {ij} $$ y ~ g_1 + g_2 + x_1 + x_2 fa lo stesso trucco? I valori F ottenuti contro x_1 e x_2 testano rispettivamente $ \ gamma_1 = 0 $ e $ \ gamma_2 = 0 $?
• Non sono sicuro di come mi sia sfuggito. Sì. …. e se vuoi testare entrambi contemporaneamente, adattali sia con che senza e passa gli oggetti lm adattati a anova (tu ‘ vedremo presto se le dai nellordine sbagliato perché alcune SS saranno negative se lo fai)

Consiglio di ottenere e leggere Discovering Statistics utilizzando R per Field. Ha una bella sezione su ANCOVA.

Per eseguire ANCOVA in R carica i seguenti pacchetti:

car compute.es effects ggplot2 multcomp pastecs WRS 

Se stai usando lm o aov (io uso aov) assicurati di impostare i contrasti utilizzando i “contrasti “prima di eseguire aov o lm. R utilizza contrasti non ortogonali per impostazione predefinita che possono rovinare tutto in un ANCOVA. Se desideri impostare contrasti ortogonali, utilizza:

contrasts(dataname$factorvariable)=contr.poly(# of levels, i.e. 3) 

quindi esegui il modello come

model.1=aov(dv~covariate+factorvariable, data=dataname) 

Per visualizzare il modello usa:

Anova(model.1, type="III") 

Assicurati di usare la “A” maiuscola Anova qui e non anova. Questo darà risultati usando il tipo III SS.

summary.lm(model.1) fornirà un altro riepilogo e includerà lR-sq. produzione.

posth=glht(model.1, linfct=mcp(factorvariable="Tukey")) ##gives the post-hoc Tukey analysis summary(posth) ##shows the output in a nice format. 

Se vuoi testare lomogeneità delle pendenze di regressione puoi anche includere un termine di interazione per IV e covariata. Sarebbe:

model=aov(dv~covariate+IV+covariate:IV, data=dataname) 

Se il termine di interazione è significativo, allora non hai omogeneità.

Commenti

• Perché i contrasti non ortogonali rovinano tutto?
• Per rispondere alla domanda precedente sul ” perché i contrasti non ortogonali confondono tutto a posto “. La risposta è che R per impostazione predefinita è non ortogonale (cioè differenza tra le medie), il che può causare problemi se si desidera vedere il contributo di ogni IV separatamente. Quando specifichiamo contrasti ortogonali diciamo a R che vogliamo che la SS per gli IV ‘ sia completamente partizionata e non si sovrapponga. In questo modo possiamo vedere la variazione attribuita a ciascun predittore in modo pulito e chiaro. Se non specifichi, R utilizza un approccio più liberale al contrasto.
• Perché linteresse per il tipo III SS?

Ecco una documentazione complementare http://goo.gl/yxUZ1R della procedura suggerita da @Butorovich. Inoltre, la mia osservazione è che quando la covariata è binaria, luso di summary (lm.object) darebbe la stessa stima IV generata da Anova (lm.object, type = “III”).

Commenti

• Non è ‘ chiaro che questa dovrebbe essere una risposta. È? In tal caso, modificare per chiarire. Se è una domanda, chiedila facendo clic su ASK QUESTION in alto & chiedendola lì. Allora possiamo aiutarti adeguatamente.
• Daccordo. Il messaggio è stato rivisto come una risposta (complementare) alla precedente.

Usiamo la regressione analisi per creare modelli che descrivono leffetto della variazione nelle variabili predittori sulla variabile di risposta. A volte, se abbiamo una variabile categoriale con valori come Sì / No o Maschio / Femmina ecc., La semplice analisi di regressione fornisce più risultati per ogni valore della variabile categoriale. In tale scenario, possiamo studiare leffetto della variabile categoriale utilizzandola insieme alla variabile predittore e confrontando le linee di regressione per ogni livello della variabile categoriale. Tale analisi viene definita Analisi della covarianza, chiamata anche ANCOVA.

Esempio
Considera i dati incorporati R imposta mtcars. In esso osserviamo che il campo am rappresenta il tipo di trasmissione (automatica o manuale). È una variabile categoriale con valori 0 e 1. Da essa può dipendere anche il valore di miglia per gallone (mpg) oltre al valore della potenza in cavalli (hp). Studiamo leffetto del valore di am sulla regressione tra mpg e hp. Viene fatto utilizzando la funzione aov() seguita dalla funzione anova() per confrontare le regressioni multiple.

Dati di input
Crea un data frame contenente i campi mpg, hp e am dal set di dati mtcars. Qui prendiamo mpg come variabile di risposta, hp come variabile predittore e am come la variabile categoriale.

input <- mtcars[,c("am","mpg","hp")] head(input) 

Quando eseguiamo il codice precedente, produce il seguente risultato:

 am mpg hp Mazda RX4 1 21.0 110 Mazda RX4 Wag 1 21.0 110 Datsun 710 1 22.8 93 Hornet 4 Drive 0 21.4 110 Hornet Sportabout 0 18.7 175 Valiant 0 18.1 105 

ANCOVA Analysis
Creiamo un modello di regressione prendendo hp come variabile predittore e mpg come variabile di risposta tenendo conto dellinterazione tra am e hp.

Modello con interazione tra variabile categoriale e variabile predittore

Crea modello di regressione1

result1 <- aov(mpg~hp*am,data=mtcars) summary(result1) 

Quando eseguiamo il codice precedente, produce il seguente risultato:

 Df Sum Sq Mean Sq F value Pr(>F) hp 1 678.4 678.4 77.391 1.50e-09 *** am 1 202.2 202.2 23.072 4.75e-05 *** hp:am 1 0.0 0.0 0.001 0.981 Residuals 28 245.4 8.8 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Questo risultato mostra che sia la potenza in cavalli che il tipo di trasmissione hanno un effetto significativo sulle miglia per gallone poiché il valore p in entrambi i casi è inferiore a 0,05. Ma linterazione tra queste due variabili non è significativa in quanto il valore p è superiore a 0,05.

Modello senza interazione tra variabile categoriale e variabile predittiva

Crea il modello di regressione2

result2 <- aov(mpg~hp+am,data=mtcars) summary(result2) 

Quando eseguiamo il codice precedente, produce il seguente risultato:

 Df Sum Sq Mean Sq F value Pr(>F) hp 1 678.4 678.4 80.15 7.63e-10 *** am 1 202.2 202.2 23.89 3.46e-05 *** Residuals 29 245.4 8.5 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Questo risultato mostra che sia la potenza in cavalli che il tipo di trasmissione hanno un effetto significativo sulle miglia per gallone poiché il valore p in entrambi i casi è inferiore a 0,05.

Confronto di due modelli
Ora possiamo confrontare i due modelli per concludere se linterazione delle variabili è veramente insignificante. Per questo utilizziamo la funzione anova().

 anova(result1,result2) Model 1: mpg ~ hp * am Model 2: mpg ~ hp + am Res.Df RSS Df Sum of Sq F Pr(>F) 1 28 245.43 2 29 245.44 -1 -0.0052515 6e-04 0.9806 

Poiché il valore p è maggiore di 0,05, concludiamo che linterazione tra potenza del motore e tipo di trasmissione non è significativa. Quindi il chilometraggio per gallone dipenderà in modo simile dalla potenza del motore dellauto sia in modalità di trasmissione automatica che manuale.

Commenti

• Quindi quale è venuto prima, questa risposta o questo post sul punto tutorial? tutorialspoint.com/r/r_analysis_of_covariance.htm . Questa risposta dovrebbe essere cancellata come plagio? O i tutorial puntavano solo a copiare da qui?