1 Ejercicio 1. Tabla de olor

1.1 1. Planteamiento y datos

Se compara el puntaje de olor de tres cafés, con seis repeticiones cada uno (18 observaciones en total). Los puntajes van de 1 a 5, siendo 5 el mejor. Se utiliza un diseño completamente al azar y un nivel de significancia de 0,05.

H0: las medias del puntaje de olor son iguales en los tres cafés (μ1 = μ2 = μ3).

H1: al menos un café tiene una media distinta.

Repetición | Café 1 | Café 2 | Café 3 |

|—:|—:|—:|—:|

1 | 2 | 5 | 1 |
2 | 4 | 3 | 3 |
3 | 3 | 4 | 3 |
4 | 2 | 3 | 3 |
5 | 4 | 3 | 3 |
6 | 5 | 3 | 3 |
olor <- c(2, 4, 3, 2, 4, 5,    # Café 1 

          5, 3, 4, 3, 3, 3,   # Café 2 

          1, 3, 3, 3, 3, 3)   # Café 3 

cafe <- factor(rep(c("Café 1", "Café 2", "Café 3"), each = 6), 

               levels = c("Café 1", "Café 2", "Café 3")) 

datos <- data.frame(cafe, olor) 

str(datos) 
## 'data.frame':    18 obs. of  2 variables:
##  $ cafe: Factor w/ 3 levels "Café 1","Café 2",..: 1 1 1 1 1 1 2 2 2 2 ...
##  $ olor: num  2 4 3 2 4 5 5 3 4 3 ...

1.2 2. Estadística descriptiva y gráficos

aggregate(olor ~ cafe, data = datos, 

          FUN = function(x) c(n = length(x), media = mean(x), 

                              mediana = median(x), de = sd(x))) 
boxplot(olor ~ cafe, data = datos, col = c("#F4B183", "#9DC3E6", "#A9D18E"), 

        ylab = "Puntaje de olor", xlab = "Café") 

stripchart(olor ~ cafe, data = datos, vertical = TRUE, method = "jitter", 

           pch = 16, ylab = "Puntaje de olor", xlab = "Café") 

points(1:3, tapply(datos$olor, datos$cafe, mean), pch = 95, cex = 3)

1.3 3. ANOVA y comparaciones múltiples

modelo <- aov(olor ~ cafe, data = datos) 

summary(modelo) 
##             Df Sum Sq Mean Sq F value Pr(>F)
## cafe         2  2.333  1.1667   1.235  0.319
## Residuals   15 14.167  0.9444
TukeyHSD(modelo) 
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = olor ~ cafe, data = datos)
## 
## $cafe
##                     diff       lwr       upr     p adj
## Café 2-Café 1  0.1666667 -1.290731 1.6240646 0.9526769
## Café 3-Café 1 -0.6666667 -2.124065 0.7907313 0.4778615
## Café 3-Café 2 -0.8333333 -2.290731 0.6240646 0.3254797
library(agricolae) 

LSD.test(modelo, "cafe", console = TRUE) 
## 
## Study: modelo ~ "cafe"
## 
## LSD t Test for olor 
## 
## Mean Square Error:  0.9444444 
## 
## cafe,  means and individual ( 95 %) CI
## 
##            olor       std r       se      LCL      UCL Min Max  Q25 Q50  Q75
## Café 1 3.333333 1.2110601 6 0.396746 2.487689 4.178977   2   5 2.25 3.5 4.00
## Café 2 3.500000 0.8366600 6 0.396746 2.654356 4.345644   3   5 3.00 3.0 3.75
## Café 3 2.666667 0.8164966 6 0.396746 1.821023 3.512311   1   3 3.00 3.0 3.00
## 
## Alpha: 0.05 ; DF Error: 15
## Critical Value of t: 2.13145 
## 
## least Significant Difference: 1.195921 
## 
## Treatments with the same letter are not significantly different.
## 
##            olor groups
## Café 2 3.500000      a
## Café 1 3.333333      a
## Café 3 2.666667      a

1.4 4. Verificación de supuestos del ANOVA

shapiro.test(residuals(modelo))                # normalidad 
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(modelo)
## W = 0.94489, p-value = 0.3505
bartlett.test(olor ~ cafe, data = datos)        # homogeneidad de varianzas 
## 
##  Bartlett test of homogeneity of variances
## 
## data:  olor by cafe
## Bartlett's K-squared = 0.9536, df = 2, p-value = 0.6208
qqnorm(residuals(modelo)); qqline(residuals(modelo)) 

plot(fitted(modelo), residuals(modelo)); abline(h = 0) 

# Prueba sobre los valores absolutos de los residuos 

anova(lm(abs(residuals(modelo)) ~ datos$cafe))

1.5 5. Prueba de Kruskal-Wallis

datos$rango <- rank(datos$olor)                # rank() promedia empates 

# automáticamente 

tapply(datos$rango, datos$cafe, sum)             # R_i 
## Café 1 Café 2 Café 3 
##   61.0   66.5   43.5
tapply(datos$rango, datos$cafe, mean)            # rango medio 
##   Café 1   Café 2   Café 3 
## 10.16667 11.08333  7.25000
kw <- kruskal.test(olor ~ cafe, data = datos) 

kw 
## 
##  Kruskal-Wallis rank sum test
## 
## data:  olor by cafe
## Kruskal-Wallis chi-squared = 2.0487, df = 2, p-value = 0.359
N <- nrow(datos) 

epsilon2 <- kw$statistic / (N - 1) 

epsilon2 
## Kruskal-Wallis chi-squared 
##                  0.1205096

1.6 6. Comparaciones múltiples: prueba de Dunn

Ejecute install.packages("FSA") una sola vez en la consola de R si aún no tiene instalado el paquete.

library(FSA) 

dunnTest(olor ~ cafe, data = datos, method = "holm") 
##        Comparison          Z   P.unadj     P.adj
## 1 Café 1 - Café 2 -0.3277254 0.7431193 0.7431193
## 2 Café 1 - Café 3  1.0427626 0.2970583 0.5941165
## 3 Café 2 - Café 3  1.3704880 0.1705346 0.5116039
dunn <- function(y, g, ajuste = "holm") { 

  r <- rank(y); N <- length(y) 

  Ri <- tapply(r, g, mean); ni <- tapply(r, g, length) 

  t <- table(r); C <- sum(t^3 - t) / (12 * (N - 1))    # corrección por empates 

  par <- combn(levels(g), 2) 

  z <- apply(par, 2, function(p) { 

    (Ri[p[1]] - Ri[p[2]]) / 

      sqrt((N * (N + 1) / 12 - C) * (1 / ni[p[1]] + 1 / ni[p[2]])) 

  }) 

  p <- 2 * pnorm(-abs(z)) 

  data.frame(Comparacion = paste(par[1, ], "-", par[2, ]), 

             Z = round(z, 4), p.sin.ajuste = signif(p, 4), 

             p.ajustado = signif(p.adjust(p, ajuste), 4)) 

} 

dunn(datos$olor, datos$cafe, "holm") 
pairwise.wilcox.test(datos$olor, datos$cafe, p.adjust.method = "bonferroni") 
## 
##  Pairwise comparisons using Wilcoxon rank sum exact test 
## 
## data:  datos$olor and datos$cafe 
## 
##        Café 1 Café 2
## Café 2 1.00   -     
## Café 3 1.00   0.82  
## 
## P value adjustment method: bonferroni

1.7 7. Redacción de resultados

El puntaje de olor se analiza mediante ANOVA de un factor y se contrasta con Kruskal-Wallis. Las comparaciones múltiples permiten revisar qué pares de cafés presentan diferencias. Interprete los valores p usando α = 0,05 y las salidas generadas por R.

2 Ejercicio 2. Tabla de sabor

2.1 1. Planteamiento y datos

Se compara el puntaje de sabor de tres cafés, con seis repeticiones cada uno (18 observaciones en total). Los puntajes van de 1 a 5, siendo 5 el mejor. Se utiliza un diseño completamente al azar y un nivel de significancia de 0,05.

H0: las medias del puntaje de sabor son iguales en los tres cafés (μ1 = μ2 = μ3).

H1: al menos un café tiene una media distinta.

Repetición | Café 1 | Café 2 | Café 3 |

|—:|—:|—:|—:|

1 | 4 | 5 | 3 |
2 | 4 | 2 | 2 |
3 | 3 | 4 | 2 |
4 | 2 | 4 | 2 |
5 | 4 | 2 | 3 |
6 | 4 | 2 | 2 |
sabor <- c(4, 4, 3, 2, 4, 4,   # Café 1 

           5, 2, 4, 4, 2, 2,   # Café 2 

           3, 2, 2, 2, 3, 2)   # Café 3 

cafe <- factor(rep(c("Café 1", "Café 2", "Café 3"), each = 6), 

               levels = c("Café 1", "Café 2", "Café 3")) 

datos <- data.frame(cafe, sabor) 

str(datos) 
## 'data.frame':    18 obs. of  2 variables:
##  $ cafe : Factor w/ 3 levels "Café 1","Café 2",..: 1 1 1 1 1 1 2 2 2 2 ...
##  $ sabor: num  4 4 3 2 4 4 5 2 4 4 ...

2.2 2. Estadística descriptiva y gráficos

aggregate(sabor ~ cafe, data = datos, 

          FUN = function(x) c(n = length(x), media = mean(x), 

                              mediana = median(x), de = sd(x))) 
boxplot(sabor ~ cafe, data = datos, col = c("#F4B183", "#9DC3E6", "#A9D18E"), 

        ylab = "Puntaje de sabor", xlab = "Café") 

stripchart(sabor ~ cafe, data = datos, vertical = TRUE, method = "jitter", 

           pch = 16, ylab = "Puntaje de sabor", xlab = "Café") 

points(1:3, tapply(datos$sabor, datos$cafe, mean), pch = 95, cex = 3)

2.3 3. ANOVA y comparaciones múltiples

modelo <- aov(sabor ~ cafe, data = datos) 

summary(modelo) 
##             Df Sum Sq Mean Sq F value Pr(>F)
## cafe         2  4.333  2.1667   2.378  0.127
## Residuals   15 13.667  0.9111
TukeyHSD(modelo) 
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = sabor ~ cafe, data = datos)
## 
## $cafe
##                     diff       lwr       upr     p adj
## Café 2-Café 1 -0.3333333 -1.764781 1.0981148 0.8196529
## Café 3-Café 1 -1.1666667 -2.598115 0.2647815 0.1198502
## Café 3-Café 2 -0.8333333 -2.264781 0.5981148 0.3133752
library(agricolae) 

LSD.test(modelo, "cafe", console = TRUE) 
## 
## Study: modelo ~ "cafe"
## 
## LSD t Test for sabor 
## 
## Mean Square Error:  0.9111111 
## 
## cafe,  means and individual ( 95 %) CI
## 
##           sabor       std r        se      LCL      UCL Min Max  Q25 Q50  Q75
## Café 1 3.500000 0.8366600 6 0.3896817 2.669413 4.330587   2   4 3.25   4 4.00
## Café 2 3.166667 1.3291601 6 0.3896817 2.336080 3.997254   2   5 2.00   3 4.00
## Café 3 2.333333 0.5163978 6 0.3896817 1.502746 3.163920   2   3 2.00   2 2.75
## 
## Alpha: 0.05 ; DF Error: 15
## Critical Value of t: 2.13145 
## 
## least Significant Difference: 1.174627 
## 
## Treatments with the same letter are not significantly different.
## 
##           sabor groups
## Café 1 3.500000      a
## Café 2 3.166667      a
## Café 3 2.333333      a

2.4 4. Verificación de supuestos del ANOVA

shapiro.test(residuals(modelo))                # normalidad 
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(modelo)
## W = 0.93289, p-value = 0.2182
bartlett.test(sabor ~ cafe, data = datos)       # homogeneidad de varianzas 
## 
##  Bartlett test of homogeneity of variances
## 
## data:  sabor by cafe
## Bartlett's K-squared = 3.8115, df = 2, p-value = 0.1487
qqnorm(residuals(modelo)); qqline(residuals(modelo)) 

plot(fitted(modelo), residuals(modelo)); abline(h = 0) 

# Prueba sobre los valores absolutos de los residuos 

anova(lm(abs(residuals(modelo)) ~ datos$cafe))

2.5 5. Prueba de Kruskal-Wallis

datos$rango <- rank(datos$sabor)                 # rank() promedia empates 

# automáticamente 

tapply(datos$rango, datos$cafe, sum)             # R_i 
## Café 1 Café 2 Café 3 
##   72.5   60.5   38.0
tapply(datos$rango, datos$cafe, mean)            # rango medio 
##    Café 1    Café 2    Café 3 
## 12.083333 10.083333  6.333333
kw <- kruskal.test(sabor ~ cafe, data = datos) 

kw 
## 
##  Kruskal-Wallis rank sum test
## 
## data:  sabor by cafe
## Kruskal-Wallis chi-squared = 4.1093, df = 2, p-value = 0.1281
N <- nrow(datos) 

epsilon2 <- kw$statistic / (N - 1) 

epsilon2 
## Kruskal-Wallis chi-squared 
##                  0.2417258

2.6 6. Comparaciones múltiples: prueba de Dunn

library(FSA) 

dunnTest(sabor ~ cafe, data = datos, method = "holm") 
##        Comparison         Z    P.unadj     P.adj
## 1 Café 1 - Café 2 0.6944563 0.48739612 0.4873961
## 2 Café 1 - Café 3 1.9965618 0.04587281 0.1376184
## 3 Café 2 - Café 3 1.3021055 0.19288032 0.3857606
dunn <- function(y, g, ajuste = "holm") { 

  r <- rank(y); N <- length(y) 

  Ri <- tapply(r, g, mean); ni <- tapply(r, g, length) 

  t <- table(r); C <- sum(t^3 - t) / (12 * (N - 1))    # corrección por empates 

  par <- combn(levels(g), 2) 

  z <- apply(par, 2, function(p) { 

    (Ri[p[1]] - Ri[p[2]]) / 

      sqrt((N * (N + 1) / 12 - C) * (1 / ni[p[1]] + 1 / ni[p[2]])) 

  }) 

  p <- 2 * pnorm(-abs(z)) 

  data.frame(Comparacion = paste(par[1, ], "-", par[2, ]), 

             Z = round(z, 4), p.sin.ajuste = signif(p, 4), 

             p.ajustado = signif(p.adjust(p, ajuste), 4)) 

} 

dunn(datos$sabor, datos$cafe, "holm") 
pairwise.wilcox.test(datos$sabor, datos$cafe, p.adjust.method = "bonferroni") 
## 
##  Pairwise comparisons using Wilcoxon rank sum exact test 
## 
## data:  datos$sabor and datos$cafe 
## 
##        Café 1 Café 2
## Café 2 1.00   -     
## Café 3 0.14   0.95  
## 
## P value adjustment method: bonferroni

2.7 7. Resultados

El puntaje de sabor se analiza mediante ANOVA de un factor y se contrasta con Kruskal-Wallis. Las comparaciones múltiples permiten revisar qué pares de cafés presentan diferencias. Interprete los valores p usando α = 0,05 y las salidas generadas por R.

3 Conclusión general

Compare los resultados obtenidos para olor y sabor. Recuerde que un valor p menor que 0,05 indica evidencia estadística de diferencias entre los cafés; un valor p igual o mayor que 0,05 no permite rechazar la hipótesis nula.