Interpretación de los gráficos de residuos para detectar heteroscedasticidad

Question

Estoy intentando medir el efecto de diferentes antibióticos en el crecimiento de E. coli y S.aureus . Como mis datos para ambas bacterias no son normales, he optado por hacer una prueba de Kruskal-Wallis. Pero uno de los supuestos de la prueba de Kruskal-Wallis es que no debe haber heteroscedasticidad. Para comprobar la heteroscedasticidad he hecho algunos gráficos de residuos, pero me cuesta interpretarlos.

Para mi E. coli datos, hice el siguiente gráfico:

¿Demuestra esto que mi E. coli ¿los datos son heteroscedásticos?

Además, cuando hice un gráfico residual para mi S.aureus he obtenido el siguiente resultado. Este parece sugerir que el S.aureus los datos no son heteroscedásticos. ¿Estoy en lo cierto?

Edita: Mi variable de respuesta es el diámetro de las zonas de inhibición del crecimiento (en mm). Por ejemplo, la siguiente tabla muestra mis datos para E. coli .

Hice un anova de Welch en mi E. coli datos pero no obtuve ningún valor p. ¿Puede alguien explicarme qué estoy haciendo mal?

One-way analysis of means (not assuming equal variances)

data:  values and ind
F = NaN, num df = 17, denom df = NaN, p-value = NA
```

Answer 1

En realidad, la prueba de Kruskal-Wallis no tiene un supuesto de homocedasticidad. Los supuestos para la prueba dados por Conover † son los siguientes.

1. Todas las muestras son aleatorias de sus respectivas poblaciones.

2. Además de la independencia dentro de cada muestra, existe independencia mutua entre las distintas muestras.

3. La escala de medición es, como mínimo, ordinal.

4. O bien el k las funciones de distribución de la población son idénticas, o bien algunas de las poblaciones tienden a dar valores mayores que otras poblaciones.

† Conover, Practical Nonparametric Statistics, 3rd.

---

Edita: Dado que OP ha proporcionado algunos datos, puedo añadir más a esta respuesta.

Como se ha discutido en algunos comentarios, los datos que tiene OP son por naturaleza continuos (diámetro en mm), pero tienen poca precisión, por lo que actúan más como una variable discreta. En este caso, probablemente me vería tentado a utilizar un modelo para valores discretos. (¿Quizás una regresión binomial negativa?) Pero como los datos son inusuales en el sentido de que algunos grupos no tienen varianza, y otros tienen considerablemente más, podría sentirme mejor con el enfoque no paramétrico.

A continuación se realiza una prueba de Kruskal-Wallis y, después, una prueba post-hoc de Dunn (1964). No hice p correcciones de valor aquí para múltiples pruebas, aunque es posible que desee.

Además, al final, hay una función que informa de las estadísticas del tamaño del efecto por pares, incluidas las de Vargha y Delaney A . Esta estadística del tamaño del efecto está relacionada con la probabilidad de que una observación en un grupo sea mayor que una observación en el otro grupo. Por ejemplo, comparando C30 con el control, vda=0, y de hecho cada observación en C30 es mayor que cada observación en el control, pero estos dos tratamientos no son significativamente diferentes en la prueba de Dunn. Depende de usted cómo utilice esta información.

if(!require(tidyr)){install.packages("tidyr")}
if(!require(FSA)){install.packages("FSA")}
if(!require(PMCMRplus)){install.packages("PMCMRplus")}
if(!require(rcompanion)){install.packages("rcompanion")}

### Import data and translate it to long format

Data = read.table(header=T, text="
Control KF30 P10  APR15 x510 TE30   AK30  C30   TOB10  CN10   N3     efiran Bleach Soap ETOH70 ETOH96 Detergent Mouthwash
6       23  6    15    11    11     18    7     17     16     17     16     6      6     6     9     6          15
6       10  6    15    12    22     18    7     16     16     17     15     6      6     6     7     6          6
6       11  6    15    12    22     18    8     15     15     17     14     6      6     6     6     6          6
6       10  6    17    12    16     21    8     18     20     19     19     6      6     6     6     6          6
6       12  6    15    11    23     17    8     16     16     15     16     7      6     7     7     6          6
6       12  6    20    10    22     18    8     16     16     16     19     6      6     6     7     6          6
6       11  6    16    12    18     19    8     15     12     18     14     12     6     6     7     6          6
6       12  6    16    13    21     20    9     17     15     17     13     12     7     6     7     6          6
6       12  6    15    12    16     20    8     16     15     16     30     10     6     6     10    6          6
")

library(tidyr)

DataLong <- gather(Data, Treatment, Diameter, Control:Mouthwash, factor_key=TRUE)

str(DataLong)

   ### 'data.frame':    162 obs. of  2 variables:
   ###  $ Treatment: Factor w/ 18 levels "Control","KF30",..: 1 1 1 1 1 1 1 1 1 2 ...
   ###  $ Diameter : int  6 6 6 6 6 6 6 6 6 23 ...

### Plot and summarize data

library(FSA)

Summarize(Diameter ~ Treatment, data=DataLong)

   ###    Treatment n      mean        sd min Q1 median Q3 max
   ### 1    Control 9  6.000000 0.0000000   6  6      6  6   6
   ### 2       KF30 9 12.555556 4.0034707  10 11     12 12  23
   ### 3        P10 9  6.000000 0.0000000   6  6      6  6   6
   ### 4      APR15 9 16.000000 1.6583124  15 15     15 16  20
   ### 5       X510 9 11.666667 0.8660254  10 11     12 12  13
   ### 6       TE30 9 19.000000 4.0311289  11 16     21 22  23
   ### 7       AK30 9 18.777778 1.3017083  17 18     18 20  21
   ### 8        C30 9  7.888889 0.6009252   7  8      8  8   9
   ### 9      TOB10 9 16.222222 0.9718253  15 16     16 17  18
   ### 10      CN10 9 15.666667 2.0615528  12 15     16 16  20
   ### 11        N3 9 16.888889 1.1666667  15 16     17 17  19
   ### 12    efiran 9 17.333333 5.1961524  13 14     16 19  30
   ### 13    Bleach 9  7.888889 2.6666667   6  6      6 10  12
   ### 14      Soap 9  6.111111 0.3333333   6  6      6  6   7
   ### 15    ETOH70 9  6.111111 0.3333333   6  6      6  6   7
   ### 16    ETOH96 9  7.333333 1.3228757   6  7      7  7  10
   ### 17 Detergent 9  6.000000 0.0000000   6  6      6  6   6
   ### 18 Mouthwash 9  7.000000 3.0000000   6  6      6  6  15

plot(Diameter ~ Treatment, data=DataLong)

   ### (Plot)

### Kruskal-Wallis and Dunn tests

kruskal.test(Diameter ~ Treatment, data=DataLong)

   ### Kruskal-Wallis rank sum test
   ### 
   ### Kruskal-Wallis chi-squared = 142.47, df = 17, p-value < 2.2e-16

### Order groups by median

Sum = Summarize(Diameter ~ Treatment, data=DataLong)

Sum2 = Sum[order(Sum$median),]

DataLong$Treatment = factor(DataLong$Treatment, levels=Sum2$Treatment)    

plot(Diameter ~ Treatment, data=DataLong)

   ### (Plot with groups ordered by medians)

library(PMCMRplus)

Dunn = kwAllPairsDunnTest(Diameter ~ Treatment, data = DataLong, p.adjust.method = "none")

Dunn

   ### (Large output table)

summaryGroup(Dunn)

   ###           median Q25 Q75 n Sig. group
   ### Control        6   6   6 9          a
   ### P10            6   6   6 9          a
   ### Bleach         6   6  10 9         ab
   ### Soap           6   6   6 9          a
   ### ETOH70         6   6   6 9          a
   ### Detergent      6   6   6 9          a
   ### Mouthwash      6   6   6 9          a
   ### ETOH96         7   7   7 9        abc
   ### C30            8   8   8 9        abc
   ### KF30          12  11  12 9         cd
   ### x510          12  11  12 9        bcd
   ### APR15         15  15  16 9         de
   ### TOB10         16  16  17 9         de
   ### CN10          16  15  16 9         de
   ### efiran        16  14  19 9         de
   ### N3            17  16  17 9         de
   ### AK30          18  18  20 9          e
   ### TE30          21  16  22 9          e

### Report Vargha and Delaney's A for pairwise comparisons

multiVDA(Diameter ~ Treatment, data=DataLong)$pairs

   ### (Large table) 

   ### Note with vda, 0.50 = no effect; 0 or 1 = stochastic dominance

Answer 2

Parece que en ambos gráficos, sufren de heteroscedasticidad, aunque este último está en mejor forma. Bueno, el test de Kruskal Wallis no asume homoscedasticidad, sin embargo es mejor no usarlo, porque asume que las distribuciones entre las poblaciones son idénticas, así que tener heteroscedasticidad no ayuda. En primer lugar, debería probar la prueba de Bartlett para comprobar la heteroscedasticidad de las muestras y, a continuación, la anova de Welch como alternativa. Intente estudiar aquí como material de consulta. Algunos documentos realmente buenos deberían serlo:

Prueba de Bartlett

Prueba de Welch .

La última fuente que proporcioné trata sobre cómo la prueba de Welch es robusta al error de tipo I, que es el principal problema que causa la heteroscedasticidad.