Inteligencia Analítica de Datos con R

Análisis de Varianza (ANOVA)

Msc. Roberto Trespalacios

Universidad Tecnológica de Bolivar

6/8/23

Análisis de varianza (ANOVA - ANalysis Of VAriance)

  • El Análisis de la Varianza (ANOVA - ANalysis Of VAriance) que constituye, sin lugar a dudas, una de las herramientas más valiosas de la Inferencia Estadística.

  • Desarrollado hacia 1930 por R.A. FISHER, el Anova es una de las técnicas básicas para el estudio de observaciones que dependen de varios factores.

  • Es una herramienta fundamental en el análisis de los modelos de Regresión Lineal y de Diseño de Experimentos.

Análisis de varianza (ANOVA - ANalysis Of VAriance)

Aspectos teóricos: ANOVA de un Factor

  • El ANOVA de un factor es el modelo básico de análisis de la varianza. No puede considerarse estrictamente un modelo multivariante, al involucrar únicamente a dos variables, la variable respuesta y el factor explicativo.

  • A continuación presentaremos los aspectos teóricos del Análisis de Varianza de un factor

Sean \[\begin{align*} Y = & \text{ variable respuesta, cuantitativa continua}\\ X= & \text{ variable explicativa, categórica con I categoría} \end{align*}\]

Aspectos teóricos: ANOVA de un Factor

Se dispone de una muestra de n observaciones de la variable Y, que podemos considerar estructuradas de la siguiente forma:

donde, en general:

  • \(y_{jh}\) representa la observación \(h\)-ésima (con \(h = 1,2,...,n_j\)) en la columna o categoría \(j\)
  • \(n_1, n_2, ..., n_I\) los tamaños de muestra respectivos en las \(I\) categorías de la variable \(X\).

Aspectos teóricos: ANOVA de un Factor

  • Una variable respuesta, \(Y\), de tipo numérico con observaciones
  • Una variable predictora, \(F\), de tipo categórico con \(I\) grupos o niveles distintos de tamaños muestrales \(n_1, n_2, \ldots , n_I\), de forma que \(n = n_1 + n_2 + \ldots + n_I\) y el vector de observaciones de la respuesta se puede escribir como:

\[\begin{align*} &y_{11},\dots,y_{1n_1}\\ &y_{21},\dots,y_{2n_2}\\ &y_{I1},\dots,y_{In_I} \end{align*}\]

donde el primer subíndice indica el nivel del factor y el segundo la posición dentro del conjunto de datos de dicho nivel del factor.

Descomposición de la variabilidad

  • El análisis de la varianza resolverá la identificación de la situación a través de la descomposición de la variabilidad de la variable respuesta en componentes que permitirán construir las pruebas de hipótesis adecuadas.
  • Estamos interesados en contrastar hipótesis para el efecto de un solo factor, que se presenta con \(I\) niveles, sobre la variable respuesta.
  • Para ello realizamos el siguiente contraste:

\[\begin{align*} H_0:& \mu_1=\mu_2=\dots=\mu_I\\ H_1:& \mu_i\neq\mu_j, \text{ para algún } i\neq j \end{align*}\]

Es decir, contrastamos dos hipótesis:

  • H0: no hay diferencia en las medias de los \(I\) tratamientos.
  • H1: de que al menos una media difiere de otra.

Descomposición de la variabilidad

  • Todo este planteamiento se puede formalizar de manera general para cualquier experimento unifactorial.

Contexto del problema

Supongamos un factor con I niveles y para el nivel i-ésimo se obtienen \(n_i\) observaciones de la variable respuesta. Entonces podemos postular el siguiente modelo

\[y_{ij} = \mu + \beta_i + e_{ij}, \qquad i = 1,\dots, I; j = 1, \dots, n_i\]

donde:

  • \(y_{ij}\): es la variable aleatoria que representa la observación \(j\)-ésima del \(i\)-ésimo tratamiento (Variable respuesta).

  • \(\mu\): Es un efecto constante, común a todos los niveles del factor, denominado media global.

\(\beta_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, llamado efecto del tratamiento \(i\)-ésimo.

  • \(e_{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.
  • \(e_{ij}\) Reciben el nombre de perturbaciones o error experimental.

Descomposición de la variabilidad

  • Nuestro objetivo es estimar el efecto de los tratamientos 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 entre sí.

  • En este modelo se distinguen dos situaciones según la selección de los tratamientos:

    • Modelo de efectos fijos.
    • Modelo de efectos aleatorios.

Modelo de efectos fijos

  • En el modelo de efectos fijos el experimentador decide qué niveles concretos se van a considerar
  • Las conclusiones que se obtengan sólo son aplicables a esos niveles, no pudiéndose hacer extensivas a otros niveles no incluidos en el estudio.

Modelo de efectos aleatorios

  • En el modelo de efectos aleatorios, los niveles del factor se seleccionan al azar; es decir los niveles estudiados son una muestra aleatoria de una población de niveles.
  • Las conclusiones que se obtengan se generalizan a todos los posibles niveles del factor, hayan sido explícitamente considerados en el estudio o no.

Supuestos del modelo y tamaños muestrales

Supuestos

El ANOVA parte de algunos supuestos o hipótesis que han de cumplirse:

  • Normalidad de los datos.
  • Independencia de las observaciones.
  • La distribución de los residuales debe ser normal.
  • Homocedasticidad: homogeneidad de las varianzas.

Tamaños muestrales de los tratamientos

En cuanto a los tamaños muestrales de los tratamientos, los modelos se clasifican en:

  • modelo equilibrado o balanceado si todas las muestras son del mismo tamaño \(n_i = n\).
  • modelo no-equilibrado o no-balanceado si los tamaños muestrales ni son distintos.

Descomposición de la variabilidad de la variable respuesta

Para nuestro modelo,

\[y_{ij} = \mu + \beta_i + e_{ij}, i = 1,\dots, I; j = 1, \dots, n_i\]

Consideremos los siguientes estadísticos muestrales obtenidos a partir de los datos:

\[\begin{align*} \bar{y}_{..}&=\frac{\sum_{ij} y_{ij}}{n} \text{, media global de la variable respuesta.}\\ \bar{y}_{i.}&=\frac{\sum_{j} y_{ij}}{n_i}\text{, media de la variable respuesta en la categoría } i \text{ del factor, } i=1,2,..., I \end{align*}\]

Contraste de hipótesis

El contraste de hipótesis planteado anteriormente, supongamos que tenemos los niveles \(A, B, C, D\) de un factor \(F\), la idea es probar que:

\[\begin{align*} H_0:& \mu_A=\mu_B=\mu_C=\mu_D, \\ H_1:& \mu_i\neq\mu_j, \text{ para algún } i\neq j \end{align*}\]

Descomposición de la variabilidad de la variable respuesta

En el contraste de hipótesis planteado anteriormente: \(H_0: \mu_1=\mu_2=\dots=\mu_I\) vs. $H_1: _i_j, i$, la variabilidad se descompone así:

\[SS_T = SS_{\beta} + SS_e\]

donde

  • \(SS_T\): es la suma total de cuadrados o variabilidad total de la respuesta \(Y\): \[\displaystyle SS_T=\sum_{i=1}^I\sum_{j=1}^{n_i} (y_{ij}-\bar{y}_{..})^2\]

  • \(SS_{\beta}\) es la suma de cuadrados entre tratamientos o variabilidad explicada: \[\displaystyle SS_{\beta}=\sum_{i=1}^In_i (\bar{y}_{i.}-\bar{y}_{..})^2\]

  • \(SS_e\) es la suma dentro de los tratamientos, variabilidad no explicada o residual: \[\displaystyle SS_e=\sum_{ij}(y_{ij}-\bar{y}_{i.})^2\]

Tabla ANOVA

  • Los grados de libertad representan el número de fuentes independientes de variación para esa medida de variabilidad.

  • Los cuadrados medios o medias cuadráticas representan las variabilidades promediadas por los grados de libertad que las producen.

  • El estadístico para el contraste formulado es:

\[F=\frac{MS_{\beta}}{MS_e} \sim F_{(I-1, n-I)}\]

Ejemplo 1: Aplicación a la administración

Se realizó un estudio en Silicon Valley, para comprobar si en condiciones determinadas, como el nivel de preparación académica, influye significativamente en el salario de un trabajador. Para ello se toman 5 individuos homogéneos (edad, salud, raza, etc) con niveles de preparación académica como sigue: 5 con licenciatura y 5 con maestría y 5 con doctorado. Sus salarios (miles de dólares) y se recogen en una tabla obteniéndose los siguientes datos:

Licenciatura Maestría Doctorado
180 172 163
173 158 170
175 167 158
182 160 162
181 175 170
157 162 182

Solución

Vemos que es un problema de un factor,con tres niveles y cada nivel con 5 replicas. El modelo es el siguiente:

  • \(y_{ij} = \mu + \beta_i + e_{ij}, i = 1,\dots, 3; j = 1, \dots, 5\)
  • Debe verificarse :

Encontremos cada término de la suma. Primero, la Suma de Cuadrados Total: \(SS_T\):

  • La media general es: \(\bar{y}_{..}= 169.73\)

por lo tanto, la Suma de Cuadrados Total es:

\[\begin{align*} \displaystyle SS_T&=\sum_{i=1}^I\sum_{j=1}^{n_i} (y_{ij}-\bar{y}_{..})^2 \\ & =(180-169.73)^2 +\cdots+ (181-169.73)^2 + \\ & =(172-169.73)^2 +\cdots+ (175-169.73)^2 +\\ & =(163-169.73)^2 +\cdots+ (170-169.73)^2 = 936.9 \end{align*}\]

Esta \(SS_T\) tiene un número de grados de libertad que es igual al número total de datos menos uno: \[gl_{total} = 15 - 1 = 14\]

Solución (continuación)

Encontremos Suma de Cuadrados del Factor \(SS_{\beta}\). Para esto, calculemos las medias para cada nivel del factor(Licenciatura =1, Maestría= 2 Doctorado =3).

\[\begin{align*} \bar{y}_{.1}& = 178.2\\ \bar{y}_{.2}& = 166.4\\ \bar{y}_{.3}& = 164.6 \end{align*}\]

por lo tanto, la suma de Cuadrados del Factor es:

\[\begin{align*} \displaystyle SS_{\beta}=&\sum_{i=1}^In_i (\bar{y}_{i.}-\bar{y}_{..})^2\\ =& 5(178.2-169.73)^2+5(166.4-169.73)^2+5(164.6-169.73)^2\\ =& 545.7 \end{align*}\]

  • Recordemos que son tres niveles del factor, con 5 repeticiones para cada nivel del factor.

  • \(SS_{\beta}\) tiene un número de grados de libertad que es igual al número de niveles del factor menos uno: \[gl_{factor} = 3 - 1 = 2\]

Solución (continuación)

Encontremos Suma de Cuadrados Residual(errores): \(SS_{e}\). Recordemos que las medias de cada nivel Licenciatura, Maestría, Doctorado, son respectivamente: \[\bar{y}_{.1} = 178.2, \qquad \bar{y}_{.2} = 166.4, \qquad \bar{y}_{.3} = 164.6\].

por lo tanto, realizando los cálculos, tenemos que:

\[\begin{align*} \displaystyle SS_e=\sum_{ij}(y_{ij}-\bar{y}_{i.})^2 \\ & =(180-178.2)^2 +\cdots+ (181-178.2)^2 + \\ & =(172-166.4)^2 +\cdots+ (175-166.4)^2 +\\ & =(163-164.6)^2 +\cdots+ (170-164.6)^2 = 391.2 \end{align*}\]

Los grados de libertad asociados a la \(SS_e\) se obtienen por diferencia entre los \(gl_{total} =14\) y los \(gl_{factor}= 2\), luego,

\[gl_{resid} = 15 - 3 = 12\]

## Solución (continuación)

Podemos verificar, que \[\begin{align*} SStotal &= SSfactor + SSerror\\ 936.9& = 545.7 + 391.2 \end{align*}\]

Por lo tanto, la tabla ANOVA es:

  • Los grados de libertad representan el número de fuentes independientes de variación para esa medida de variabilidad.
  • Los cuadrados medios o medias cuadráticas representan las variabilidades promediadas por los grados de libertad que las producen.

Solución (continuación)

De la tabla anterior, podemos observar que el estadístico para el contraste formulado es:

\[F=\frac{MS_{\beta}}{MS_e}=\frac{272.85}{27.94} \sim F_{(I-1, n-I)}=F_{(3-1, 15-3)}\]

  • Si vemos en la tabla F de Fisher para:
    • grados de libertad de numerador = 2,
    • grados de libertad de denominador = 27
    • nivel de significancia = 0.05

\[F= 9.76 > 3.8 = F_{(2, 13),0.05}\]

  • Por lo tanto, se rechaza la hipótesis nula de que el efecto del de la raza sobre la media de la estatura no es significativo,
  • Es decir, rechazamos que: \(\mu_A=\mu_B=\mu_C=\mu_D\)

Por lo tanto, podemos afirmar que existe al menos una diferencia significativas entre las medias de alguna(s) de las estaturas para cada raza.

Solución en R

  1. Lea los datos nivel_salario

  2. ¿Se verifican los supuestos del modelo de análisis de la varianza?

  3. ¿Existe evidencia de que el nivel de educación influye en el salario?

  4. ¿Una mayor educación es favorable?. ¿Qué grado de educación cree que tiene mayores salarios?

  5. Leemos los datos del computador.

nivel_salario <- read_csv("C:/Users/Roberto Trespalacio/Documents/sem_3.3/datos/nivel_salario3.csv")

veamos el encabezado de los datos

head(nivel_salario)

Estamos ante un análisis balanceado ya que todas las muestras tienen el mismo tamaño (hay 5 individuos de cada nivel de educación).

table(nivel_salario)
  1. ¿Se verifican los supuestos del modelo de análisis de la varianza?

Para poder aplicar esta técnica estadística se han de verificar las siguientes condiciones:

  • Normalidad en la distribución del salario por nivel de educación.
  • Homocedasticidad o igualdad de varianzas poblacionales del salario por nivel de educación.
  • Independencia de las observaciones, que se consigue extrayendo de las poblaciones muestras aleatorias independientes. Esta hipótesis suele ser asumida como consecuencia del diseño del experimento, pues los distintos niveles del factor son asignados al azar (en este caso concreto, los productos son elegidos aleatoriamente). Por ello, solo comprobaremos los otros dos supuestos.

Solución en R

Boxplot por nivel

Antes de verificar si se cumplen las condiciones necesarias para llevar a cabo un ANOVA, exploraremos los datos mediante un gráfico de caja:

boxplot(nivel_salario$salario ~ nivel_salario$educacion)  

El gráfico parece sugerir que no hay grandes desviaciones respecto a la normalidad aunque podrían existir diferencias en la variabilidad. En cualquier caso realizaremos tests de hipótesis para confirmarlo estadísticamente.

Normalidad

Mediante el contraste de normalidad de Shapiro-Wilk, se puede comprobar si una variable sigue una distribución normal. Se contrasta la hipótesis nula de que la distribución de la que proceden los datos es normal para cada nivel de educación:

Licenciatura

shapiro.test(nivel_salario$salario[nivel_salario$educacion == "licenciatura"])

Maestría

shapiro.test(nivel_salario$salario[nivel_salario$educacion == "maestria"])

Doctorado

shapiro.test(nivel_salario$salario[nivel_salario$educacion == "doctorado"])

También puede llegarse a este resultado mediante la función shapiro.test combinada con tapply. Ésta última función permite aplicar el test de Shapiro-Wilk sobre la variable estatura en cada uno de los niveles del factor raza:

tapply(nivel_salario$salario, nivel_salario$educacion, shapiro.test)

Por tanto, en cada uno de los grupos, es decir, para cada uno de los niveles de educación, no rechazamos la hipótesis nula de normalidad al 5% (al ser el p-valor superior a 0.05 en todos los casos). En la muestra no hay evidencia para rechazar que los datos en cada grupo provienen de una distribución normal.

Homocedasticidad

De manera descriptiva, una primera aproximación para estudiar cómo son las varianzas poblacionales de los grupos es calcular las desviaciones típicas muestrales del salario para cada nivel de edcuación. Para ello utilizamos la función sd combinada con tapply:

tapply(nivel_salario$salario, nivel_salario$educacion, sd)
  • Se observa cierta diferencia en la variabilidad en los salarios por grupos, pero para contrastar si esa diferencia es estadísticamente significativa emplearemos un test de homogeneidad de varianzas.
  • Uno de los test disponibles es el de Bartlett, que se puede aplicar solo cuando las observaciones provienen de poblaciones normales.
  • La hipótesis nula es que las varianzas son idénticas en todos los grupos y, por tanto, se cumple el requisito de homocedasticidad.
  • La función bartlett.test se puede usar de la siguiente forma, incluyendo en primer lugar la variable dependiente y después el factor.
bartlett.test(nivel_salario$salario, nivel_salario$educacion)
  • Como el p-valor es superior a 0.05, no rechazaríamos la hipótesis nula de igualdad de varianzas al 5%, y, por tanto, no hay evidencia en contra del requisito de homocedasticidad.

  • Si el factor únicamente tuviera dos niveles podríamos utilizar también la función var.test, que está indicada para el contraste de varianzas en dos muestras provenientes de poblaciones normales, como ya se vio en la anterior clase.

Contraste de hipótesis

  1. ¿Existe evidencia de que el nivel de educación influye en en el salario?

La hipótesis que se desea contrastar es si el salario es igual en todos los niveles de educación frente a la alternativa de que hay diferencias,

  • \(H_0: \mu_1=\mu_2=\mu_3\) todas las medias son iguales
  • \(H_1:\) al menos una de las medias difiere del resto

siendo \(\mu_i\) el salario medio con el método i, i=Licenciatura, Maestria, Doctorado.

  • Podriamos realizar todos los contrastes, pero al efectuar un mayor número de contrastes aumenta la probabilidad de error de tipo I del conjunto (la probabilidad de detectar una diferencia significativa en alguno de los test cuando no existen diferencias), algo que hay que evitar.
  • El diseño de un análisis de la varianza (siempre que se cumplan las condiciones para ello) permite eludir este problema.

Contraste de hipótesis

Otra forma de escribirlo es:

\(H_0: \mu_{Licenciatura}=\mu_{Maestria}\), \(H_0: \mu_{Licenciatura}=\mu_{Doctorado}\) y \(H0: \mu_{Maestria}=\mu_{Doctorado}\)

  • Para realizar un ANOVA, R ofrece varias funciones, siendo la más sencilla oneway.test.
  • Esta función permite llevar a cabo el contraste F del ANOVA de un factor cuando se verifica el supuesto de homoscedasticidad (especificando el argumento var.equal = TRUE)
  • También el contraste de Welch en caso contrario (con var.equal = FALSE, que es la opción por defecto). Como en los datos no se han detectado evidencias contra la igualdad de varianzas,
oneway.test(nivel_salario$salario ~ nivel_salario$educacion, var.equal = TRUE)
# También oneway.test(salario ~ educacion, var.equal = TRUE, nivel_salario)

Al ser el p-valor inferior a 0.01, se puede concluir que hay diferencias en el salario medio según la educación (a todos los niveles de significación usuales).

La misma conclusión se obtiene con la función aov,

aov(nivel_salario$salario ~ nivel_salario$educacion)
# También aov(salario ~ educacion, nivel_salario)

Contraste de hipótesis

  • Aunque no proporciona la tabla completa del test ANOVA, en particular el p-valor, si no se usa seguidamente la función summary.
  • Así, repetimos el comando aov, pero ahora almacenando su resultado en el objeto ANOVA1, para después aplicar a dicho objeto la función summary:
ANOVA1 <- aov(salario ~ educacion, nivel_salario)
summary(ANOVA1)

Observe que para poder utilizar la función aov la variable dependiente debe ser un factor. Se puede comprobar de qué tipo es con class:

class(nivel_salario$educacion)
  • Si no fuera un factor, para conseguir que lo sea, bastaría con emplear la función factor, escribiendo factor(nivel_salario$educacion).

  • Se pueden calcular las medias para cada método, mediante la función tapply,

tapply(nivel_salario$salario, nivel_salario$educacion, mean)

Gráfico de medias

Podemos igualmente, representar las medias con el comando plotMeans, que está en el paquete RcmdrMisc. Para ello, hay que instalar previamente el paquete con install.packages(“RcmdrMisc”) y cargarlo:

library(RcmdrMisc)
plotMeans(nivel_salario$salario, 
          nivel_salario$educacion, 
          error.bars = "conf.int", 
          xlab = "Educacion",
          ylab = "Salario",
          main = "Gráfico de medias por grupos"
          )
  • En el gráfico resultante, se tienen los intervalos de confianza al 95% para los salarios medios de cada nivel de educación.

  • Conclusión: Existe evidencia en los datos de que el salario depende del nivel de educación a un tamaño del test del 5%.

    • Esto es, la educación influye significativamente en el salario de los individuos.

Análisis de diferencia de medias “Análisis post-hoc”

  1. ¿Una mayor educaciónes favorable para el sueldo?. ¿Qué nivel de educación cree que tiene mayores sueldos?
  • Como se han detectado diferencias significativas en los salarios por educación, el paso siguiente es determinar cuáles son niveles de educación son diferentes.

  • Esto se denomina “análisis post-hoc” y consiste en realizar todos los contrastes de igualdad de pares de medias controlando la probabilidad global de error tipo I (comparaciones múltiples).

Para ello se puede utilizar, por ejemplo, el método de Bonferroni, incluido en la función pairwise.t.test:

pairwise.t.test(nivel_salario$salario, nivel_salario$educacion, p.adj = "bonferroni")

Concluciones del análisis post-hoc

  • En la tabla aparecen tres p-valores que permiten establecer qué parejas de medias son significativamente diferentes.
  • Al 5% el salario medio de las personas con licenciatura es significativamente distinta a las personas con doctorado (p-valor = 0.00069 < 0.05)
  • Al 5% el salario medio de las personas con maestría es significativamente diferente a las personas con doctorado (p-valor = 0.00293 < 0.05)
  • Al 5%, no hay diferencias significativas entre el salario medio de las personas con maestría y el salario medio de las personas con licenciatura (p-valor = 1.00000 > 0.05).

Método de Tukey para el análisis post-hoc

También puede emplearse el método de Tukey con la función TukeyHSD.

TukeyHSD(ANOVA1, "educacion")
  • Este método nos ofrece intervalos de confianza (por defecto al 95%) para la diferencia de las medias poblacionales entre cada par de métodos.
  • Las diferencias de medias para las que el intervalo de confianza correspondiente no contiene el valor 0 son estadísticamente significativas.
  • Si el intervalo de confianza contiene el valor 0, no hay diferencia significativa la estatura media.
  • A la misma conclusión se llega observando el p-valor.
  • Esta función incluye, además, un gráfico específico para representar los intervalos de confianza.
plot(TukeyHSD(ANOVA1, "educacion"))

Conclusión

Conclusión: Al 5% solo se aprecian diferencias significativas entre los niveles de educación Licenciatura-Doctorado y Maestria-Doctorado.

ANOVA de dos Factores

El ANOVA de dos factores es el segundo modelo básico de análisis de la varianza. No puede un modelo multivariante, al involucrar únicamente a dos variables, la variable respuesta y dos factores explicativos.

  • A continuación presentaremos los aspectos teóricos del Análisis de Varianza de dos factores

Sean

\[\begin{align*} Y = & \text{ variable respuesta, cuantitativa continua}\\ X_i= & \text{ variables explicativas, categórica con II categoría y varios niveles cada una} \end{align*}\]

Aspectos generales del Anova de dos factores

Este diseño de anova que permite estudiar simultáneamente los efectos de dos fuentes de variación.

  • En el ejemplo 1, en el que se estudiaron los salarios de acuerdo al nivele educativo, se podría plantear que, quizás, otro factor, por ejemplo el sexo y ver si es diferente para los hombres y las mujeres, en cuyo caso, y si el número de hombres y mujeres en cada muestra no fuera el mismo, podría ocurrir que una parte del efecto atribuido al nivel educactivo fuera debido al sexo.

  • En cualquier caso, el investigador puede estar interesado en estudiar si hay, o no, diferencia en los salarios según el sexo.

  • En un anova de dos vías se clasifica a los individuos de acuerdo a dos factores (o vías) para estudiar simultáneamente sus efectos.

  • Si pensamos en el ejemplo 1, si tuvieramos 10 hombres y 10 mujeres, entonces, el primer factor(nive educactivo) tiene \(a\) niveles y el segundo(sexo) tiene \(b\), se tendrán \(ab\) muestras o unidades experimentales, cada una con \(n\) individuos o repeticiones.

Aspectos teóricos: ANOVA de dos Factores

Se dispone de una muestra de n observaciones de la variable Y, que podemos considerar estructuradas de la siguiente forma:

donde, en general, \(y_{jh}\) representa la observación \(h\)-ésima (con \(h = 1,2,...,n_j\)) en la columna o categoría \(j\), y \(n_1, n_2, ..., n_I\) los tamaños de muestra respectivos en las \(I\) categorías de la variable \(X\).

Aspectos teóricos: ANOVA de dos Factores

  • Una observación individual se representa como:

\(Y_{ijk}\), donde \(i=1,2...,a\),\(j=1,2...,b\) y \(j=1,2...,n\),

  • El primer subíndice indica el nivel del primer factor, el segundo el nivel del segundo factor y el tercero la observación dentro de la muestra. Los factores pueden ser ambos de efectos fijos (se habla entonces de modelo I), de efectos aleatorios (modelo II) o uno de efectos fijos y el otro de efectos aleatorios (modelo mixto). El modelo matemático de este análisis es:

\[\begin{align*} Y_{ijk} = & \mu + \alpha_i+\beta_j+(\alpha\beta)_{ij}+e_{ijk} \text{ Modelo I}\\ Y_{ijk} = & \mu + A_i+B_j+(AB)_{ij}+e_{ijk} \text{ Modelo II}\\ Y_{ijk} = & \mu + \alpha_i+\ beta_j+(\alpha B)_{ij}+ e_{ijk} \text{ Modelo III (mixto)} \end{align*}\]

donde \(\mu\) es la media global, ai o \(Ai\) el efecto del nivel \(i\) del factor 1, \(B_j\) el efecto del nivel \(j\) del factor 2 y \(e_{ijk}\) los errores aleatorios alrededor de las medias, que también se asume que están normalmente distribuidas, son independientes y tienen media 0 y varianza \(\sigma^2\).

Aspectos teóricos: ANOVA de dos Factores

  • A las condiciones de muestreo aleatorio, normalidad e independencia, este modelo añade la de aditividad de los efectos de los factores.

  • A los términos \((\alpha\beta)_{ij}\), \((AB)_{ij}\), \((\alpha B)_{ij}\), se les denomina interacción entre ambos factores y representan el hecho de que el efecto de un determinado nivel de un factor sea diferente para cada nivel del otro factor.

Ejemplo ilustrativo

Para entender mejor este concepto de interacción veamos un ejemplo sencillo sobre un anova de dos factores, cada uno con dos niveles: supóngase un estudio para analizar el efecto de un somnífero teniendo en cuenta el sexo de los sujetos. Se eligen al azar dos grupos de hombres y otros dos de mujeres. A un grupo de hombres y otro de mujeres se les suministra un placebo y a los otros grupos el somnífero. Se mide el efecto por el tiempo que los sujetos tardan en dormirse desde el suministro de la píldora.

Ejemplo ilustrativo

Se trata de un anova de dos factores (sexo y fármaco) fijos, cada uno con dos niveles (hombre y mujer para el sexo y somnífero y placebo para el fármaco). Los dos tipos de resultados posibles se esquematizan en la figura

  • En la figura A se observa que las mujeres tardan más en dormirse, tanto en el grupo tratado como en el grupo placebo (hay un efecto del sexo) y que los tratados con placebo tardan más en dormirse que los tratados con somnífero en ambos sexos (hay un efecto del tratamiento). Ambos efectos son fácilmente observables.

  • Sin embargo en la figura B es difícil cuantificar el efecto del somnífero pues es distinto en ambos sexos y, simétricamente, es difícil cuantificar el efecto del sexo pues es distinto en ambos grupos de tratamiento. En este caso, se dice que existe interacción.

Ejemplo ilustrativo

Podría, incluso, darse el caso de que se invirtieran los efectos de un factor para los distintos niveles del otro, es decir, que las mujeres se durmieran antes con el somnífero y los hombres antes con el placebo.

  • La interacción indica, por tanto, que los efectos de ambos factores no son aditivos: cuando se dan juntos, su efecto no es la suma de los efectos que tienen cuando están por separado, por lo que, si en un determinado estudio se encuentra interacción entre dos factores, no tiene sentido estimar los efectos de los factores por separado.

  • A la interacción positiva, es decir, cuando el efecto de los factores actuando juntos es mayor que la suma de efectos actuando por separado, en Biología se le denomina sinergia o potenciación y a la interacción negativa inhibición.

  • En el ejemplo que estamos viendo, se diría que el ser mujer inhibe el efecto del somnífero, o que el ser hombre lo potencia (según el sexo que se tome como referencia).

Descomposición de la varianza

Suma de cuadrados total

La suma de cuadrados total en un anova de 2 factores, es: \[\displaystyle SS_T=\sum_{i=1}^a\sum_{j=1}^{b} \sum_{k=1}^n(y_{ijk}-\bar{y}_{...})^2\]

(donde para representar las medias se ha usado la convención habitual de poner un punto (.) en el lugar del subíndice con respecto al que se ha sumado) que dividida por sus grados de libertad, \(abn - 1\), estima la varianza \(s^2\) en el supuesto de que las \(ab\) muestras provengan de una única población.

Se puede demostrar que

\[SS_T = SS_{\alpha}+SS_{\beta} + SS_{\alpha\beta}+ SS_e\]

Supuestos del modelo

El ANOVA parte de algunos supuestos o hipótesis que han de cumplirse:

  • Normalidad de los datos.
  • Independencia de las observaciones.
  • La distribución de los residuales debe ser normal.
  • Homocedasticidad: homogeneidad de las varianzas.

Tabla ANOVA de dos Factores

A continuación presentamos la tabla anova de dos factores.

Contraste de hipótesis-Modelo I

Los contrastes de hipótesis planteados son los siguientes

  • No existe interacción (MSAB/MSE) \[H_0: (\alpha\beta)_{ij}=0,\] para \(i=1,2,...,a\) y \(j=1,2,...,b\)

  • No existe efecto del primer factor, es decir, diferencias entre niveles del primer factor (MSA/MSE) \[H_0: \mu_1=\mu_2=\dots=\mu_a\]

  • No existe efecto del segundo factor ( MSB/MSE) \[H_0: \mu_1=\mu_2=\dots=\mu_b\]

  • Si se rechaza la primera hipótesis de no interacción, no tiene sentido contrastar las siguientes. En este caso lo que está indicado es realizar un análisis de una vía entre las ab combinaciones de tratamientos para encontrar la mejor combinación de los mismos.

Contraste de hipótesis-Modelo II (Ejemplo)

  • La interacción se contrasta, como en el modelo I, con MSAB/MSE, si se rechaza la hipótesis nula se contrastarían cada uno de los factores con MSA/MSAB y MSB/MSAB.

  • A continuación mostraremos un ejemplo del Modelo II

Ejercicio

A partir de la siguiente tabla de un anova de 2 factores, responda las preguntas.

  • Realizar los contrastes adecuados e interprete.
  • Escriba el modelo.

Ejemplo en R

Vamos a análisar los datos nivel_salario pero con los factores educación y sexo. En este caso se trata de un test ANOVA de dos factores, al tener una variable dependiente (salario) y dos variables cualitativas o factores (educación y sexo). Con más de un factor ya no puede usarse la función oneway.test, por lo que en este ejemplo trabajaremos con aov.

Es un diseño balanceado o equilibrado, puesto que las muestras de los grupos, determinados por las combinaciones de los niveles de los dos factores, tienen el mismo tamaño (18 observaciones).

¿Qué vamos a trabajar?

  1. Lea los datos nivel_salario
  2. ¿Se verifican los supuestos del modelo de análisis de la varianza para el factor (variable) salario?
  3. ¿existe evidencia de que el efecto en salario es distinto según la educación o el sexo?
  4. ¿Qué nivel de educación tiene mejor sueldo?

Solución

  1. Lea los datos
nivel_salario <- read_csv("C:/Users/Roberto Trespalacio/Documents/sem_3.3/datos/nivel_salario3.csv")
head(nivel_salario)

Estamos ante un análisis balanceado ya que todas las muestras tienen el mismo tamaño (hay 6 individuos en cada nivel educativo).

table(nivel_salario$educacion)

Supuestos del modelo

  1. ¿Se verifican los supuestos del modelo de análisis de la varianza?

Para poder aplicar esta técnica estadística se han de verificar las siguientes condiciones:

  • Normalidad en los salarios por tipo de educación y sexo.
  • Homocedasticidad o igualdad de varianzas poblacionales del salario por educacion y sexo.
  • Independencia de las observaciones, que se consigue extrayendo de las poblaciones muestras aleatorias independientes. Esta hipótesis suele ser asumida como consecuencia del diseño del experimento, pues los distintos niveles del factor son asignados al azar (en este caso concreto, los individuos son elegidos aleatoriamente). Por ello, solo comprobaremos los otros dos supuestos.

Exploración de los datos

  • Antes de verificar si se cumplen las condiciones necesarias para llevar a cabo un ANOVA, exploraremos los datos mediante un gráfico de caja:
boxplot(nivel_salario$salario ~ nivel_salario$educacion)  

# tambien puede escribirse
boxplot(salario ~ educacion, nivel_salario)  

boxplot(nivel_salario$salario ~ nivel_salario$sexo)  

# tambien puede escribirse
boxplot(salario ~ sexo, nivel_salario)  

Interpretación

Normalidad

  • Mediante el contraste de normalidad de Shapiro-Wilk, se puede comprobar si una variable sigue una distribución normal.
  • Se contrasta la hipótesis nula de que la distribución de la que proceden los datos es normal para cada método:
tapply(nivel_salario$salario, nivel_salario$educacion, shapiro.test)
  • En cada uno de los grupos, es decir, para cada uno de los métodos, no rechazamos la hipótesis nula de normalidad al 5% (al ser el p-valor superior a 0.05 en todos los casos).
  • En la muestra no hay evidencia para rechazar que los datos en cada grupo provienen de una distribución normal.

Normalidad para factor educación

tapply(nivel_salario$salario, nivel_salario$educacion, shapiro.test)

Normalidad para el factor sexo

tapply(nivel_salario$salario, nivel_salario$sexo, shapiro.test)

Homocedasticidad

Homocedasticidad para el factor educacion

tapply(nivel_salario$salario, nivel_salario$educacion, sd)

Homocedasticidad para el factor sexo

tapply(nivel_salario$salario, nivel_salario$sexo, sd)
  • Se observa cierta diferencia en la variabilidad de tiempos por grupos, pero para contrastar si esa diferencia es estadísticamente significativa emplearemos un test de homogeneidad de varianzas. - Uno de los test disponibles es el de Bartlett, que se puede aplicar solo cuando las observaciones provienen de poblaciones normales.
  • La hipótesis nula es que las varianzas son idénticas en todos los grupos y, por tanto, se cumple el requisito de homocedasticidad.
  • La función bartlett.test se puede usar de la siguiente forma, incluyendo en primer lugar la variable dependiente y después el factor,
bartlett.test(nivel_salario$salario, nivel_salario$educacion)

# tambien puede escribirse así
bartlett.test(salario ~ educacion, nivel_salario)
bartlett.test(nivel_salario$salario, nivel_salario$sexo)

# tambien puede escribirse así
bartlett.test(salario ~ sexo, nivel_salario)
  • Como el p-valor es superior a 0.05, no rechazaríamos la hipótesis nula de igualdad de varianzas al 5%, y, por tanto, no hay evidencia en contra del requisito de homocedasticidad.

  • Si el factor únicamente tuviera dos niveles podríamos utilizar también la función var.test, que está indicada para el contraste de varianzas en dos muestras provenientes de poblaciones normales, como ya se vio en la anterior clase.

Efectos de los factores

  1. ¿existe evidencia de que hay diferencia de salario según el nivel de educación?

1. No existe interacción (MSAB/MSE)

\[H_0: (\alpha\beta)_{ij}=0,\] para \(i=1,2,...,a\) y \(j=1,2,...,b\)

Es decir, si hay interacción entre los niveles de educación y sexo

2. No existe efecto del primer factor, es decir, diferencias entre niveles del primer factor (MSA/MSE)

\[H_0: \mu_1=\mu_2=\dots=\mu_a\]

Es decir, si hay diferencias entre salario para los niveles de educación

3. No existe efecto del segundo factor (MSB/MSE)

\[H_0: \mu_1=\mu_2=\dots=\mu_b\]

Es decir, si hay diferencias entre salario para homres y mujeres.

Efectos de los factores

  • Si se rechaza la primera hipótesis de no interacción, no tiene sentido contrastar las siguientes. En este caso lo que está indicado es realizar un análisis de una vía entre las ab combinaciones de tratamientos para encontrar la mejor combinación de los mismos.

  • En este caso se trata de un test ANOVA de dos factores, al tener una variable dependiente (salario) y dos variables cualitativas o factores (Educación y Sexo). Usaremos la función de R, aov.

  • Es un diseño balanceado o equilibrado, puesto que las muestras de los seis grupos, determinados por las combinaciones de los niveles de los dos factores, tienen el mismo tamaño (18 observaciones). Como ya lo comprobamos.

Interacción entre el educación y sexo

En primer lugar, estudiamos si existe interacción entre el sexo y la dieta, esto es, si el efecto del nivel de educación en el salario depende del sexo. Para incluir la interacción en el modelo se usa el signo * entre los factores dentro de la función aov:

anova2 <- aov(salario ~ educacion * sexo , nivel_salario)
summary(anova2)
  • Observamos que el factor educación, nos dice que hay diferencias significativas en al menos dos medias del factor (p-valor = 0.00137)
  • El factor sexo, nos dice que no hay diferencias significativas en las medias del factor (p-valor = 0.56059)
  • Vemos también, que no existe interacción significativa entre los factores (educación y sexo) (p-valor = 0.60693).
interaction.plot(nivel_salario$educacion,
                 nivel_salario$sexo, 
                 nivel_salario$salario,  
                 col = c("red", "blue"), 
                 ylab = "Salario", 
                 xlab = "Educación",
                 trace.label = "Sexo")
  • Si representamos el gráfico de los salarios por nivel de educación y sexo, se tiene que las líneas correspondientes a los distintos niveles del segundo factor en el comando interaction.plot son más o menos “paralelas”(con una insersección que al perrecer no es significativa), indicando la ausencia de interacción significativa.

¿Y si no hay intercacción?

Una vez comprobada la ausencia de interacción, procedemos a realizar el ANOVA de dos factores considerando solo los efectos principales. Para ello, usamos el signo + para separar la lista de factores, en lugar del signo * empleado anteriormente. Así:

anova3 <- aov(salario ~ educacion + sexo , nivel_salario)   
summary(anova3) 

Observamos diferencias significativas entres los niveles de educación, mientras que en el factor sexo no hay vemos diferencias.

boxplot(salario ~ educacion, nivel_salario)
boxplot(salario ~ sexo, nivel_salario)

Análisis post-hoc

TukeyHSD(anova3, "educacion")
TukeyHSD(anova3, "sexo")

Vemos que no muestra una diferencia significativa respecto al sexo.

Ejercicio: ANOVA de dos factores con efectos fijos

Una forma de medir el estado físico de una persona es medir su porcentaje de grasa corporal. El porcentaje de grasa corporal promedio varía con la edad, pero según ciertas pautas, el intervalo normal para hombres es del 15-20 % de grasa corporal, y para mujeres, del 20-25 %. Los datos de muestra vienen de un grupo de hombres y mujeres que hicieron ejercicio en un gimnasio tres veces por semana durante un año. Luego, su entrenador medıa la grasa corporal.

Dieta
Paciente A B C
Hombre 13 16 20
18 14 19
18 25 16
Mujer 22 15 15
13 22 14
19 21 28

En este caso se trata de un test ANOVA de dos factores, al tener una variable dependiente (porcentage de grasa) y dos variables cualitativas o factores (dieta y sexo). Con más de un factor ya no puede usarse la función oneway.test, por lo que en este ejemplo trabajaremos con aov.

Es un diseño balanceado o equilibrado, puesto que las muestras de los seis grupos, determinados por las combinaciones de los niveles de los dos factores, tienen el mismo tamaño (21 observaciones). Para comprobarlo

  1. Lea los datos grasa
  2. ¿Se verifican los supuestos del modelo de análisis de la varianza para el factor (variable) grasa?
  3. ¿existe evidencia de que el efecto en el porcentaje de grasa es distinto según el sexo?
  4. ¿Qué dieta cree que tiene mejores resultados?