¿Cómo se hace un ANOVA SS de tipo III en R con códigos de contraste?

Question

Por favor, proporcione el código de R que permite realizar un ANOVA entre sujetos con contrastes de -3, -1, 1, 3. Entiendo que hay un debate sobre el tipo de suma de cuadrados (SS) apropiado para dicho análisis. Sin embargo, como el tipo de SS por defecto utilizado en SAS y SPSS (Tipo III) se considera el estándar en mi área. Por lo tanto, me gustaría que los resultados de este análisis coincidieran perfectamente con los generados por esos programas estadísticos. Para ser aceptada una respuesta debe llamar directamente a aov(), pero otras respuestas pueden ser votadas (especialmente si son fáciles de entender/usar).

sample.data <- data.frame(IV=rep(1:4,each=20),DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))

Editar: Tenga en cuenta que el contraste que solicito no es un simple contraste lineal o polinómico, sino que se trata de un contraste derivado de una predicción teórica, es decir, el tipo de contrastes analizados por Rosenthal y Rosnow.

Answer 1

Las sumas de cuadrados de tipo III para el ANOVA están disponibles a través del Anova() de la función coche paquete.

La codificación del contraste puede hacerse de varias maneras, utilizando C() El contr.* (como indica @nico), o directamente el contrasts() función/argumento. Esto se detalla en §6.2 (pp. 144-151) de Estadística Aplicada Moderna con S (Springer, 2002, 4ª edición). Tenga en cuenta que aov() es sólo una función que envuelve al lm() función. Es interesante cuando se quiere controlar el término de error del modelo (como en un diseño dentro de un sujeto), pero por lo demás ambos dan los mismos resultados (y sea cual sea la forma de ajustar el modelo, se pueden obtener resúmenes tipo ANOVA o LM con summary.aov o summary.lm ).

No tengo SPSS para comparar los dos resultados, pero algo así como

> library(car)
> sample.data <- data.frame(IV=factor(rep(1:4,each=20)),
                            DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))
> Anova(lm1 <- lm(DV ~ IV, data=sample.data, 
                  contrasts=list(IV=contr.poly)), type="III")
Anova Table (Type III tests)

Response: DV
            Sum Sq Df F value    Pr(>F)    
(Intercept)  18.08  1  21.815  1.27e-05 ***
IV          567.05  3 228.046 < 2.2e-16 ***
Residuals    62.99 76                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

vale la pena probar en primera instancia.

Acerca de la codificación de factores en R frente a SAS: R considera la línea de base o el nivel de referencia como el primer nivel en el orden lexicográfico, mientras que SAS considera el último. Por lo tanto, para obtener resultados comparables, hay que utilizar contr.SAS() ou à relevel() su factor R.

Answer 2

Esto puede parecer un poco de autopromoción (y supongo que lo es). Pero he desarrollado una lsmeans para R (disponible en CRAN) que está diseñado para manejar exactamente este tipo de situaciones. Así es como funciona para su ejemplo:

> sample.data <- data.frame(IV=rep(1:4,each=20),DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))
> sample.aov <- aov(DV ~ factor(IV), data = sample.data)

> library("lsmeans")
> (sample.lsm <- lsmeans(sample.aov, "IV"))
 IV    lsmean        SE df   lower.CL  upper.CL
  1 -3.009669 0.2237448 76 -3.4552957 -2.564043
  2 -3.046072 0.2237448 76 -3.4916980 -2.600445
  3  1.147080 0.2237448 76  0.7014539  1.592707
  4  3.049153 0.2237448 76  2.6035264  3.494779

> contrast(sample.lsm, list(mycon = c(-3,-1,1,3)))
 contrast estimate       SE df t.ratio p.value
 mycon    22.36962 1.000617 76  22.356  <.0001

Si lo desea, puede especificar otros contrastes en la lista. Para este ejemplo, obtendrá los mismos resultados con el contraste polinómico lineal incorporado:

> con <- contrast(sample.lsm, "poly")
> con
 contrast   estimate        SE df t.ratio p.value
 linear    22.369618 1.0006172 76  22.356  <.0001
 quadratic  1.938475 0.4474896 76   4.332  <.0001
 cubic     -6.520633 1.0006172 76  -6.517  <.0001

Para confirmarlo, observe que el "poly" le indica que llame a poly.lsmc que produce estos resultados:

> poly.lsmc(1:4)
  linear quadratic cubic
1     -3         1    -1
2     -1        -1     3
3      1        -1    -3
4      3         1     1

Si desea realizar una prueba conjunta de varios contrastes, utilice la función test con la función joint = TRUE . Por ejemplo,

> test(con, joint = TRUE)

Esto producirá una prueba de "tipo III". A diferencia de car::Anova() lo hará correctamente independientemente de la codificación de contraste utilizada en la etapa de ajuste del modelo. Esto se debe a que las funciones lineales que se prueban se especifican directamente y no implícitamente a través de la reducción del modelo. Una característica adicional es que se detecta un caso en el que los contrastes que se prueban son linealmente dependientes, y se produce la estadística de prueba y los grados de libertad correctos.

Answer 3

Quizás quieras echar un vistazo a esta entrada del blog:

Obtención de los mismos resultados de ANOVA en R que en SPSS - las dificultades con las sumas de cuadrados de Tipo II y Tipo III

( Spoiler: añadir options(contrasts=c("contr.sum", "contr.poly")) al principio de su guión)

Answer 4

Cuando se hacen contrastes, se hace una combinación lineal específica y declarada de las medias de las celdas en el contexto del término de error apropiado. Como tal, el concepto de "Tipo de SS" no tiene sentido con los contrastes. Cada contraste es esencialmente el primer efecto utilizando un SS de tipo I. El "Tipo de SS" tiene que ver con lo que se parcializa o se contabiliza por los otros términos. En el caso de los contrastes, no se excluye ni se tiene en cuenta nada. El contraste se mantiene por sí mismo.

Answer 5

El hecho de que las pruebas de tipo III se utilicen en su lugar de trabajo es la más débil de las razones para seguir utilizándolas. SAS ha hecho un gran daño a la estadística en este sentido. La exégesis de Bill Venables, a la que se hace referencia más arriba, es un gran recurso al respecto. Simplemente diga no al tipo III; se basa en una noción defectuosa de equilibrio y tiene menor potencia debido a la ponderación tonta de las celdas en el caso de desequilibrio.

Una forma más natural y menos propensa a errores para obtener contrastes generales, y poder describir lo que se hizo, es la que proporciona el programa R rms paquete contrast.rms función. Los contrastes pueden ser muy complejos, pero para el usuario son muy sencillos porque se expresan en términos de diferencias de valores predictivos. Se admiten pruebas y contrastes simultáneos. Esto maneja efectos de regresión no lineales, efectos de interacción no lineales, contrastes parciales, todo tipo de cosas.