Cuando se comparan tres o más grupos los resultados del test F indican que existe (o no) alguna diferencia entre las medias de los grupos comparados. Sin embargo, usuamente el investigador está interesado en comparaciones específicas entre algunos grupos. Por ejemplo, en el contexto de un estudio con 4 condiciones experimentales podemos plantear las siguientes hipótesis:
Caso 1: \[H_o=\mu_1-\mu_2=0\] Caso 2: \[H_o=\frac{1}{3}(\mu_1+\mu_2+\mu_3)=\mu_4\]
En el primer caso, nuestro interés es comparar el grupo 1 y el grupo 2 entre ellos. En el segundo caso, nuestro interés es comparar el grupo 4 con el promedio de los grupos 1 a 3. En el contexto de ANOVA este tipo de hipótesis específicas se puede analizar mediante lo que se conoce como contrastes o comparaciones. Los contrastes son simples cuando se comparan medias individuales (Caso 1) o complejos la comparación involucra agrupar medias de alguna manera (Caso 2). Un contraste es una combinación lineal de medias que suma cero. Específicamente, esto se consigue asignándole un coeficiente (c) a cada una de las medias que refleje las comparaciones que queremos hacer
Caso 1: \[H_o:1\mu_1-1\mu_2+0\mu_3+0\mu_4=\mu_1-\mu_2=0\] Caso 2: \[H_o:\frac{1}{3}\mu_1+\frac{1}{3}\mu_2+\frac{1}{3}\mu_3-1\mu_4\]
Aunque existe un número infinito de coeficientes que sumados dan cero, una restricción adicional que recomiendan 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.) (p. 182) es que los coeficientes positivos sumen 1 y los negativos -1, lo que se cumple en ambos ejemplos.
Los contrastes se representan usando la letra griega psi (\(\psi\)) y se pueden representar conjuntamente como una sumatoria de la forma, \[\psi=\sum_{j=1}^{a}c_j\mu_j\] Que se puede estimar mediante \[\hat{\psi}=\sum_{j=1}^{a}c_j\overline{Y}_j\] Por último, el valor de F para un contraste se obtiene mediante
\[F=\frac{\hat{\psi}^2/\sum_{j=1}^{a}(c_j^2/n_j)}{MS_W}\] Que se usa para evaluar \[H_0:\psi=\sum_{j=1}^{a}c_j\mu_j=0\] Esto es, cuando Ho es cierta la suma de las medias poblacionales (\(\mu_j\)) ponderadas por sus respectivos coeficientes es igual a cero.
Veamos un ejemplo de cómo usar contrastes con la data del restaurant.
sujeto<-paste0("s",1:18) # paste adds a space between the characters
condicion<-c(rep("Música en vivo",6),rep("Música envasada",6),rep("Sin música",6))
satisfaccion <- c(3,4,1,2,0,2,6,2,4,3,4,5,7,7,5,5,8,4)
df <-data.frame(sujeto,condicion,satisfaccion)
df## sujeto condicion satisfaccion
## 1 s1 Música en vivo 3
## 2 s2 Música en vivo 4
## 3 s3 Música en vivo 1
## 4 s4 Música en vivo 2
## 5 s5 Música en vivo 0
## 6 s6 Música en vivo 2
## 7 s7 Música envasada 6
## 8 s8 Música envasada 2
## 9 s9 Música envasada 4
## 10 s10 Música envasada 3
## 11 s11 Música envasada 4
## 12 s12 Música envasada 5
## 13 s13 Sin música 7
## 14 s14 Sin música 7
## 15 s15 Sin música 5
## 16 s16 Sin música 5
## 17 s17 Sin música 8
## 18 s18 Sin música 4Hay distintas opciones para ingresar los contrastes. Primero, usando la fórmula anterior podemos obtener un valor para F mediante el siguiente código
# Primero, separamos la data en tres grupos para
# obtener los estadísticos descriptivos que usaremos
# más adelante
Group.1 <- df[df$condicion=="Música en vivo", "satisfaccion"]
Group.2 <- df[df$condicion=="Música envasada", "satisfaccion"]
Group.3 <- df[df$condicion=="Sin música", "satisfaccion"]
# Descriptivos
Y.bar1 <- mean(Group.1)
s2.1 <- var(Group.1)
n.1 <- length(Group.1)
Y.bar2 <- mean(Group.2)
s2.2 <- var(Group.2)
n.2 <- length(Group.2)
Y.bar3 <- mean(Group.3)
s2.3 <- var(Group.3)
n.3 <- length(Group.3)A continuación ingresamos los contrastes para evaluar la diferencia entre el grupo con Música en vivo y Música envasada. Al grupo sin música le asignamos un coeficiente de cero. El valor de \(\psi\) en este caso es -2, que resulta de sutraer 2 (la media del grupo con Música en vivo) de 4 (la media del grupo con Música envasada).
#vector coeficientes
c.weights <- c(1,-1, 0)
#vector de medias de los grupos
Ybar <- c(Y.bar1, Y.bar2, Y.bar3)
#vector de tamaños de los grupos
n <- c(n.1, n.2, n.3)
# Obtenemos psi sumando el producto de
# cada media y su respectivo coeficiente
Psi <- sum(c.weights*Ybar)
# Creamos un objeto para extraer algunos valores
# de esta lista
ANOVA.Object <-aov(satisfaccion ~ as.factor(condicion), data=df)
ANOVA_out <- anova(ANOVA.Object)
# Así extraemos MSw
MS.W <- ANOVA_out$"Mean Sq"[2]
# Y combinamos esta información para obtener F
# de acuerdo a la fórmula de más arriba
F <- (Psi^2/sum(c.weights^2/n))/MS.W
F## [1] 5.625Los grados de libertad del modelo restringido (bajo Ho) corresponde a N - (a-1). La razón es que al asumir \(\mu_1-\mu_2=0\) sólo un parámetro representa la diferencia entre ambos grupos. El resto de los grupos tiene cada uno su propio parámetro. Los grados de libertad del modelo saturado es igual a N - a. Por lo tanto, los grados de libertad para el denominador de F es, \[[df_R]-[df_F]=[N-(a-1)]-[(N-a)]\] \[df_R-df_F=N-a+1-N+a\] \[df_R-df_F=1\] El valor crítico para F en este caso es 4.54
qf(.95,1,15)## [1] 4.543077Luego, dado que \(F_obs = 5.625 > F_crit = 4.543\) se rechaza la hipótesis nula. En otras palabra, el valor estimado de \(\hat{\psi}=-2\) es estadísticamente significativamente menor que cero, lo que quiere decir que existe menos satisfacción cuando hay un grupo con música en vivo en comparación con un grupo con música envasada
Vamos a evaluar dos contrastes. Primero, replicando los resultados obtenidos en el paso anterior tenemos que
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#Requiere el paquete emmeans
#install.packages("emmeans")
require(emmeans)## Loading required package: emmeans## Warning: package 'emmeans' was built under R version 3.5.3# Primero creamos un objeto del tipo 'emmGrid'
# que usaremos luego para evaluar diferencias
# entre grupos.
# specs especifica el factor
aov1.emm = emmeans(aov1, specs=~condicion)
# Para comparar el grupo musica en vivo
# vs el grupo música envasada
Contrasts.1 = list(Musica_en_vivo_vivo_vs_Musica_Envasada = c(1,-1,0))
contrast(aov1.emm, Contrasts.1)## contrast estimate SE df t.ratio p.value
## Musica_en_vivo_vivo_vs_Musica_Envasada -2 0.843 15 -2.372 0.0315Aunque en este caso se usa el estadístico t, el valor p (0.0315) asociado al valor t observado (-2.372) es el mismo que el valor p para F(1,15)=5.625. Esto es así porque \(t^2=F\).
1-pf(5.625,1,15)## [1] 0.03151799Si queremos evaluar más de un contraste, hay que hacer una corrección o ajuste a los valores p.
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#Requiere el paquete emmeans
#install.packages("emmeans")
require(emmeans)
aov1.emm = emmeans(aov1, specs =~condicion)
# Agreguemos un contraste para la diferencia entre la
# condición sin música y las otras dos condiciones
# combinadas
Contrasts.2 = list(Musica_en_vivo_vivo_vs_Musica_Envasada = c(-1,1,0),
Todos_vs_Sin_Musica = c(1/2,1/2,-1))
contrast(aov1.emm, Contrasts.2)## contrast estimate SE df t.ratio p.value
## Musica_en_vivo_vivo_vs_Musica_Envasada 2 0.843 15 2.372 0.0315
## Todos_vs_Sin_Musica -3 0.730 15 -4.108 0.0009Es importante notar que los valores p asociados a los valores t observados requieren un ajuste ya que se realizó más de una prueba t. En este caso, debemos usar un valor p ajustado igual a \(\alpha/C\), donde C representa el número de contrastes. en este caso, C = 2, por lo que debemos evaluar t contra un valor \(\alpha=.05/2=0.025\). De acuerdo a este ajuste sólo la diferencia entre el grupo sin música y las otras dos condiciones es estadísticamente significativa. En este caso, esto tiene sentido porque la conclusión del estudio es simplemente que no es bueno usar música, ni envasada ni no envasada.
El auste de Bonferroni también se puede solicitar directamente a través de la función contrast especificando la opción “bonferroni” en ajustes (adjust). En este caso, el ajuste se ve reflejado directamente en el valor p ajustado.
contrast(aov1.emm,Contrasts.2,adjust = "bonferroni")## contrast estimate SE df t.ratio p.value
## Musica_en_vivo_vivo_vs_Musica_Envasada 2 0.843 15 2.372 0.0630
## Todos_vs_Sin_Musica -3 0.730 15 -4.108 0.0019
##
## P value adjustment: bonferroni method for 2 testsCuando más de un contraste se aplica a un conjunto de datos (como en el caso del Contraste.2) hay que considerar en qué medida dichos contrastes añaden información redundante o sobrelapada. Cuando los contrastes son linealmente independientes, los valores obtenidos de \(\psi\) no pueden ser reproducidos mediante una combinación lineal de los otros contrastes. Por ejemplo:
\[\psi_1=\mu_1-\mu_2\] \[\psi_2=\mu_1-\mu_3\] \[\psi_3=\frac{1}{2}(\mu_1+\mu_2)-\mu_3\]
En este caso toda la información de \(\psi_3\) puede ser obtenida combinando \(\psi_1\) y \(\psi_2\) ya que \(\psi_3=\)_2-_1$
Por otro lado (y tal vez más importante) se dice que dos contrastes son ortogonales cuando los tamaños de muestra de los grupos es la misma y la suma de sus productos es igual a cero. \[\sum c_{1j}c_{2j}=0\] En el caso de los contrastes que usamos en nuestro ejemplo, podemos decir que ambos contrastes son ortogonales ya que,
\[(-1)(1/2)+(1)(1/2)+(0)(-1)=0\] Cuando esta condición no se cumple existe un solapamiento entre la varianza explicada por cada contraste. Adicionalmente, cuando los tamaños de los grupos son desiguales los contrastes necesariamente serán no-ortogonales. Veremos como resolver este problema más adelante, pero en general, hay que usar lo que se conoce como suma de cuadrados Tipo I o Tipo III. En el caso de la SS Tipo I, cada contraste se ingresa secuencialmente de modo que cada nuevo contraste captura la proporción de varianza en la VD no explicada por el(los) contraste(s) anterior(es). La SS Tipo III sólo considera la varianza en la VD no explicada por ningún otro factor (o contraste en este caso). Por lo tanto, la SS Tipo III es más conservadora que la SS Tipo I (salvo que n por celda sean iguales o muy parecidos).