Esta no es una pregunta estrictamente estadística: puedo leer todos los libros de texto sobre los supuestos de ANOVA. Estoy tratando de descubrir cómo los analistas de trabajo reales manejan los datos que no cumplen con los supuestos. Revisé muchas preguntas en este sitio buscando respuestas y sigo encontrando publicaciones sobre cuándo no usar ANOVA (en un contexto matemático abstracto e idealizado) o cómo hacer algunas de las cosas que describo a continuación en R. I Realmente estoy tratando de averiguar qué decisiones la gente realmente toma y por qué.

Estoy ejecutando análisis de datos agrupados de árboles (árboles reales, no árboles estadísticos) en cuatro grupos. Tengo datos de aproximadamente 35 atributos para cada árbol y estoy revisando cada atributo para determinar si los grupos difieren significativamente en ese atributo. Sin embargo, en un par de casos, los supuestos ANOVA se violan ligeramente porque las variaciones no son iguales (según una prueba de Levene, usando alfa = .05).

Tal como lo veo, mis opciones son: 1. Transformar los datos y ver si cambia el p-val de Levene. 2. Use una prueba no paramétrica como un Wilcoxon (si es así, ¿cuál?). 3. Realice algún tipo de corrección al resultado ANOVA, como un Bonferroni (¿no estoy seguro de si existe algo como esto?). He probado las dos primeras opciones y obtuve resultados ligeramente diferentes: en algunos casos, un enfoque es significativo y el otro no. Tengo miedo de caer en la trampa de pesca de valor p, y estoy buscando consejos que me ayuden a justificar qué enfoque usar.

También he leído algunas cosas que sugieren que la heterocedasticidad no es realmente un problema tan grande para ANOVA a menos que los medios y las variaciones estén correlacionados (es decir, ambos aumentan juntos), por lo que tal vez pueda ignorar el resultado de Levene a menos que vea un patrón como ese? Si es así, ¿hay alguna prueba para esto?

Finalmente, debo agregar que estoy haciendo este análisis para su publicación en una revista revisada por pares, por lo que cualquier enfoque que adopte debe pasar a los revisores. Entonces, si alguien puede proporcionar enlaces a ejemplos similares publicados, eso sería fantástico.

Estoy tratando de descubrir cómo los analistas de trabajo reales manejan los datos que no cumplen con los supuestos.

Depende de mis necesidades, qué supuestos se violan, de qué manera, qué tan mal, cuánto afecta la inferencia y, a veces, el tamaño de la muestra.

Estoy ejecutando análisis de datos agrupados de árboles en cuatro grupos. Tengo datos de aproximadamente 35 atributos para cada árbol y estoy revisando cada atributo para determinar si los grupos difieren significativamente en ese atributo. Sin embargo, en un par de casos, los supuestos ANOVA se violan ligeramente porque las variaciones no son iguales (según una prueba de Levene, usando alfa = .05).

1) Si los tamaños de muestra son iguales, no tiene mucho problema. ANOVA es bastante (nivel) robusto a diferentes variaciones si las n son iguales.

2) probar la igualdad de la varianza antes de decidir si asumir que se recomienda en una serie de estudios. Si tiene alguna duda real de que estarán cerca de ser iguales, es mejor simplemente asumir que son desiguales.

Algunas referencias:

Zimmerman, DW (2004),
"Una nota sobre pruebas preliminares de igualdad de varianzas".
Br. J. Math. Stat. Psychol , Mayo ; 57 (Pt 1): 173-81.
http://www.ncbi.nlm.nih.gov/pubmed/15171807

Henrik da tres referencias aquí

3) Lo que importa es el tamaño del efecto, en lugar de si su muestra es lo suficientemente grande como para decirle que son significativamente diferentes. Entonces, en muestras grandes, una pequeña diferencia en la varianza se mostrará como altamente significativa por la prueba de Levene, pero esencialmente no tendrá consecuencias en su impacto. Si las muestras son grandes y el tamaño del efecto (la relación de las variaciones o las diferencias en las variaciones) está bastante cerca de lo que debería ser, entonces el valor p no tiene ninguna consecuencia. (Por otro lado, en muestras pequeñas, un buen valor p grande es de poca comodidad. De cualquier manera, la prueba no responde la pregunta correcta).

Tenga en cuenta que hay un ajuste de tipo Welch-Satterthwaite para la estimación del error estándar residual y df en ANOVA, tal como lo hay en las pruebas t de dos muestras.

2. Use una prueba no paramétrica como un Wilcoxon (si es así, ¿cuál?).

Si está interesado en las alternativas de cambio de ubicación, todavía está asumiendo una extensión constante. Si está interesado en alternativas mucho más generales, tal vez podría considerarlo; la muestra k equivalente a una prueba de Wilcoxon es una prueba de Kruskal-Wallis.

Hacer algún tipo de corrección al resultado ANOVA

Vea mi sugerencia anterior de considerar Welch-Satterthwaite, eso es un 'tipo de corrección'.

(Alternativamente, puede lanzar su ANOVA como un conjunto de pruebas t tipo Welch por pares, en cuyo caso es probable que desee ver un Bonferroni o algo similar)

También he leído algunas cosas que sugieren que la heterocedasticidad no es realmente un problema tan grande para ANOVA a menos que los medios y las variaciones estén correlacionados (es decir, ambos aumentan juntos)

Tendrías que citar algo así. Habiendo examinado varias situaciones con pruebas t, no creo que sea claramente cierto, por lo que me gustaría ver por qué piensan así; quizás la situación está restringida de alguna manera. Sin embargo, sería bueno si fuera así, porque los modelos lineales generalizados pueden ayudar con esa situación.

Finalmente, debo agregar que estoy haciendo este análisis para su publicación en una revista revisada por pares, por lo que cualquier enfoque que adopte debe pasar a los revisores.

Es muy difícil predecir qué podría satisfacer a sus revisores. La mayoría de nosotros no trabajamos con árboles.

En realidad, no es muy difícil manejar la heterocedasticidad en modelos lineales simples (p. Ej., Modelos ANOVA de una o dos vías).

Robustez de ANOVA

Primero, como otros han notado, el ANOVA es increíblemente robusto a las desviaciones del supuesto de variaciones iguales, especialmente si tiene datos aproximadamente equilibrados (igual número de observaciones en cada grupo). Las pruebas preliminares sobre variaciones iguales, por otro lado, no lo son (aunque la prueba de Levene es mucho mejor que la prueba F que se enseña comúnmente en los libros de texto). Como lo expresó George Box:

¡Hacer la prueba preliminar sobre las variaciones es como embarcarse en un bote de remos para averiguar si las condiciones son lo suficientemente tranquilas para que un transatlántico salga del puerto!

Aunque el ANOVA es muy robusto, ya que es muy fácil tener en cuenta la heterocedaticidad, hay pocas razones para no hacerlo.

Pruebas no paramétricas.

Si está realmente interesado en las diferencias de medias , las pruebas no paramétricas (p. Ej., La prueba de Kruskal-Wallis) realmente no sirven de nada. Lo hacen diferencias entre los grupos de prueba, pero lo hacen no en las diferencias generales de ensayo en los medios.

Datos de ejemplo

Generemos un ejemplo simple de datos donde uno quisiera usar ANOVA, pero donde la suposición de varianzas iguales no es cierta.

set.seed(1232)
pop = data.frame(group=c("A","B","C"),
                 mean=c(1,2,5),
                 sd=c(1,3,4))
d = do.call(rbind, rep(list(pop),13))
d$x = rnorm(nrow(d), d$mean, d$sd)

Tenemos tres grupos, con diferencias (claras) tanto en medias como en variaciones:

stripchart(x ~ group, data=d)

ANOVA

No es sorprendente que un ANOVA normal maneje esto bastante bien:

> mod.aov = aov(x ~ group, data=d)
> summary(mod.aov)
            Df Sum Sq Mean Sq F value  Pr(>F)    
group        2  199.4   99.69   13.01 5.6e-05 ***
Residuals   36  275.9    7.66                    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Entonces, ¿qué grupos difieren? Usemos el método HSD de Tukey:

> TukeyHSD(mod.aov)
  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = x ~ group, data = d)

$group
        diff        lwr      upr     p adj
B-A 1.736692 -0.9173128 4.390698 0.2589215
C-A 5.422838  2.7688327 8.076843 0.0000447
C-B 3.686146  1.0321403 6.340151 0.0046867

Con un valor P de 0.26, no podemos reclamar ninguna diferencia (en medias) entre el grupo A y B. E incluso si no tomáramos en cuenta que hicimos tres comparaciones, no obtendríamos un P bajo . valor ( P  = 0.12):

> summary.lm(mod.aov)
[…]
Coefficients:
            Estimate Std. Error t value  Pr(>|t|)    
(Intercept)   0.5098     0.7678   0.664     0.511    
groupB        1.7367     1.0858   1.599     0.118    
groupC        5.4228     1.0858   4.994 0.0000153 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.768 on 36 degrees of freedom

¿Porqué es eso? Sobre la base de la trama, no es una diferencia bastante clara. La razón es que ANOVA asume variaciones iguales en cada grupo y estima una desviación estándar común de 2.77 (mostrada como 'Error estándar residual' en la summary.lmtabla, o puede obtenerla tomando la raíz cuadrada del cuadrado medio residual (7.66) en la tabla ANOVA).

Pero el grupo A tiene una desviación estándar (población) de 1, y esta sobreestimación de 2.77 hace que sea (innecesariamente) difícil obtener resultados estadísticamente significativos, es decir, tenemos una prueba con (demasiado) baja potencia.

'ANOVA' con variaciones desiguales

Entonces, ¿cómo ajustar un modelo adecuado, uno que tenga en cuenta las diferencias en las variaciones? Es fácil en R:

> oneway.test(x ~ group, data=d, var.equal=FALSE)
    One-way analysis of means (not assuming equal variances)

data:  x and group
F = 12.7127, num df = 2.000, denom df = 19.055, p-value = 0.0003107

Por lo tanto, si desea ejecutar un simple 'ANOVA' unidireccional en R sin asumir variaciones iguales, use esta función. Básicamente es una extensión de (Welch) t.test()para dos muestras con variaciones desiguales.

Por desgracia, no funciona con TukeyHSD()(o la mayoría de las otras funciones que se utilizan en aovlos objetos), por lo que incluso si estamos seguros de que no hay diferencias entre los grupos, no sabemos donde están.

Modelando la heteroscedasticidad

La mejor solución es modelar las variaciones explícitamente. Y es muy fácil en R:

> library(nlme)
> mod.gls = gls(x ~ group, data=d,
                weights=varIdent(form= ~ 1 | group))
> anova(mod.gls)
Denom. DF: 36 
            numDF  F-value p-value
(Intercept)     1 16.57316  0.0002
group           2 13.15743  0.0001

Todavía diferencias significativas, por supuesto. Pero ahora las diferencias entre el grupo A y B también se han vuelto estadísticamente significativas ( P  = 0.025):

> summary(mod.gls)
Generalized least squares fit by REML
  Model: x ~ group
  […]
Variance function:
 Structure: Different standard
            deviations per stratum
 Formula: ~1 | group 
 Parameter estimates:
       A        B        C 
1.000000 2.444532 3.913382 

Coefficients:
               Value Std.Error  t-value p-value
(Intercept) 0.509768 0.2816667 1.809829  0.0787
groupB      1.736692 0.7439273 2.334492  0.0253
groupC      5.422838 1.1376880 4.766542  0.0000
[…]
Residual standard error: 1.015564 
Degrees of freedom: 39 total; 36 residual

¡Así que usar un modelo apropiado ayuda! También tenga en cuenta que obtenemos estimaciones de las desviaciones estándar (relativas). La desviación estándar estimada para el grupo A se puede encontrar en la parte inferior de los resultados, 1.02. La desviación estándar estimada del grupo B es 2.44 veces esto, o 2.48, y la desviación estándar estimada del grupo C es similar 3.97 (escriba intervals(mod.gls)para obtener intervalos de confianza para las desviaciones estándar relativas de los grupos B y C).

Corrección para pruebas múltiples

Sin embargo, realmente deberíamos corregir las pruebas múltiples. Esto es fácil usando la biblioteca 'multcomp'. Desafortunadamente, no tiene soporte incorporado para objetos 'gls', por lo que primero tendremos que agregar algunas funciones auxiliares:

model.matrix.gls <- function(object, ...)
    model.matrix(terms(object), data = getData(object), ...)
model.frame.gls <- function(object, ...)
  model.frame(formula(object), data = getData(object), ...)
terms.gls <- function(object, ...)
  terms(model.frame(object),...)

Ahora manos a la obra:

> library(multcomp)
> mod.gls.mc = glht(mod.gls, linfct = mcp(group = "Tukey"))
> summary(mod.gls.mc)
[…]
Linear Hypotheses:
           Estimate Std. Error z value Pr(>|z|)    
B - A == 0   1.7367     0.7439   2.334   0.0480 *  
C - A == 0   5.4228     1.1377   4.767   <0.001 ***
C - B == 0   3.6861     1.2996   2.836   0.0118 *  

Todavía diferencia estadísticamente significativa entre el grupo A y el grupo B! ☺ E incluso podemos obtener intervalos de confianza (simultáneos) para las diferencias entre las medias grupales:

> confint(mod.gls.mc)
[…]
Linear Hypotheses:
           Estimate lwr     upr    
B - A == 0 1.73669  0.01014 3.46324
C - A == 0 5.42284  2.78242 8.06325
C - B == 0 3.68615  0.66984 6.70245

Usando un modelo aproximadamente (aquí exactamente) correcto, ¡podemos confiar en estos resultados!

Tenga en cuenta que para este ejemplo simple, los datos para el grupo C realmente no agregan ninguna información sobre las diferencias entre el grupo A y B, ya que modelamos medias separadas y desviaciones estándar para cada grupo. Podríamos haber utilizado las pruebas t por pares corregidas para comparaciones múltiples:

> pairwise.t.test(d$x, d$group, pool.sd=FALSE)
    Pairwise comparisons using t tests with non-pooled SD 

data:  d$x and d$group 

  A       B      
B 0.03301 -      
C 0.00098 0.02032

P value adjustment method: holm 

Sin embargo, para modelos más complicados, por ejemplo, modelos de dos vías o modelos lineales con muchos predictores, usar GLS (mínimos cuadrados generalizados) y modelar explícitamente las funciones de varianza es la mejor solución.

Y la función de varianza no necesita ser simplemente una constante diferente en cada grupo; podemos imponerle estructura. Por ejemplo, podemos modelar la varianza como una potencia de la media de cada grupo (y, por lo tanto, solo necesitamos estimar un parámetro, el exponente), o tal vez como el logaritmo de uno de los predictores del modelo. Todo esto es muy fácil con GLS (y gls()en R).

Los mínimos cuadrados generalizados son, en mi humilde opinión, una técnica de modelado estadístico muy poco utilizada. En lugar de preocuparse por las desviaciones de los supuestos del modelo, ¡ modele esas desviaciones!

1. De hecho, puede haber alguna transformación de sus datos que produzca una distribución aceptablemente normal. Por supuesto, ahora su inferencia se trata de los datos transformados, no de los datos no transformados.

2. Suponiendo que está hablando de un ANOVA de una vía, la prueba de Kruskal-Wallis es un análogo no paramétrico apropiado al ANOVA de una vía. La prueba de Dunn (no la prueba de suma de rango de variedad de jardín) es quizás la prueba no paramétrica más común apropiada para comparaciones múltiples por pares post hoc , aunque hay otras pruebas como la prueba de Conover-Iman (estrictamente más poderosa que la prueba de Dunn después del rechazo) del kruskal-Wallis), y la prueba Dwass-Steele-Crichtlow-Fligner.

3. $\alpha$

El ANOVA se basa en una relación entre las variaciones dentro del grupo y entre grupos. No estoy completamente seguro de lo que quiere decir con heterocedasticidad en este contexto, pero si se refiere a variaciones desiguales entre grupos, me parece que fundamentalmente romperá la lógica de la hipótesis nula de la prueba.

Una simple consulta de Google Scholar para la "prueba de Dunn" junto con un término general de su disciplina debería devolver muchos ejemplos publicados.

Referencias

Conover, WJ e Iman, RL (1979). En procedimientos de comparaciones múltiples . Informe técnico LA-7677-MS, Laboratorio científico de Los Alamos.

Crichtlow, DE y Fligner, MA (1991). En comparaciones múltiples sin distribución en el análisis de varianza unidireccional . Comunicaciones en estadística: teoría y métodos , 20 (1): 127.

Dunn, DO (1964). Múltiples comparaciones usando sumas de rango . Technometrics , 6 (3): 241–252.

Me parece que está haciendo el trabajo de pies y está haciendo todo lo posible, pero le preocupa que sus esfuerzos no sean lo suficientemente buenos como para pasar su trabajo por encima de los revisores. Muy un problema del mundo real. Creo que todos los investigadores luchan con análisis que parecen estar al límite o incluso francamente violar supuestos de vez en cuando. Después de todo, hay millones de artículos que evalúan, por ejemplo, los efectos del tratamiento en 3 grupos pequeños de ratones con algo así como 6 a 7 ratones en cada grupo. ¡Cómo saber si los supuestos de Anova se cumplen en un documento así!

He revisado una gran cantidad de artículos, especialmente en el campo de la fisiopatología cardiovascular, y en realidad nunca me siento 100% seguro de si puedo confiar en los datos o no en un artículo que leo. Pero para mí como revisor, en realidad tienden a pensar que los problemas pueden surgir en tan muchos niveles en la ciencia que es probable que haya mucho sentido cavar demasiado profundo en las estadísticas - después de todo, todo el conjunto de datos podría fabricarse y yo nunca en un millón de años poder decirlo. En consecuencia, siempre habrá un elemento de confianza en este campo de trabajo, que los investigadores nunca deben abusar.

La sugerencia más real que le daría es que debe pensar todo detenidamente antes de enviar y asegurarse de que podrá responder con sinceridad cualquier pregunta que le hagan los revisores. Siempre y cuando hayas hecho todo lo posible, tus intenciones son honestas y duermes bien por la noche, creo que deberías estar bien.