Ejemplos_ANOVA
Constanza Contreras
4/27/2020
Impacto ubicación estante en las ventas
A consumer-goods manufacturer conducted an experiment to assess the impact of retail shelf location on sales. Two factors where controlled in the experiment: (1) the height in the shelf; and (2) the facings in the shelf. Height was varied on three levels: knee high, waist high and eye level. Facings has two levels: half-width and full width. For each configuration, the retailer collected three sales observations. (The data of the experiment is available in shelf_test.dta, a Stata dataset. HEIGHT and FACINGS are numerical variables with labels. The other variables with extension Str and Num are string and numerical variables).
library(foreign)
library(ggplot2)
shelf <- read.dta("shelf_test.dta")
head(shelf, n=5)## sales HEIGHT FACINGS heightStr facingsStr heightNum facingsNum
## 1 70 knee half knee half 2 2
## 2 75 knee half knee half 2 2
## 3 79 knee half knee half 2 2
## 4 91 knee full knee full 2 1
## 5 90 knee full knee full 2 1Nos piden analizar el impacto de la ubicación de la estantería en las ventas. Es decir, si existe una relación entre la variable \(HEIGHT\) y la variable \(sales\). Para esto, visualizemos los datos con un gráfico de boxplot.
ggplot(data = shelf, aes(x = HEIGHT, y = sales, color = HEIGHT)) +
geom_boxplot() +
theme_bw() #theme_bw : white background and gray grid lines
En el experimento para evaluar el impacto, se posicionaron los productos a tres niveles: a la altura de los ojos (\(eye\)), a la altura de la cintura (\(waist\)) y a la altura de las rodillas (\(knee\)). En el gráfico apreciamos como distribuyen las ventas en cada una de estas categorías. Visualmente existen diferencias en las ventas en posicionar los productos a estos 3 distintos niveles, siendo al nivel de la cintura donde se aprecia el valor mas alto de ventas. ¿Qué deberíamos esperar al comparar las medias de estos tres grupos?
Vamos a realizar un ANOVA de una vía: compararemos las medias de las ventas en estos 3 grupos (categorías). Por lo cual en este caso \(I=3\). Y para cada categoría tenemos 6 observaciones. Entonces \(J=6\). Recordemos que el modelo de ANOVA es el siguiente: \[ Y_{ij} = \mu + \alpha_{i} + \epsilon_{ij}\] Donde \(Y_{ij}\) es la observación de la persona \(j\) en el grupo \(i\), \(\mu\) es la media poblacional de donde provienen todas las observaciones, \(\alpha_{i}\) es el efecto que existe por pertenecer al grupo \(i\) y en \(\epsilon_{ij}\) dejamos todo aquello que no podemos observar, toda la variabilidad entre las observaciones que no es explicada por los factores de nuestro modelo (es un error aleatorio). A modo de ejemplo, en nuestros datos: $ Y_{knee,1} = 70$, \(Y_{knee,2}=75\).
anova_shelf <- aov(shelf$sales~shelf$HEIGHT)
summary(anova_shelf)## Df Sum Sq Mean Sq F value Pr(>F)
## shelf$HEIGHT 2 316.3 158.17 3.695 0.0496 *
## Residuals 15 642.2 42.81
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1¿Cómo interpretamos estos resultados? En la fila \(shelf\$HEIGHT\) se despliegan los resultados ENTRE GRUPOS (lo que corresponde al SSB). En efecto, el valor de SSB es de 316.3. Recordemos que el estimador insesgado para la varianza intergrupos es \(\frac{SSB}{(I-1)}\). En este caso, son 3 grupos, por eso en la columna \(Df\) (grados de libertad) aparece \(I-1=2\). Entonces, el estimador insesgado de la varianza intergrupo nos queda en \(316.3/2=158.17\). Este es el valor que se reporta en la columna “Mean Sq”.
Con respecto a la segunda fila, la que dice “Residuals”, esta corresponde a los resultados INTRA GRUPO. La suma de cuadrados intagrupos (SSW) es de 642.2. El estimador insesgado de la varianza intragrupo es \(\frac{SSW}{I(J-1)}\). En este caso \(I(J-1)=3*(6-1)=15\), que es justamente el valor calculado en los grados de libertad (columna \(Df\)). Luego, el valor del estimador nos queda en \(642.2/15=42.81\) (el valor reportado en la columna “Mean Sq”). Con estos valores, se construye el estadístico del test, el cual es: \[ F = \frac{SSB/(I-1)}{SSW/I(J-1)} = 3.695\]
Este valor es reportado en la columna “F value”. En la columna “Pr(>F)” se entrega el p-valor del test de hipótesis, que en este caso es de 0.0496. Por lo cual, a un nivel de significancia de 5% se rechaza la hipótesis nula del ANOVA.
¿Qué concluimos al rechazar ANOVA? Que AL MENOS UNA DE LAS MEDIAS es distinta de las otras significativamente. En base al gráfico Boxplot, ¿Cuál cree usted que es el grupo distinto?.
Verifiquemos con los test de diferencia de medias entre los 3 grupos:
- Comparamos \(waist\) con \(knee\):
t.test(x = shelf[shelf$HEIGHT == "waist",]$sales,
y = shelf[shelf$HEIGHT == "knee",]$sales,
alternative = "two.sided", mu = 0, var.equal = TRUE, conf.level = 0.95)##
## Two Sample t-test
##
## data: shelf[shelf$HEIGHT == "waist", ]$sales and shelf[shelf$HEIGHT == "knee", ]$sales
## t = 2.4492, df = 10, p-value = 0.0343
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 0.8724937 18.4608396
## sample estimates:
## mean of x mean of y
## 91.66667 82.00000- Comparamos \(waist\) con \(eye\):
t.test(x = shelf[shelf$HEIGHT == "waist",]$sales,
y = shelf[shelf$HEIGHT == "eye",]$sales,
alternative = "two.sided", mu = 0, var.equal = TRUE, conf.level = 0.95)##
## Two Sample t-test
##
## data: shelf[shelf$HEIGHT == "waist", ]$sales and shelf[shelf$HEIGHT == "eye", ]$sales
## t = 2.6103, df = 10, p-value = 0.02603
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 1.146854 14.519813
## sample estimates:
## mean of x mean of y
## 91.66667 83.833333.Comparamos \(eye\) con \(knee\):
t.test(x = shelf[shelf$HEIGHT == "eye",]$sales,
y = shelf[shelf$HEIGHT == "knee",]$sales,
alternative = "two.sided", mu = 0, var.equal = TRUE, conf.level = 0.95)##
## Two Sample t-test
##
## data: shelf[shelf$HEIGHT == "eye", ]$sales and shelf[shelf$HEIGHT == "knee", ]$sales
## t = 0.42941, df = 10, p-value = 0.6767
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -7.679483 11.346150
## sample estimates:
## mean of x mean of y
## 83.83333 82.00000Observamos que los dos primeros test rechazan la hipótesis nula a una significancia del 5% pero no el último test. Lo cual coincide con nuestras apreciaciones del gráfico de BoxPlot.
ANOVA de 2 vías
Ahora, agreguemos la información que nos da la variable \(facing\). Para analizar si hay un efecto con esta variable además del efecto con \(HEIGHT\) podemos realizar un ANOVA de 2 vías. Antes de modelar, veamos visualmente como se comporta esta variable en cada categoría con respecto a las ventas:
ggplot(data = shelf, aes(x = HEIGHT, y = sales, colour = FACINGS, group = FACINGS)) +
stat_summary(fun.y = mean, geom = "point") + stat_summary(fun.y = mean,
geom = "line") + labs(y = "mean (result)") + theme_bw()
Este gráfico se conoce como “Interaction plot” o gráfico de interacción. Es para visualizar si existe un efecto en la interacción de 2 variables. Recordemos que el ANOVA de 2 vías agrega efectos por cada subgrupo y además podemos incoporar un efecto de la interacción. ¿Cómo interpretamos este gráfico? Uno esperaría que, si no existe un efecto de la variable \(facings\) en las categorías de \(HEIGHT\) entonces las curvas debiesen ser paralelas, pero en este caso, vemos que claramente no lo son. Este gráfico existe como función en R. Si no queremos usar ggplot podemos usar la siguiente función.
interaction.plot(shelf$HEIGHT, shelf$FACINGS, shelf$sales, fun = mean,
type = c("l"), legend = TRUE,
trace.label = "FACINGS",
fixed = FALSE,
xlab = "HEIGHT",
ylab = "SALES")
Realicemos ahora un ANOVA de 2 vías. Recordemos el modelo de ANOVA de 2 vías: \[ Y_{ijk} = \mu + \alpha_{i} + \beta_{j} + \delta_{ij} + \epsilon_{ijk}\] Donde \(Y_{ijk}\) es la observación k del grupo \(j\) y grupo \(i\). A modo de ejemplo: \(Y_{knee, half,1} = 70\), \(Y_{knee,full,1} = 90\). \(\alpha_{i}\) es el efecto por pertenecer al grupo \(i\), \(\beta_{j}\) es el efecto por pertenecer al grupo \(j\) y \(\delta_{ij}\) es el efecto de la interacción de pertenecer a ambos grupos (ser grupo \(i\) y ser grupo \(j\)). Extendiendo el ejemplo que hemos usado en clases, I representa a los distintos remedios que probamos en el experimento y J al género (Hombre, mujer). Luego, \(\delta_{remedio1,mujer}\) es el efecto que existe por ser mujer Y tomar el medicamento 1.
En el caso de ANOVA de 2 vías, tengo varias hipotesis nula y construyo el estadístico según el efecto que quiera evaluar. Dicho de otra forma, la suma total de cuadrados (SST) se descompone como:
\[ SST = SSB_{I} + SSB_{J} + SSB_{IJ} + SSW\] Donde \(SSB_{I}\) es la variabilidad entre grupos \(I\) (como lo vimos en el ANOVA de 1 vía), luego aparecen los otros términos por el efecto que agregamos de este nuevo grupo \(J\): \(SSB_{J}\) es la variabilidad entre grupos \(J\) y \(SSB_{IJ}\) es la variabilidad entre los grupos que se forman por pertenecer a ambos grupos a la vez (a cierto grupo \(i\) y cierto grupo \(j\)). Finalmente, \(SSW\) es la variabilidad intra-grupos que viene dada por el error aleatorio.
Para mayor desglose de como se obtienen estos términos revisen el capítulo de ANOVA del libro de John Rice. Esta en el Syllabus cuál capítulo es.
Para realizar ANOVA de 2 vías, solo necesitamos agregar esta nueva variable a la fórmula. Primero realizaremos el modelo sin considerar la interacción de las dos variables (sin el término \(\delta_{ij}\):
anova2_shelf_sin <- aov(shelf$sales~shelf$HEIGHT + shelf$FACINGS)
summary(anova2_shelf_sin)## Df Sum Sq Mean Sq F value Pr(>F)
## shelf$HEIGHT 2 316.3 158.2 13.70 0.000506 ***
## shelf$FACINGS 1 480.5 480.5 41.61 1.52e-05 ***
## Residuals 14 161.7 11.5
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Aquí la interpretación de los resultados es análoga a la de 1 vía. En este caso la fila correspondiente a \(shelf\$HEIGHT\) nos indican los cálculos para el efecto de \(HEIGHT\), lo que sería el \(SSB_{I}\) en nuestra descomposición de términos. \(shelf\$facing\) corresponde a los cálculos sobre el efecto de \(facing\), lo que en nuestra descomposión de términos sería \(SSB_{J}\). Y Residuals corresponde a la variabilidad intragrupo (SSW).
Para medir el efecto de la altura en las ventas tenemos que comparar la varianza intergrupo de \(HEIGHT\) con la varianza intragrupo. Esto ya lo hicimos en el ANOVA de 1 vía. Nuestro estadístico es: \[ F_{HEIGHT} = \frac{SSB_{I}/(I-1)}{SSW/N-4} = \frac{316.3/2}{161.7/14}=\frac{158.2}{11.5}=13.7\] Para este caso que no se incluye la interacción, es más fácil ver los grados de libertad en la varianza intragrupo como:dimensión del problema - parámetros estimados. Tenemos 18 observaciones y se estiman 4 promedios (ocupamos 4 promedios): \(\bar{Y}_{...}\),\(\bar{Y}_{i..}\), \(\bar{Y}_{.j.}\) y \(\bar{Y}_{ij.}\). Por eso la tabla en la columna \(Df\) arroja el valor 14 en la fila residuals: 18 - 4 (total de nuestras observaciones menos los 4 promedios que debe estimar).
Para medir el efecto de la posición del producto (facing) en las ventas tenemos que comparar la varianza intergrupo de \(facing\) con la varianza intragrupo.
Nuestro estadístico es: \[ F_{facing} = \frac{SSB_{J}/(J-1)}{SSW/N-4} = \frac{480.5/1}{161.7/14}=\frac{480.5}{11.5}=41.61\] En ambos casos vemos que se rechaza la hipótesis nula a una significancia del 5%: existe una diferencia significativa entre las medias de ventas en los distintos niveles de \(HEIGHT\) y los distintos niveles de \(facing\).
Ahora, consideremos la interacción entre ambas variables. Para agregar la interacción a la fórmula se deben multiplicar las variables:
anova2_shelf <- aov(shelf$sales~shelf$HEIGHT + shelf$FACINGS + shelf$HEIGHT*shelf$FACINGS)
summary(anova2_shelf)## Df Sum Sq Mean Sq F value Pr(>F)
## shelf$HEIGHT 2 316.3 158.2 18.019 0.000243 ***
## shelf$FACINGS 1 480.5 480.5 54.741 8.29e-06 ***
## shelf$HEIGHT:shelf$FACINGS 2 56.3 28.2 3.209 0.076502 .
## Residuals 12 105.3 8.8
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Nuestro estadístico para la interacción es: \[ F_{height*facing} = \frac{SSB_{IJ}/(I-1)(J-1)}{SSW/IJ(K-1)} = \frac{56.3/2}{105.3/12}=\frac{28.2}{8.8}=3.209\] Para la interacción, debemos mirar los valores de la tabla que aparecen en la fila \(shelf\$HEIGHT:shelf\$FACINGS\).
Observamos que ambas categorías son significativas por si solas (a una significancia de 0,1%), y que la interacción también lo es pero a una significancia del 10%. Al 5% de significancia no podemos rechazar la hipótesis nula de que no existe un efecto en la interacción de \(HEIGHT\) y \(FACING\).