ANOVA

Comparación de modelos

Gonzalo J. Muñoz

2019-04-03

Para entender la lógica de ANOVA usaré la perspectiva de comparación de modelos que es una manera general para determinar la adecuación de un modelo lineal (Maxwell, Delaney, and Kelley 2018Maxwell, S. E., H. D. Delaney, and K. Kelley. 2018. Designing Experiments and Analyzing Data: A Model Comparison Perspective. 3rd ed. New York, NY: Routledge.) Volvamos al ejemplo inicial que usamos en la Introducción.

Nivel de satisfacción con el restaurant para los distintos tipos de música.

sujeto condicion satisfaccion
s1 Música en vivo 3
s2 Música en vivo 4
s3 Música en vivo 1
s4 Música en vivo 2
s5 Música en vivo 0
s6 Música en vivo 2
s7 Música envasada 6
s8 Música envasada 2
s9 Música envasada 4
s10 Música envasada 3
s11 Música envasada 4
s12 Música envasada 5
s13 Sin música 7
s14 Sin música 7
s15 Sin música 5
s16 Sin música 5
s17 Sin música 8
s18 Sin música 4

Imaginemos que queremos predecir el valor s1 pero no sabemos aun a cuál de las condiciones experimentales ha sido expuesto. Nuestra mejor predicción para s1 es el promedio del grupo o gran media. De hecho, la gran media es nuestra mejor predicción para todos los sujetos de la muestra en el sentido que cualquier otro valor distinto de la media nos entregaría un error de predicción mayor.

## Loading required package: ggplot2

El error de predicción para este modelo lo podemos obteniendo la suma de las diferencia al cuadrado de cada puntaje respecto de la gran media.1 Si no elevamos al cuadrado estas diferencias la suma de estos valores sería cero.

\[\sum_{i=1}^{n}e^2_i = \sum_{i=1}^{n}(Y_i-\overline{Y})^2\]

Estos es, réstele a cada valor observado (1…n) la media muestral, eleve el resultado al cuadrado, y sume los valores para cada caso. Para nuestro caso, la suma de cuadrados para este modelo es 80.

EL modelo anterior tiene la restricción de basarse únicamente en la media del grupo. Sin embargo, como sabemos, existe un factor que está asociado a la satisfacción de quienes asisten al restaurant y que nos va a ayudar a explicar mejor (en teoría) los reesultados que observamos. Es decir, podemos generar un modelo alternativo que incluya el efecto del tratamiento sobre la variable dependiente. Nuestra presunción es que incluir este factor en el modelo nos ayudará a reducir el error de predicción del modelo inicial.

De manera análoga a lo que hicimos en el paso anterior, vamos ahora a determiar el error de predicción de cada sujeto pero ahora en relación a la media de la condición experimental a la que pertenece.

\[\sum_j\sum_i e^2_{ij} = \sum_{j=1}^j\sum_{i=1}^{n_j}(Y_{ij}-\overline{Y_j})^2\]

donde j es un índice que representa el número de grupo y i el sujeto. En este caso particular tenemos 3 grupos (con 6 personas por grupo) por lo que el error de predicción corresponde a la suma del error de predicción para cada grupo. Es decir, el error de predicción se obtiene para cada sujeto en relación a la media de su grupo. En este caso ese valor es 32.

El hecho de que el error de predicción del modelo inicial (i.e., 80) sea mayor que el error de predicción del modelo alternativo (i.e., 32) sugiere que el modelo alternativo es preferible en el sentido que usar este modelo reduce error de predicción del modelo inicial. Le llamatemos al modelo inicial modelo restringido del inglés restricted model y al modelo alternativo modelo saturado (del inglés full model. Una manera de estandarizar la manera en que comparamos modelos es usando una proporción entre el incremento en el error asociado al modelo restringido en relación al modelo saturado.

Reordemos que el modelo restringido es más simple y por lo tanto es esperable que haya más error asociado a este modelo. El modelo saturado es más complejo pero tiende a minimizar el error de predicción. Conceptualmente, la comparación de modelos se vería de esta forma:

\[incrementoproporcional del error = \frac{incremento en el error}{mínimo error}\]

Sin embargo, para que la comparación sea “justa” hay que incorporar el hecho de que el modelo restringido se basa en sólo un parámetro (la gran media) mientras que el modelo alternativo utiliza tres parámetros (la media de cada grupo). Si llevamos esto al extremo, podríamos tomar como valor predicho el valor observado de cada individuo y tendríamos un modelo que se ajusta perfectamente a los datos pero que sería poco informativo. Como en ciencia nos interesa tener modelos simples, vamos a “castigar” al modelo alternativo por su relativa complejidad.

En general, cada parámetro adicional que utiliza un modelo supone la pérdida de un grado de libertad (gl). Los gl se definene como el número de observaciones independientes de un estudio menos el número de parámetros estimados. Para nuestro ejemplo, el modelo restringido tiene \(N-1=18-1=17\) gl y el modelo saturado \(N-3=18-3=15\).

Integrando la noción de grados de libertad para comparar modelos tenemos

\[F=\frac{(E_R - E_S)/(gl_R-gl_S)}{E_S/gl_S}\]

Nótese que esta versión revisada de la fórmula anterior da como resultado un estadístico, que se denomina F, y que se utiliza para determinar la adecuación relativa de ambos modelos.

Reemplazando con los datos del ejemplo anterior tenemos que el valor (o razón) F en este caso es:

\[F=\frac{(80 - 32)/(17-15)}{32/15}=\frac{48/2}{32/15}=\frac{24}{2.1333}=11.25\] Para determinar si el valor de F es estadísticamente debemos determinar el valor crítico de F para \(\alpha\) 0.05 con gl 2 y 15. En R es posible obtener el valor crítico mediante qf(P,glb,glw)

qf(.95,2,15)
## [1] 3.68232

Como \(F_{obs} = 11.25 > F_{crit} = 3.68\) se rechaza la hipótesis nula (más sobre esto más adelante)

El resultado de comparar ambos modelos es equivalente al modelo que evaluamos en la introducción mediante la función aov2 ¡BAM!

aov1<-aov(satisfaccion~condicion,data=df)
summary(aov1)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## condicion    2     48  24.000   11.25 0.00104 **
## Residuals   15     32   2.133                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Matemáticamente, podemos reducir la fórmula anterior para obtener la razón F se a la siguiente expresión:

\[F=\frac{[\sum_{j=1}^an_j(\overline{Y}_j-\overline{Y})^2]/(a-1)}{\sum_{j}\sum_{i=1}^{n_j}(Y_{ij}-\overline{Y}_j)^2/(N-a)}=\frac{MS_{B}}{MS_{W}}\] donde a corresponde al número de grupos y N al número de observaciones o sujetos.

El numerador en la fórmula anterior se conoce tradicionalmente como media cuadrática inter-grupos (mean square between o \(MS_{B}\)) y el denominador media cuadrática intra-grupos (mean square within o \(MS_{W}\)). En esta concepción, la razón F es básicamente una comparación de la varianza entre medias de grupos (\(MS_{B}\)) y la varianza intra-grupo o varianza error (\(MS_{within}\)). Cuando \(MS_{B]\) es mayor que (\(MS_{W}\)) diríamos que las diferencias entre las medias grupales son mayores que las diferencias que uno esperaría “por azar”. Más adelante veremos que es posible asignarle una probabilidad a este valor F y tomar en relación a esa probabilidad una decisión de rechazo o no rechazo de la hipótesis nula.

References