Diseños balanceados con efectos fijos
Gerardo Martín
9/3/2021
1 Introducción
Como vimos anteriormente los diseños balanceados son aquellos que cuentan con la misma cantidad de unidades observacionales y experimentales por tratamiento y bloque. Aquí comenzaremos a ver cómo se ajusta e interpreta un ANOVA de dos vías en R.
Suponiendo que el diseño experimental únicamente contempla dos criterios de clasificación \(A\) y \(B\), podemos ajustar un ANOVA de dos vías utilizando aov(y ~ A * B, data = datos). Recordemos que A * B en R quiere decir A + B + A : B.
2 Ejemplo de ANOVA de dos vías para un diseño balanceado
Para este caso utilizaremos una base de datos reales que se distribuye por default en R llamada ToothGrowth. Esta base proviene de un estudio del efecto de la vitamina C sobre el crecimiento de los dientes de los conejillos de indias. Los conejillos de indias, como los humanso, no sintetizan vitamina C, por lo que son muy susceptiles de padecer escorbuto, una enfermedad metabólica por acumulación de radicales libres cuyos primeros síntomas son la pérdida de piezas dentales.
Para ver la base de datos podemos escribir el nombre del objeto que de la base de datos en la consola de R:
str(ToothGrowth)## 'data.frame': 60 obs. of 3 variables:
## $ len : num 4.2 11.5 7.3 5.8 6.4 10 11.2 11.2 5.2 7 ...
## $ supp: Factor w/ 2 levels "OJ","VC": 2 2 2 2 2 2 2 2 2 2 ...
## $ dose: num 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 ...Y podemos ver que consta de tres columnas y 60 filas (observaciones), la primera len contiene las mediciones de los dientes, supp es el tratamiento de suplementación con vitamina C (VC) y jugo de naranja (OJ), y la tercera dose contiene la dosis de 0.5, 1.0 y 2.0. Esta última es contínua pero para hacer ANOVA la transformaremos en categórica.
Comenzaremos por almancenar ToothGrowth en otro objeto y transformar dose en un factor:
datos <- ToothGrowth
datos$dose <- as.factor(paste0("D", datos$dose)) #Transformando dose en un factor
str(datos)## 'data.frame': 60 obs. of 3 variables:
## $ len : num 4.2 11.5 7.3 5.8 6.4 10 11.2 11.2 5.2 7 ...
## $ supp: Factor w/ 2 levels "OJ","VC": 2 2 2 2 2 2 2 2 2 2 ...
## $ dose: Factor w/ 3 levels "D0.5","D1","D2": 1 1 1 1 1 1 1 1 1 1 ...Podemos ahora ver si ToothGrowth tiene datos de un estudio balanceado, verifigando que todos los niveles de los factores supp y dose tengan el mismo número de observaciones. Para ello haremos una tabulación cruzada:
table(datos$supp)##
## OJ VC
## 30 30table(datos$dose)##
## D0.5 D1 D2
## 20 20 20Entonces, lo que nos queda por ver es si len se presta para hacer ANOVA:
Figure 2.1: Gráfico de cajas de la longitud de dientes en ToothGrowth
Es posible que haya grupos que se desvíen de la normalidad, pero podemos ver, después de ajustado el modelo, si esto es un problema.
2.1 Ajustando el modelo
Para comenzar hay que entender la hipótesis nula que estaremos probando. El primer criterio de clasificación es la suplementación con jugo de naranja y con vitamina C, el segundo es la dosis de cada uno. De modo que trataremos de rechazar la hipótesis nula de que el suplemento y su dosis no afectan el crecimiento de los dientes.
Comenzaremos por ajustar un modelo donde hay interacción entre supp y dose. Para mantener la interpretación sencilla lo haremos con aov:
modelo.1 <- aov(len ~ supp * dose, data = datos)
summary(modelo.1)## Df Sum Sq Mean Sq F value Pr(>F)
## supp 1 205.4 205.4 15.572 0.000231 ***
## dose 2 2426.4 1213.2 92.000 < 2e-16 ***
## supp:dose 2 108.3 54.2 4.107 0.021860 *
## Residuals 54 712.1 13.2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Podemos ver que todos los componentes del modelo afectan significativamente la media de los tratamientos, con lo que de entrada podemos rechazar \(H_0\). Sin embargo, si vemos con cuidado la suma de cuadrados (columna Sum Sq), el factor dosis explica una gran proporción de la variación, pues contiene la suma de cuadrados con el mayor valor, y tiene la \(P\) más pequeña de todas. Para ver si los otros factores son importantes en realidad, podemos ajustar otro modelo quitando el término o factor que contribuyó con la menor suma de cuadrados que fue supp:dose:
modelo.2 <- aov(len ~ supp + dose, data = datos)
summary(modelo.2)## Df Sum Sq Mean Sq F value Pr(>F)
## supp 1 205.4 205.4 14.02 0.000429 ***
## dose 2 2426.4 1213.2 82.81 < 2e-16 ***
## Residuals 56 820.4 14.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Aquí vemos que la suma de cuadrados de los residuales aumentó un poco y bajó la significancia de supp. Esta diferencia entre modelos sugiere que el término de interacción es importante, aunque la dosis es el principal factor que afecta el crecimiento de los dientes. Chequemos ahora los residuales de ambos modelos.
2.1.1 Invirtiendo el orden de los factores
Como se mencionó en Introducción al modelo ANOVA de dos vías, el orden de especificación de los factores de primer orden (antes de las interacciones especificaddas con :), afecta los resultados, por lo que ahora ajustaremos el modelo haciendo dose * supp:
modelo.3 <- aov(len ~ dose * supp, data = datos)
summary(modelo.3)## Df Sum Sq Mean Sq F value Pr(>F)
## dose 2 2426.4 1213.2 92.000 < 2e-16 ***
## supp 1 205.4 205.4 15.572 0.000231 ***
## dose:supp 2 108.3 54.2 4.107 0.021860 *
## Residuals 54 712.1 13.2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Y efectivamente, los valores de \(P\) para supp son diferentes que cuando especificamos el modelo como supp * dose, pues la base de comparación para supp ahora es dose. También podemos apreciar que \(P\) para supp : dose es igual, y esto se debe a que tiene la misma base de comparación. El valor idéndico de \(P\) para dose se debe a que el número más pequeño que R (la computadora) puede mostrar es \(2 \times 10^{-16}\).
2.1.2 Diagnóstico de modelo.1
par(mfrow = c(2,2))
plot(modelo.1)
Figure 2.2: Gráficas de diagnóstico de modelo.1.
A juzgar por estas gráficas es posible que el modelo sea confiable para rechazar \(H_0\), veamos sin embargo si modelo.2 tiene otras cualidades.
2.1.3 Diagnósticos de modelo.2
par(mfrow = c(2,2))
plot(modelo.2)
Figure 2.3: Gráficas de diagnóstico de modelo.2.
Para este modelo podemos que que los residuales no están centrados en cero y que la línea roja tampoco es recta o completamente horizontal. Esto corrobora que modelo.1 es una mejor opción para rechazar \(H_0\) y que sí podemos concluir que la dosis del suplemento afecta el crecimiento de los dientes.