5 ANOVA
5.1 Introducción
El análisis de la varianza fue desarrollado por Fisher en 1925 con el objetivo de comparar entre sí varios grupos o tratamientos. Normalmente es de interés estudiar cómo se diferencian estos grupos (o tratamientos), teniendo en cuenta la incidencia de otros factores cualitativos o cuantitativos (factores ambientales), cuya influencia es eliminada mediante una adecuada descomposición de la variabilidad total de un experimento en componentes independientes que puedan atribuirse a distintas causas.
Por ejemplo, una compañía agricola que emplea diversos fertilizantes desea comprobar si éstos tienen efectos diferentes sobre el rendimiento de la planta de jitomate en donde se usan estos fertilizantes. La compañía tiene a su disposición 5 tipos de fertilizantes y para comparar su eficacia fumiga, con cada uno de los fertilizantes, un cierto número de parcelas de terreno de la misma calidad y de igual superficie. Cabe mencionar que el número de parcelas a las que les corresponde cada fertilizante no ncesariamnete deben ser las mismas, es decir el fertilizante 1 se aplica a 8 parcelas mientras que el fertilizante 2 se aplicó sólo en 6. Al recoger la cosecha se mide el rendimiento de la semilla en peso de jitomates por unidad de superficie. Este peso de jitomates nos permitirá saber si los fertilizantes están teniendo efecto sobre el rendimiento. En este ejemplo se esta tratando de controlar los efectos de los otros factores al tener parcelas de igual superficie, en condiciones lo más homogéneas posible y con la planta derivada de la misma semilla, de esta maner se espera que de observarse diferencias en las cosechas se puedan atrubuir solamente al fertilizante.
Esencialmente el análisis de varianza determina si la discrepancia entre las medias de los tratamientos es mayor de lo que podría esperarse razonablemente de la discrepancia existente dentro de los tratamientos. En este tipo de modelos, la variable respuesta se considera del tipo continuo, mientras que las variables experimentales o factores son variables categóricas o categorizadas en niveles.
Similar al proceso de inferencia, para llegar a conclusiones estadísticas correctas es necesario observar el resultado tras la repetición del experimento en varias unidades experimentales para cada una de las diversas condiciones que indica el diseño pero lo más homogéneas posibles dentro de cada una. Esto reducirá la variabilidad y, por tanto, aumentará la capacidad estadística de detectar cambios o identificar variables influyentes. Ya que con una variabilidad muy grande respecto al error experimental no se pueden detectar diferencias entre tratamientos.
En general, en todo análisis de varianza es necesario considerar tres etapas:
Diseño del experimento a fin de obtener observaciones de una variable Y, combinando adecuadamente los factores incidentes.
Planteamiento de hipótesis, cálculo de sumas de cuadrados (residuales, de desviación de la hipótesis, etc.) y obtención de los cocientes F. Esta parte del análisis se formula mediante la teoría de los modelos lineales.
Toma de decisiones e interpretación de los resultados. Planteamiento “a posteriori” de nuevas hipótesis.
5.2 Modelo de una vía
Supongamos que una variable \(Y\) ha sido observada bajo \(k\) condiciones experimentales distintas. Puede ser que las observaciones provengan de \(k\) poblaciones, o bien tratarse de réplicas para cada uno de los \(k\) niveles de un factor.
Indiquemos por \(y_{ij}\) la réplica \(j (j = 1, . . . , n_i)\) en la población o nivel \(i (i = 1, . . . , k)\), donde \(n_i\) es el número de réplicas en la población o nivel \(i\). Entonce el conjunto de datos se puede ver como:
\[ \begin{align*} Nivel 1 \quad y_{11}, y_{12}, . . . , y_{1n_1}\\ Nivel 2 \quad y_{21}, y_{22}, . . . , y_{2n_2}\\ .\\ .\\ .\\ Nivel k \quad y_{k1}, y_{k2}, . . . , y_{kn_k} \end{align*} \]
Con estos datos se pueden calcular las siguientes medias:
Media en la población o nivel \(i\):
\[y_{i.} =\frac{1}{n_i}\sum_{j=1}^{n_i}y_{ij}\]
Media general:
\[\bar{y} = y_{\cdot\cdot} =\frac{1}{n}\sum_{i=1}^{k}\sum_{j=1}^{n_i}y_{ij}\]
donde \(n=\sum_{i=1}^{k}n_i\) es el número total de observaciones.
El objetivo es estimar los efectos de los distintos \(k\) niveles y contrastar la hipótesis de que todos los niveles del factor producen el mismo efecto, frente a la alternativa de que al menos dos difieren significativamente entre sí.
Es importante mencionar que con este diseño de los datos es posible proponer distintos modelos lineales que en teoría o en la formulación matemática se adaptan al diseño de los datos antes mencionado. Por ejemplo se puede proponer el siguiente modelo lineal.
\[y_{ij}=\mu_i+\epsilon_{ij},\ i=1,...,k;\ j=1,...,n_i\]
siendo \((\mu_1,\mu_2,...\mu_k)'\) el vector de parámetros y \(\epsilon_{ij}\) siendo las variables aleatorias que agrupan información proveniente de la combinación del factor y la observación (o lo no considerado en \(\mu_i\)). En este modelo podemos considerar la matriz \(X\) definida por una matriz identidad de \(k\)x\(k\) y entonce nos lleva a un modelo de la forma \[Y=\mu X+\epsilon\]
Esta es la representación matricial del modelo de regresión lineal simple ya conocido. Y por lo tanto sabemos que un estimador del vector de parámetros \(\mu\) es \(\hat\mu=y_{i.}\)
Asimismo podriamos expresar \(\mu_{i}\) como la suma de dos términos \(\mu\), común a todas las observaciones, y \(\tau_{i}\) específica a cada nivel, entonces el modelo propuesto sería:
\[y_{ij}=\mu+\tau_i+\epsilon_{ij},\ i=1,...,k;\ j=1,...,n_i\]
donde
\(y_{ij}\) es la variable aleatoria que representa la observación j-ésima del i-ésimo tratamiento (nivel i-ésimo del factor).
\(\mu\) es un efecto constante, común a todos los niveles, denominado media global.
\(\tau_i\) es la parte de \(y_{ij}\) debida a la acción del nivel i-ésimo, que será común a todos los elementos sometidos a ese nivel del factor. Los efectos de los factores \(\tau_{i}\) son las desviaciones de la media de cada nivel con respecto a la media general, por esta razón se debe verificar la relación \(\sum_{i=1}^k\tau_{i}=0\).
\(\epsilon_{ij}\) son variables aleatorias que engloban un conjunto de factores, cada uno de los cuales influye en la respuesta sólo en pequeña magnitud pero que de forma conjunta debe tenerse en cuenta en la especificación y tratamiento del modelo; es decir, las perturbaciones o error experimental pueden interpretarse como las variaciones causadas por todos los factores no analizados y que dentro del mismo tratamiento variarán de unos elementos a otros.
Para este modelo lineal supondremos que los errores \(\epsilon_{ij}\):
Tienen media cero
Tienen varianza constante (hipótesis de homocedasticidad)
Son estadísticamente independientes entre sí
Su distribución es normal
Analogamente al modelo anterior, se trata de contrastar si todos los niveles del factor producen el mismo efecto versus al menos un nivel del factor produce un efecto distinto.
En ambos casos si la hipótesis de que todos los nivels producen el mismo efecto es cierta, entonces todos los niveles tienen la misma media, \(\mu\), y la pertenencia a un grupo u otro es irrelevante, y podemos considerar todas las observaciones como provenientes de una única población.
5.2.1 Estimación de parámetros
La hipótesis de normalidad sobre los términos de error conlleva el hecho de que las variables \(y_{ij}\) sean normales e independientes, por lo que es inmediato construir la función de verosimilitud asociada a la muestra \(y=(y_{11},...,y_{1n_1},y_{21},...,y_{2n_2},...,y_{k1},...,y_{kn_k})\)
\[L(\mu,\tau_i,\sigma^2)=(2\pi\sigma^2)^{n/2}exp \left( -\frac{1}{2\sigma^2}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\mu-\tau_i]^2 \right)\]
Los estimadores máximo-verosímiles para los parámetros \(\mu,\tau_{i},\sigma^2\) son los valores para los cuales la función de verosimilitud alcanza su máximo. Por conveniencia, en vez de maximizar la función de verosimilitud, se maximiza el logaritmo de esta. En este caso,
\[ln(L(\mu,\tau_i,\sigma^2))=-\frac{n}{2}ln(2\pi)-\frac{n}{2}ln(\sigma^2)-\frac{1} {2\sigma^2}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\mu-\tau_i]^2\]
Las derivadas parciales respecto de los parámetros del modelo son:
\[\frac{\partial lnL}{\partial \mu}= \frac{1} {2}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\mu-\tau_i] \]
\[\frac{\partial lnL}{\partial \tau_i}= \frac{1} {\sigma^2}\sum_{j=1}^{n_i}[y_{ij}-\mu-\tau_i] \quad i=1,...,k \]
\[\frac{\partial lnL}{\partial \sigma^2}= \frac{N}{2\sigma^2}+\frac{1} {2(\sigma^2)^2}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\mu-\tau_i]^2\]
Igualando a cero y despejando los parámteros obtenemos los siguientes estimadores máximo verosimiles:
\[\hat\mu=\bar y_{..}\]
\[\hat\tau_i=\bar y_{i.}-\bar y{..}\]
\[\hat\sigma^2=\frac{1}{n}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\hat\mu-\hat\tau_i]^2=\frac{1}{n}\sum_{i=1}^{k}\sum_{j=1}^{n_i}[y_{ij}-\bar y_{i.}]^2=\frac{1}{n}\sum_{i=1}^{k}n_is_i^2\]
donde \(s_i^2\) es la varianza muestral del i-ésimo nivel.
5.2.2 Residuos
Los residuos se definen como las diferencias entre los valores observados \(y_{ij}\) y los valores previstos por el modelo \(\hat y_{ij}\) y los denotamos por
\[e_{ij}=y_{ij}-\hat y_{ij}=y_{ij}-\hat\mu-\hat \tau_{i}=y_{ij}-\bar y_{i.}\]
Los residuos son los estimadores de los errores aleatorios \(\epsilon_{ij}=y_{ij}-\mu-\tau_i\), los cuales son variables aleatorias no observables. Y es posible verificar que la suma de los residuos es cero.
Por otro lado con esta definición de los residuos, tenemos que el estimador máximo-verosimil \(\hat \sigma^2\) se puede escribir en términos de la suma de los residuos.
\[\hat \sigma^2=\frac{\sum_{ij}e_{ij}^2}{n}\]
5.2.3 Propiedades de los estimadores máximo verosímiles
Recordaremos algunas de las propiedades de los estimadores máximo verosímiles de \(\mu, \tau_{i}, \sigma^2\), estas propiedades son de utilidad para la descomposición de la varianza y las pruebas de hipótesis que se verán mas adelante.
5.2.3.1 Propiedades de \(\hat \mu\)
\(\hat\mu\) es un estimador insesgado de \(\mu\)
La varianza de \(\hat\mu\) es \(\frac {\sigma^2}{n}\)
\(\hat\mu\) tiene distribución normal, puesto que dicho estimador es combinación lineal de las variables \(y_{ij}\) y éstas son variables aleatorias independientes con distribución normal.
5.2.3.2 Propiedades de \(\hat \tau_{i}\)
\(\hat\tau_{i}\) es un estimador insesgado de \(\tau_{i}\)
La varianza de \(\hat\tau_{i}\) es \((n-n_i)\frac {\sigma^2}{nn_i}\)
\(\hat\tau_{i}\) tiene distribución normal, puesto que dicho estimador es combinación lineal de las variables \(y_{ij}\) y éstas son variables aleatorias independientes con distribución normal.
5.2.3.3 Propiedades de \(\hat \sigma^2\)
\(\hat \sigma^2\) NO es un estimador insesgado de \(\sigma^2\)
\(\frac {n\hat \sigma^2}{\sigma^2}\sim\chi^2_{(n-k)}\)
5.2.4 Descomposición de la variabilidad
Para comparar los efectos de los distintos niveles de un factor se emplea la técnica estadística denominada análisis de la varianza, abreviadamente ANOVA, que está basada en la descomposición de la variabilidad total de los datos en distintas componentes.
Se considera la siguiente identidad:
\[y_{ij}-\overline{y}=(y_{i.}-\overline{y})+(y_{ij}-y_{i.})\]
Elevando al cuadrado y sumando se obtiene:
\[\sum_{ij}(y_{ij}-\overline{y})^2=\sum_{ij}(y_{i.}-\overline{y})^2+\sum_{ij}(y_{ij}-y_{i.})^2+2\sum_{ij}(y_{i.}-\overline{y})(y_{ij}-y_{i.})\]
Considerando que
\[\sum_{ij}(y_{i.}-\overline{y})(y_{ij}-y_{i.})=\sum_{ij}(y_{ij}-y_{i.})y_{i.}-\sum_{ij}(y_{ij}-y_{i.})\overline{y}=0\]
La igualdad se da por que \(\hat y_{ij}=\hat\mu_i=y_{i.}\) entonces \(\sum_{ij}(y_{ij}-y_{i.})\) representa la suma de los residuales la cual es cero.
Por lo tanto la ecuación queda como:
\[\sum_{ij}(y_{ij}-\overline{y})^2=\sum_{ij}(y_{i.}-\overline{y})^2+\sum_{ij}(y_{ij}-y_{i.})^2\]
Generalmente se hace un cambio de notación:
\(SC_{T}\) es la suma de cuadrados total, mide la variabilidad de las observaciones del total respecto a la media global y es denotado por \(SC_{T}=\sum_{ij}(y_{ij}-\overline{y})^2\)
\(SC_{E}\) es la suma de cuadrados entre grupos, es denotado por \(SC_{E}=\sum_{ij}(y_{i.}-\overline{y})^2=\sum_{j}n_{i}(y_{i.}-\overline{y})^2\)
\(SC_{D}\) es la suma de cuadrados dentro de grupos o intragrupos, es denotado por \(SC_{D}=\sum_{ij}(y_{ij}-y_{i.})^2\)
De esta forma se tiene la siguiente igualdad:
\[SC_{T}=SC_{E}+SC_{D}\]
Además se observa que:
\(SC_{T}\) tiene \(n-1\) grados de libertad, ya que la suma \(\sum_{ij}(y_{ij}-\overline{y})^2\) es la parte central del estimador de la varianza \((n-1)S^2=\frac{\sum_{ij}(y_{ij}-\overline{y})^2}{n-1}\) en el cual al estimar media global con \(\overline{y}\) pierde un grado de libertad.
\(SC_{E}\) tiene \(k-1\) grado de libertad,ya que \(SC_{E}=\sum_{j}n_{i}(y_{i.}-\overline{y})^2\) la suma es sobre los \(k\) factores y se esta estimando la media global con \(\overline{y}\).
\(SC_{D}\) tiene \(n-k\) grados de libertad, ya que la regresión \(SC_{D}=SC_{T}-SC_{E},\) en el cual se ha visto que \(SC_{T}\) y \(SC_{E}\) tienen \(n-1\) y \(k-1\) grado de libertad, respectivamente, así los grados de libertad de \(SC_{D}=n-1-(k-1),\) entonces \(SC_{D}\) tienen n-k grados de libertad.
Sabemos que el modelo propuesto considera que los errores tienen distribución normal, por lo tanto dado que \(\epsilon_{ij} \sim N(0,\sigma^2)\):
\(\frac{SC_{D}}{\sigma^2}\sim \chi^2_{(n-k)}.\)
Bajo el supuesto de que las medias de las \(k\) poblaciones son iguales entonces \(SC_{E}/(k-1)\) es otro estimador insesgado de \(\sigma^2\) y \(\frac{SC_{E}}{\sigma^2}\sim\chi^2_{(k-1)}\)
Debido a un resultado de probabilidad se sabe que si \(x \sim \chi^2_{(n)}\) y \(y \sim \chi^2_{(m)}\) y si \(x\) es independiente a \(y\) entonces:
\[\frac{x/n}{y/m}\sim F_{(n,m)}\]
El contraste estadístico de interés en este modelo, como mencionamos al principio de esta sección, es el que tiene como hipótesis nula la igualdad de medias de los tratamientos, la cual se puede expresar de la siguiente forma:
\[H_0: \mu_1=\mu_2=...=\mu_k=\mu\]
o
\[H_0: \tau_1=\tau_2=...=\tau_k=0\]
Por lo tanto se puede aplicar la prueba \(F\) de Fisher en el análisis de varianza para probar las hipótesis de que las medias de las \(k\) poblaciones son iguales, quedando el estadístico F
\[F=\frac{\frac{SC_{E}}{k-1}}{\frac{SC_{D}}{n-k}} \sim F_{(k-1,n-k)}\]
El siguiente cuadro resume la información anterior, mejor conocida como Tabla ANOVA:
\[ \begin{array}{|c| c| c| c| c|} \hline &Grados\ de\ libertad & Suma\ de\ Cuadrados & Cuadrado\ Medio & Prueba\ F \\ \hline \hline Entre\ grupos & k-1 & SC_{E}=\sum_{j}n_{i}(y_{i.}-\overline{y})^2 & \frac{SC_{E}}{k-1} &\frac{SC_{E}/(k-1)}{SC_{D}/(n-k)} \\ \hline Dentro\ grupos & n-k & SC_{D}=\sum_{ij}(y_{ij}-y_{i.})^2& \frac{SC_{D}}{n-k} & -\\ \hline Total & n-1 & SC_{T}=\sum_{ij}(y_{ij}-\overline{y})^2 & - & - \\ \hline \end{array} \]
En esta tabla tenemos entonces la descomposición de la varianza observada en los datos. En la columna de suma de cuadrados tenemos la suma de cuadrados totales que y la descomposición de esta en la diferencia observada entre los grupos y la diferencia observada dentro de los grupos. En la tabla ANOVA también se incluye el estadístico de prueba F que se usa para concluir estadisticamente si los distintos niveles o grupos en los que están clasificadas las observaciones son significativamente diferentes o no.
Similar a lo estudiado en los modelos de regresión lineal, en este análisis de varianza también se tiene el coeficiente de determinación, el cual con frecuencia se le asocia a la proporción de la variación explicada por los grupos, ya que \(0<SC_{E}<SC_{T}\) entonces los valores del coeficiente de determinación están entre \(0<R^2<1.\)
Se define el coeficiente de determinación del modelo del como:
\[R^2=\frac{SC_{E}}{SC_{T}}=1-\frac{SC_{D}}{SC_{T}}.\]
5.2.5 Ejemplo
Se desean comparar dos medicamentos D (diurético), B (betabloquedor) con un producto inocuo P (placebo). Se tomó una muestra de 15 individuos hipertensos cuyas condiciones iniciales eran suficientemente homogéneas y se asignaron los tratamientos al azar.
El objetivo del estudio es ver cómo actuan los tres tratamientos frente a la hipertensión, concretamente si la disminuyen. Para cubrir con el objetivo, se ha elegido la variable observable “porcentaje de descenso de la presión arterial media.” Los datos obtenidos se presentan en la siguiente tabla:
| D | B | P |
|---|---|---|
| 22 | 20 | 10 |
| 28 | 28 | 5 |
| 30 | 35 | 0 |
| 15 | 19 | 14 |
| 17 | 33 | 18 |
Primero calculamos las medias por medicamento \(y_{D.}\), \(y_{B.}\) y \(y_{P.}\) asi como la media global \(y_{..}\)
x=c(22,18,30,15,17,20,28,35,19,33,10,5,0,14,18)
y=c(rep("D",5),rep("B",5),rep("P",5))
A=data.frame(x,y)
aggregate(A$x,by=list(A$y), FUN=mean)## Group.1 x
## 1 B 27.0
## 2 D 20.4
## 3 P 9.4mean(A$x)## [1] 18.93333\[y_{B.}=27.00\]
\[y_{D.}=20.40\]
\[y_{P.}=9.40\]
\[y_{..}=18.93\]
Ahora calculamos \(SC_E\) y \(SC_T\)
SCE=5*(27-18.93)^2+5*(20.4-18.93)^2+5*(9.4-18.93)^2
SCE## [1] 790.5335SCT=sum((x-mean(x))^2)
SCT## [1] 1348.933Por lo tanto \(SC_D\) se puede calcular como la diferencia entre \(SC_T\) y \(SC_E\)
SCD=SCT-SCE
SCD## [1] 558.3998La estadística \(F\) tomará entonces el valor de
Est_F=(SCE/2)/(SCD/12)
Est_F## [1] 8.494274Y su correspondiente p-value es
1-pf(Est_F,2,12)## [1] 0.005031862Aplicando las funciones de R que calculan el análisis de varianza obtenemos
anova(lm(x~y,data=A))## Analysis of Variance Table
##
## Response: x
## Df Sum Sq Mean Sq F value Pr(>F)
## y 2 790.53 395.27 8.4943 0.005032 **
## Residuals 12 558.40 46.53
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1mod <- aov(x ~ y, data = A)
mod## Call:
## aov(formula = x ~ y, data = A)
##
## Terms:
## y Residuals
## Sum of Squares 790.5333 558.4000
## Deg. of Freedom 2 12
##
## Residual standard error: 6.821535
## Estimated effects may be unbalancedDe los resultados es posible concluir que la diferencia entre los tres fármacos es claramente significativa.
5.3 Modelo de dos vías
Una variable o factor cuyo efecto sobre la respuesta no es directamente de interés pero que se introduce en el experimento para obtener comparaciones homogéneas se denomina una variable bloque. Por ejemplo, en una investigación para comparar la efectividad de varios fertilizantes (tratamientos) se puede considerar las fincas donde se prueban como un factor bloque.
El efecto de la finca sobre la producción no es de interés y el objetivo es comparar los fertilizantes eliminando el efecto de la pertenencia de una cosecha a una finca en específico. En general, se supone que no hay interacción entre la variable bloque y los factores de interés.
Supongamos que la variable respuesta está afectada por dos causas de variabilidad, es decir, por dos variables o factores cualitativos A y B, con a y b niveles respectivamente. El factor A es el factor principal, mientras que el factor B es una variable bloque. Supongamos también que tenemos unicamente una observación combinación de niveles. Eso significa tener tantas unidades experimentales por bloque como tratamientos o niveles del factor principal y que la asignación del tratamiento se hace al azar en cada bloque. Entonces, podemos representar las observaciones del siguiente modo:
| \(B_1\) | \(B_2\) | \(\dots\) | \(B_b\) | ||
|---|---|---|---|---|---|
| \(A_1\) | \(y_{11}\) | \(y_{12}\) | \(\dots\) | \(y_{1b}\) | \(y_{1.}\) |
| \(A_2\) | \(y_{21}\) | \(y_{22}\) | \(\dots\) | \(y_{2b}\) | \(y_{2.}\) |
| \(\dots\) | \(\vdots\) | ||||
| \(A_a\) | \(y_{a1}\) | \(y_{a2}\) | \(\dots\) | \(y_{ab}\) | \(y_{a.}\) |
| \(y_{.1}\) | \(y_{.2}\) | \(\dots\) | \(y_{.b}\) | \(y_{..}\) |
Con estos datos se pueden calcular las siguientes medias:
Media en el factor o nivel \(i\):
\[y_{i.} =\frac{1}{b}\sum_{j=1}^{b}y_{ij}\]
Media en el bloque \(j\):
\[y_{.j} =\frac{1}{a}\sum_{i=1}^{a}y_{ij}\]
Media general:
\[\bar{y} = y_{\cdot\cdot} =\frac{1}{ab}\sum_{i=1}^{a}\sum_{j=1}^{b}y_{ij}\]
Con esta presentación de los datos, se puede decir que A es el factor principal, mientras que B es una variable bloque. O también se puede decir que A es el factor fila y B el factor columna.
Realizando un proceso similar al que ANOVA de un factor lo que buscamos es plantear un modelo lineal que ajuste a los datos con las caracteristicas antes descritas. El primero que presentamos en un modelo en que suponemos que tanto el efecto fila como el efecto columna son aditivos.
5.3.1 Modelo aditivo
En el modeo aditivo de dos factores se asume que los efectos se suman, por lo que el modelo matemático queda de la siguiente forma:
\[y_{ij}=\mu+\alpha_i+\beta_j+\epsilon_{ij},\quad i=1,...,a;\quad j=1,...,b\]
donde:
\(\mu\) es la media general.
\(\alpha_i\) es el efecto del i-ésimo nivel del factor A.
\(\beta_j\) es el efecto del j-ésimo nivel del factor B.
Y en este caso pondremos la restricción de que la suma de los \(i\) parámetros \(\alpha\) y los \(j\) parámetros \(\beta\) sean cero. Es decir \[\sum_{i}\alpha_i=\sum_{j}\beta_{j}=0\]
Dicha restricción nos lleva a que el modelo de \(y_{ij}\) solamente depende de \(\mu,\alpha_1,...\alpha_{a-1},\beta_1,...\beta_{b-1}\) ya que \[\alpha_a=-\alpha_1-\alpha_2-...-\alpha_{a-1}\]
y
\[\beta_b=-\beta_1-\beta_2-...-\beta_{b-1}\]
Por lo tanto el numero de parámetros a estimar en este modelo serán a+b+1.
El objetivo es estimar los efectos de los distintos \(a\) niveles del factor A y contrastar la hipótesis de que todos los niveles del factor producen el mismo efecto. Analogamente se buscará que los distintos \(b\) niveles del factor B producen el mismo efecto. Con lo que se estaría concluyendo que el factor A y/o el factor B no son significativos.
Para la estimación de los parámetros de este modelo se utiliza el método de mínimos cuadrados, entonces partiendo de la ecuación \(y_{ij}-\mu-\alpha_i-\beta_j\) se suman ceros y queda la siguiente igualdad.
\[y_{ij}-\mu-\alpha_i-\beta_j=(\bar y-\mu)+(y_{i.}-\bar y-\alpha_i)+(y_{.j}-\bar y-\beta_j)+(y_{ij}-y_{i.}-y_{.j}+\bar y)\]
la cual al elevar al cuadrado y sumar sobre todos los valores obtenemos del lado izquierdo la función a minimizar y del lado derecho después de desarrollar los cuadrados y observar que los productos cruzados se anulan queda:
\[\sum_{ij}(y_{ij}-\mu-\alpha_i-\beta_j)^2=\sum_{ij}(\bar y-\mu)^2+\sum_{ij}(y_{i.}-\bar y-\alpha_i)^2+\sum_{ij}(y_{.j}-\bar y-\beta_j)^2+\sum_{ij}(y_{ij}-y_{i.}-y_{.j}+\bar y)^2\]
por lo tanto, considerando la restricción de que la suma de los parámetros \(\alpha\) y \(\beta\) deben ser cero llegamos a que los valores de \(\mu, \alpha_i, \beta_i\) que minimizan la ecuación serán.
\[\hat \mu= \bar y\]
\[\hat\alpha_i=y_{i.}-\bar y\]
\[\hat\beta_j= y_{.j}-\bar y\]
Entonces la suma de los cuadrados de los residuales estará dado por \(SC_R=\sum_{ij}(y_{ij}-y_{i.}-y_{.j}+\bar y)^2\)
5.3.1.1 Hipótesis
En este caso la hipótesis se plantean para cada factor. Siendo las hipótesis nulas las que dicen que el factor A o B según el caso no muestran efecto.
Entonces la hipótesis de que el factor principal A no es significativo (no hay efecto fila) es:
\[H_0^A: \alpha_1=\alpha_2=...=\alpha_a=0\] Y la hipótesis de que el factor B no es significativo (no hay efecto columna) es:
\[H_0^B: \beta_1=\beta_2=...=\beta_b=0\]
La prueba para contrastar las hipótesis será una F y se obtendrá de forma similar al caso de una vía con la descomposición de la suma de los cuadrados.
\[\sum_{ij}(y_{ij}-\overline{y})^2=b\sum_{ij}(y_{i.}-\overline{y})^2+a\sum_{ij}(y_{.j}-\overline{y})^2-\sum_{ij}(y_{ij}-y_{i.}-y_{.j}+\overline{y})^2\]
\[SC_T=SC_F+SC_C+SCR\]
donde \(SC_T\) es la suma de cuadrados total, \(SC_F\) es la suma de cuadrados entre filas, \(SC_C\) es la suma de cuadrados entre columnas y \(SCR\) es la suma de cuadrados residual.
Se observa que si \(H_0^B\) es cierta entonces el modelo se convertirá en un modelo de una vía \[y_{ij}=\mu +\alpha_i+\epsilon_{ij}\] Y la suma de cuadrados residual será entonces
\[SCR_H=\sum_{ij}(y_{ij}-y_{i.})^2=\sum_{ij}(y_{.j}+\overline{y})^2\sum_{ij}(y_{ij}-y_{i.}-y_{.j}+\overline{y})^2=SC_C+SCR\]
Entonces la regla de decisión para la prueba de hipótesis esta dada por el estadístico
\[F=\frac{\frac{SC_C}{b-1}}{\frac{SCR}{(a-1)(b-1)}}\]
Bajo \(H_0\) el estadítico tiene distribución \(F\) con \((a-1)\) y \((b-1)\) grados de libertad
5.3.1.2 Descomposición de la variabilidad
Para comparar los efectos de los distintos niveles de ambos factores emplearemos el análisis de la varianza (ANOVA). Similar a el caso de un factor, el análisis se basa en la descomposición de la variabilidad total de los datos en distintas componentes.
| Grados de libertad | Suma de Cuadrados | Cuadrado Medio | Prueba F | |
|---|---|---|---|---|
| Entre filas | \(a-1\) | \(SC_{F}=b\sum_{i}(y_{i.}-\overline{y})^2\) | \(\frac{SC_{F}}{a-1}\) | \(\frac{SC_{F}/(a-1)}{SCR/[(a-1)(b-1)]}\) |
| Entre columnas | \(b-1\) | \(SC_{C}=a\sum_{j}(y_{.j}-\overline{y})^2\) | \(\frac{SC_{C}}{b-1}\) | \(\frac{SC_{C}/(b-1)}{SCR/[(a-1)(b-1)]}\) |
| Residuo | \((a-1)(b-1)\) | \(SCR=\sum_{ij}(y_{ij}-y_{i.}-y{.j}+\overline{y})^2\) | \(\frac{SCR}{(a-1)(b-1)}\) | - |
| Total | \(ab-1\) | \(SC_{T}=\sum_{ij}(y_{ij}-\overline{y})^2\) | - | - |
Cuando el efecto de la variable bloque no es significativo se puede considerar el modelo más sencillo con un solo factor, prescindiendo de los bloques. Sin embargo, si hay diferencias entre los bloques, el modelo en bloques aleatorizados es mucho más eficaz en la detección de diferencias entre tratamientos.
5.3.1.3 Coeficientes de determinación parcial
El coeficiente de determinación se define como
\[R^2=1-\frac{SCR}{SC_T}=\frac{SC_F+SC_C}{SC_T}\]
De modo que los coeficientes de determinación parcial
\[R^2_F=\frac{SC_F}{SC_T}\]
y \[R^2_C=\frac{SC_C}{SC_T}\]
indican el porcentaje de la variabilidad total explicada por el factor principal y por el factor bloque, respectivamente.
5.3.1.4 Ejemplo
Retomando el ejemplo mencionado en la introducción, se desea estudiar las diferencias entre los efectos de 4 fertilizantes sobre la producción de papas, para ello se cuenta con 5 fincas, cada una de las cuales se dividió en 4 parcelas del mismo tamaño y tipo. Los fertilizantes fueron asignados al azar en las parcelas de cada finca. El rendimiento en toneladas fue
| Finca | |||||
|---|---|---|---|---|---|
| Fertilizante | 1 | 2 | 3 | 4 | 5 |
| 1 | 2.1 | 2.2 | 1.8 | 2.0 | 1.9 |
| 2 | 2.2 | 2.6 | 2.7 | 2.5 | 2.8 |
| 3 | 1.8 | 1.9 | 1.6 | 2.0 | 1.9 |
| 4 | 2.1 | 2.0 | 2.2 | 2.4 | 2.1 |
Este diseño es especialmente utilizado en experimentación agrícola. Y el objetivo es comparar \(a = 4\) tratamientos (fertilizantes en este caso) utilizando \(b = 5\) bloques (fincas) y repartiendo aleatoriamente los a tratamientos en cada uno de los bloques (los fertilizantes son asignados al azar en las parcelas de cada finca). Para una correcta aplicación de este diseño se debe considerar que hay homogeneidad dentro de cada bloque, de modo que el efecto bloque sea el mismo para todos los tratamientos.
Nos interesa saber si hay diferencias significativas entre los tratamientos \(\alpha_i\) y entre los bloques \(\beta_j\). Los resultados obtenidos son
Primero calculamos las medias por tratamiento (\(y_{1.}\), \(y_{2.}\), \(y_{3.}\) y \(y_{4.}\)), las medias por finca, asi como la media global \(y_{..}\)
x=c(2.1,2.2,1.8,2.0,1.9,2.2,2.6,2.7,2.5,2.8,1.8,1.9,1.6,2.0,1.9,2.1,2.0,2.2,2.4,2.1)
f=c(rep("1",5),rep("2",5),rep("3",5),rep("4",5))
c=c(rep(c("1","2","3","4","5"),4))
A=data.frame(x,f,c)
MF=aggregate(A$x,by=list(A$f), FUN=mean)
MC=aggregate(A$x,by=list(A$c), FUN=mean)
MG=mean(A$x)
MF## Group.1 x
## 1 1 2.00
## 2 2 2.56
## 3 3 1.84
## 4 4 2.16MC## Group.1 x
## 1 1 2.050
## 2 2 2.175
## 3 3 2.075
## 4 4 2.225
## 5 5 2.175MG## [1] 2.14Ahora calculamos \(SC_F\), \(SC_C\) y \(SC_T\)
SCF=5*sum((MF$x-MG)^2)
SCF## [1] 1.432SCC=4*sum((MC$x-MG)^2)
SCC## [1] 0.088SCT=sum((x-mean(x))^2)
SCT## [1] 1.928SCR=SCT-SCF-SCCLa estadística \(F\) para comparar los tratamientos tomará entonces el valor de
Est_F=(SCF/3)/(SCR/12)
Est_F## [1] 14.03922Y su correspondiente p-value es
1-pf(Est_F,3,12)## [1] 0.0003136709La estadística \(F\) para comparar las fincas tomará entonces el valor de
Est_F=(SCC/4)/(SCR/12)
Est_F## [1] 0.6470588Y su correspondiente p-value es
1-pf(Est_F,3,12)## [1] 0.599646Por lo tanto se concluye que si hay un efecto por tratamiento pero no por finca.
Aplicando las funciones de R que calculan el análisis de varianza obtenemos
mod <- aov(x ~ f+c, data = A)
summary(mod)## Df Sum Sq Mean Sq F value Pr(>F)
## f 3 1.432 0.4773 14.039 0.000314 ***
## c 4 0.088 0.0220 0.647 0.639572
## Residuals 12 0.408 0.0340
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1De los resultados es posible concluir que el efecto por fila si es significativo pero el efecto por columna no es significativo, lo cual en términos del problema se traduce a que si hay evidencia de una diferencia entre los fertiizantes.
5.3.2 Modelo aditivo con interacción
Este modelo es una extensión del modeo aditivo en donde adicional a suponer que los efectos de los dos factores se suman, también se considera la presencia del efecto interacción por lo que el modelo matemático queda de la siguiente forma:
\[y_{ijk}=\mu+\alpha_i+\beta_j+\gamma_{ij}+\epsilon_{ijk},\quad i=1,...,a;\quad j=1,...,b; \quad k=1,...,r\]
donde:
\(\mu\) es la media general.
\(\alpha_i\) es el efecto del i-ésimo nivel del factor A.
\(\beta_j\) es el efecto del j-ésimo nivel del factor B.
\(\gamma_{ij}\) es la interacción entre el i-ésimo nivel del factor A y el j-ésimo nivel del factor B.
Análogo a lo visto en el modelo anterior, para la determinación de los parámetros se imponen las siguientes restricciones:
\[\sum_{i}\alpha_i=\sum_{j}\beta_{j}=\sum_{i}\gamma_{ij}=\sum_{j}\gamma_{ij}=0\]
Por otra parte, se dice que un diseño es de rango completo si el número de parámetros es igual al número de condiciones experimentales, es decir, al número de filas distintas de la matriz de diseño. En un diseño que no es de rango completo hay menos parámetros que condiciones experimentales, por lo que en realidad “admitimos” que los datos se ajustan al modelo propuesto.
En general, un modelo de rango completo se ajusta intrínsecamente a los datos sin problemas. No obstante, para poder estimar todos los parámetros es necesario disponer de más de una réplica por condición experimental.
Por esta razón en ocasiones la interacción no puede ser incluida en un modelo aditivo como la parametrización presentada.Sin embargo para resolver este problema se pueden considera las siguientes transformaciones:
\[\mu=\frac{1}{ab}\sum_{ij}\eta_{ij}\]
\[\alpha_i=\frac{1}{b}\sum_j\eta_{ij}-\mu\]
\[\beta_j=\frac{1}{a}\sum_i\eta_{ij}-\mu\]
\[\gamma_{ij}=\eta_{ij}-\mu-\alpha_i-\beta_j\]
Y con esto el modelo puede ser reparametrizado de la siguiente forma:
\[y_{ijk} = \eta_{ij}+\epsilon_{ijk}\]
Los estimadores de los parámetros serán
\[\hat \mu= \bar y\]
\[\hat\alpha_i=y_{i..}-\bar y\]
\[\hat\beta_j= y_{.j.}-\bar y\]
\[\hat\gamma_{ij}=y_{ij.}-y_{i..}-y_{.j.}+\overline{y}\]
Entonces la suma de los cuadrados de los residuales estará dado por \(SC_R=\sum_{ij}(y_{ijk}-y_{ij.})^2\). Entonces \(\hat\eta_{ij}=y_{ij.}\)
5.3.2.1 Hipótesis y tabla ANOVA
En este caso las hipótesis son las siguientes:
\[H_0^A: \alpha_1=\alpha_2=...=\alpha_a=0 \quad (no\ hay\ efecto\ fila)\]
\[H_0^B: \beta_1=\beta_2=...=\beta_b=0 \quad (no\ hay\ efecto\ columna)\]
\[H_0^I: \gamma_{ij}=0 \quad\forall i,j \quad (no\ hay\ interacción)\]
La prueba para contrastar las hipótesis será una F y se obtendrá de forma similar a los casos anteriores con la descomposición de la suma de los cuadrados.Entonces la suma de cuadrados totales \(SC_T\) se descompondrá en:
\[SC_T=\sum_{ijk} (y_{ijk}-\overline y)^2=br\sum_{ij}(y_{i..}-\overline{y})^2+ar\sum_{ij}(y_{.j.}-\overline{y})^2-r\sum_{ij}(y_{ijk}-y_{i..}-y_{.j.}+\overline{y})^2+\sum_{ijk}(y_{ijk}-y_{ij.})^2\]
\[SC_T=SC_F+SC_C+SC_I+SCR\]
La tabla ANOVA muestra la descomposición de la variabilidad total de los datos en los distintos componentes considerados en este modelo (factor fila, factor columna e interacción.
| Grados de libertad | Suma de Cuadrados | Cuadrado Medio | Prueba F | |
|---|---|---|---|---|
| Entre filas | \(a-1\) | \(SC_{F}=br\sum_{i}(y_{i..}-\overline{y})^2\) | \(\frac{SC_{F}}{a-1}\) | \(\frac{SC_{F}/(a-1)}{SCR/[(a-1)(b-1)]}\) |
| Entre columnas | \(b-1\) | \(SC_{C}=ar\sum_{j}(y_{.j.}-\overline{y})^2\) | \(\frac{SC_{C}}{b-1}\) | \(\frac{SC_{C}/(b-1)}{SCR/[(a-1)(b-1)]}\) |
| Interacción | \((a-1)(b-1)\) | \(SC_{I}=r\sum_{ij}(y_{ij.}-y_{i..}-y_{.j.}+\overline{y})^2\) | \(\frac{SC_{I}}{(a-1)(b-1)}\) | \(\frac{SC_{I}/[(a-1)(b-1)]}{SCR/[ab(r-1)]}\) |
| Residuo | \(ab(r-1)\) | \(SCR=\sum_{ijk}(y_{ijk}-y_{ij.})^2\) | \(\frac{SCR}{ab(r-1)}\) | - |
| Total | \(abr-1\) | \(SC_{T}=\sum_{ijk}(y_{ijk}-\overline{y})^2\) | - | - |
Similar a lo que hemos visto en los modelos de ANOVA previos entonces cada una de las hipótesis se probará con su estadística de prueba F que se incluye en la última columna de la tabla.