Pruebas de rangos múltiples
Posterior al rechazo de la hipótesis nula \(H_o=\mu_1=\mu_2=...=\mu_a\) en el análisis de varianza, ANOVA, para un experimento de un factor. Es necesario ir a detalle y ver ¿cuáles tratamientos son estadísticamente diferentes? esta interrogante se responden probando la igualdad de todos los posibles pares de medias, para lo que se han propuesto varios métodos, conocidos como métodos de comparaciones múltiples o pruebas de rango múltiple.
Prueba LSD (Least Significant Difference) - Diferencia mínima significativa.
1. Hipótesis
Las hipótesis a plantear corresponden a:
\[H_o: \mu_i = \mu_j \\ H_1 : u_i \neq u_j \\ \forall~i \neq j ,\\ i,j = 1, 2,...,a \]
Para \(a\) tratamientos se tienen en total \(\frac{a(a-1)}{2}\) pares de medias, o comparaciones pareadas a realizar.
2. Cálculo de \(LSD\)
La diferencia mínima significativa \(LSD\) es igual a:
\[LSD = t_{\frac{\alpha}{2},~N-a} \sqrt{\frac{2MS_{error}}{n}}\]
3. Comparación \(LSD\) con valor absoluto de diferencia de medias puntuales \(|\bar{y_{i.}}-\bar{y_{j.}}|\)
\[Si,~|\bar{y_{i.}}-\bar{y_{j.}}|>LSD ~ \rightarrow Rechazo~H_o\]
Ejemplo en \(RStudio\)
Vamos a tomar como ejemplo el caso de la fibra sintética fabricada con distintos porcentajes de algodón que se presume afecta la resistencia a la tensión:
| Peso porcentual algodón | 1 | 2 | 3 | 4 | 5 | Total | Promedio |
|---|---|---|---|---|---|---|---|
| 15% | 7 | 7 | 15 | 11 | 9 | 49 | 9.80 |
| 20% | 12 | 17 | 12 | 18 | 18 | 77 | 15.40 |
| 25% | 14 | 18 | 18 | 19 | 19 | 88 | 17.60 |
| 30% | 19 | 25 | 22 | 19 | 23 | 108 | 21.60 |
| 35% | 7 | 10 | 11 | 15 | 11 | 54 | 10.80 |
| 376 | 15.04 |
Para el ejemplo realizamos ANOVA
# Leemos los datos desde un archivo de excel creando un dataframe llamado "Datos_base"
library(readxl)
Datos_base <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2024__1/Diseno_de_Experimentos_11951/Clases/7_Pruebas_Rangos_Multiples/2024_1/Datos_base_1.xlsx")
#Definimos el factor "Porcentaje de algodón" y planteamos el modelos de datos
Datos_base$Porc_Algodon<-as.factor(Datos_base$Porc_Algodon)
modelo <- lm(Datos_base$Resistencia~Datos_base$Porc_Algodon)
# Realizamos ANOVA
anova <- aov(modelo)
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## Datos_base$Porc_Algodon 4 475.8 118.94 14.76 9.13e-06 ***
## Residuals 20 161.2 8.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Prueba LSD
library (agricolae)
Grupos <- LSD.test(y=anova, trt="Datos_base$Porc_Algodon", group=TRUE, console=TRUE)##
## Study: anova ~ "Datos_base$Porc_Algodon"
##
## LSD t Test for Datos_base$Resistencia
##
## Mean Square Error: 8.06
##
## Datos_base$Porc_Algodon, means and individual ( 95 %) CI
##
## Datos_base.Resistencia std r se LCL UCL
## Tratamiento_15 9.8 3.346640 5 1.269646 7.151566 12.44843
## Tratamiento_20 15.4 3.130495 5 1.269646 12.751566 18.04843
## Tratamiento_25 17.6 2.073644 5 1.269646 14.951566 20.24843
## Tratamiento_30 21.6 2.607681 5 1.269646 18.951566 24.24843
## Tratamiento_35 10.8 2.863564 5 1.269646 8.151566 13.44843
## Min Max Q25 Q50 Q75
## Tratamiento_15 7 15 7 9 11
## Tratamiento_20 12 18 12 17 18
## Tratamiento_25 14 19 18 18 19
## Tratamiento_30 19 25 19 22 23
## Tratamiento_35 7 15 10 11 11
##
## Alpha: 0.05 ; DF Error: 20
## Critical Value of t: 2.085963
##
## least Significant Difference: 3.745452
##
## Treatments with the same letter are not significantly different.
##
## Datos_base$Resistencia groups
## Tratamiento_30 21.6 a
## Tratamiento_25 17.6 b
## Tratamiento_20 15.4 b
## Tratamiento_35 10.8 c
## Tratamiento_15 9.8 cGrupos## $statistics
## MSerror Df Mean CV t.value LSD
## 8.06 20 15.04 18.87642 2.085963 3.745452
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD none Datos_base$Porc_Algodon 5 0.05
##
## $means
## Datos_base$Resistencia std r se LCL UCL
## Tratamiento_15 9.8 3.346640 5 1.269646 7.151566 12.44843
## Tratamiento_20 15.4 3.130495 5 1.269646 12.751566 18.04843
## Tratamiento_25 17.6 2.073644 5 1.269646 14.951566 20.24843
## Tratamiento_30 21.6 2.607681 5 1.269646 18.951566 24.24843
## Tratamiento_35 10.8 2.863564 5 1.269646 8.151566 13.44843
## Min Max Q25 Q50 Q75
## Tratamiento_15 7 15 7 9 11
## Tratamiento_20 12 18 12 17 18
## Tratamiento_25 14 19 18 18 19
## Tratamiento_30 19 25 19 22 23
## Tratamiento_35 7 15 10 11 11
##
## $comparison
## NULL
##
## $groups
## Datos_base$Resistencia groups
## Tratamiento_30 21.6 a
## Tratamiento_25 17.6 b
## Tratamiento_20 15.4 b
## Tratamiento_35 10.8 c
## Tratamiento_15 9.8 c
##
## attr(,"class")
## [1] "group"Los tratamientos pertenecientes al mismo grupo suponen medias estadísticamente iguales, los tratamientos en grupos diferentes suponen medias estadísticamente diferentes
Ejemplo de clase resuelto en RStudio - Test LSD
Durante el desarrollo de las clases de diseño de experimento se ha utilizado el siguiente problema:
Un fabricante de calzado desea mejorar la calidad de las suelas, las cuales se pueden hacer con uno de los cuatro tipos de cuero A, B, C y D disponibles en el mercado. Para ello, prueba los cueros con una máquina que hace pasar los zapatos por una superficie abrasiva; la suela de éstos se desgasta al pasarla por dicha superficie. Como criterio de desgaste se usa la pérdida de peso después de un número fijo de ciclos. Se prueban en orden aleatorio 24 zapatos, seis de cada tipo de cuero. Al hacer las pruebas en orden completamente al azar se evitan sesgos y las mediciones en un tipo de cuero resultan independientes de las demás. Los datos (en miligramos) sobre el desgaste de cada tipo de cuero se muestran en la tabla siguiente:
| Tipo de cuero | 1 | 2 | 3 | 4 | 5 | 6 |
|---|---|---|---|---|---|---|
| A | 264 | 260 | 258 | 241 | 262 | 255 |
| B | 208 | 220 | 216 | 200 | 213 | 206 |
| C | 220 | 263 | 219 | 225 | 230 | 228 |
| D | 217 | 226 | 215 | 227 | 220 | 222 |
El procedimiento para solucionarlo en R, es el siguiente:
# Realizamos análisis de varianza
cuero <- c("A", "A","A","A","A","A", "B", "B","B","B","B","B","C","C","C","C","C","C","D","D","D","D","D","D")
desgaste <- c(264,260,258,241,262,255,208,220,216,200,213,206,220,263,219,225,230,228,217,226,215,227,220,222)
cuero <- as.factor(cuero)
modelo <- lm(desgaste ~ cuero)
anova <- aov(modelo)
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## cuero 3 7019 2339.8 22.75 1.18e-06 ***
## Residuals 20 2056 102.8
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova$fitted.values## 1 2 3 4 5 6 7 8
## 256.6667 256.6667 256.6667 256.6667 256.6667 256.6667 210.5000 210.5000
## 9 10 11 12 13 14 15 16
## 210.5000 210.5000 210.5000 210.5000 230.8333 230.8333 230.8333 230.8333
## 17 18 19 20 21 22 23 24
## 230.8333 230.8333 221.1667 221.1667 221.1667 221.1667 221.1667 221.1667# Realizamos prueba LSD
library(agricolae)
grupos <- LSD.test(y=anova, trt = "cuero", group=TRUE)
grupos## $statistics
## MSerror Df Mean CV t.value LSD
## 102.825 20 229.7917 4.412809 2.085963 12.21224
##
## $parameters
## test p.ajusted name.t ntr alpha
## Fisher-LSD none cuero 4 0.05
##
## $means
## desgaste std r se LCL UCL Min Max Q25 Q50 Q75
## A 256.6667 8.286535 6 4.139746 248.0313 265.3020 241 264 255.75 259.0 261.50
## B 210.5000 7.259477 6 4.139746 201.8646 219.1354 200 220 206.50 210.5 215.25
## C 230.8333 16.339115 6 4.139746 222.1980 239.4687 219 263 221.25 226.5 229.50
## D 221.1667 4.792355 6 4.139746 212.5313 229.8020 215 227 217.75 221.0 225.00
##
## $comparison
## NULL
##
## $groups
## desgaste groups
## A 256.6667 a
## C 230.8333 b
## D 221.1667 bc
## B 210.5000 c
##
## attr(,"class")
## [1] "group"Los tratamientos pertenecientes al mismo grupo suponen medias estadísticamente iguales, los tratamientos en grupos diferentes suponen medias estadísticamente diferentes
Prueba de Tukey
1. Hipótesis
Las hipótesis a plantear corresponden a:
\[H_o: \mu_i = \mu_j \\ H_1 : u_i \neq u_j \\ \forall~i \neq j ,\\ i,j = 1, 2,...,a \]
Para \(a\) tratamientos se tienen en total \(\frac{a(a-1)}{2}\) pares de medias, o comparaciones pareadas a realizar.
2. Cálculo de \(T_\alpha\)
El valor crítico \(T_\alpha\) se calcula de la siguiente manera:
\[T_\alpha = q_\alpha (a, N-a) \sqrt{\frac{MS_{error}}{n}}\]
3. Comparación \(T_\alpha\) con valor absoluto de diferencia de medias puntuales \(|\bar{y_{i.}}-\bar{y_{j.}}|\)
\[Si,~|\bar{y_{i.}}-\bar{y_{j.}}|>T_\alpha ~ \rightarrow Rechazo~H_o\]
Ejemplo en \(RStudio\)
Vamos a tomar como ejemplo el caso de la fibra sintética fabricada con distintos porcentajes de algodón que se presume afecta la resistencia a la tensión:
| Peso porcentual algodón | 1 | 2 | 3 | 4 | 5 | Total | Promedio |
|---|---|---|---|---|---|---|---|
| 15% | 7 | 7 | 15 | 11 | 9 | 49 | 9.80 |
| 20% | 12 | 17 | 12 | 18 | 18 | 77 | 15.40 |
| 25% | 14 | 18 | 18 | 19 | 19 | 88 | 17.60 |
| 30% | 19 | 25 | 22 | 19 | 23 | 108 | 21.60 |
| 35% | 7 | 10 | 11 | 15 | 11 | 54 | 10.80 |
| 376 | 15.04 |
Para el ejemplo realizamos ANOVA
# Leemos los datos desde un archivo de excel creando un dataframe llamado "Datos_base"
library(readxl)
Datos_base <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2024__1/Diseno_de_Experimentos_11951/Clases/7_Pruebas_Rangos_Multiples/2024_1/Datos_base_1.xlsx")
#Definimos el factor "Porcentaje de algodón" y planteamos el modelos de datos
Datos_base$Porc_Algodon<-as.factor(Datos_base$Porc_Algodon)
modelo <- lm(Datos_base$Resistencia~Datos_base$Porc_Algodon)
# Realizamos ANOVA
anova <- aov(modelo)
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## Datos_base$Porc_Algodon 4 475.8 118.94 14.76 9.13e-06 ***
## Residuals 20 161.2 8.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Prueba Tukey
library (agricolae)
prueba_tukey <- TukeyHSD(anova, conf.level = 0.95)
prueba_tukey## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = modelo)
##
## $`Datos_base$Porc_Algodon`
## diff lwr upr p adj
## Tratamiento_20-Tratamiento_15 5.6 0.2270417 10.9729583 0.0385024
## Tratamiento_25-Tratamiento_15 7.8 2.4270417 13.1729583 0.0025948
## Tratamiento_30-Tratamiento_15 11.8 6.4270417 17.1729583 0.0000190
## Tratamiento_35-Tratamiento_15 1.0 -4.3729583 6.3729583 0.9797709
## Tratamiento_25-Tratamiento_20 2.2 -3.1729583 7.5729583 0.7372438
## Tratamiento_30-Tratamiento_20 6.2 0.8270417 11.5729583 0.0188936
## Tratamiento_35-Tratamiento_20 -4.6 -9.9729583 0.7729583 0.1162970
## Tratamiento_30-Tratamiento_25 4.0 -1.3729583 9.3729583 0.2101089
## Tratamiento_35-Tratamiento_25 -6.8 -12.1729583 -1.4270417 0.0090646
## Tratamiento_35-Tratamiento_30 -10.8 -16.1729583 -5.4270417 0.0000624En este caso comparamos e \(diff\) que corresponde al \(T_\alpha\) con cada una de las diferencias entre medias o se compara el \(p-adj\) con el nivel de significancia \(\alpha\) para rechazar, o no, las hipótesis \(H_o\)
Ejemplo de clases resuelto en RStudio - Test de Tukey
Durante el desarrollo de las clases de diseño de experimento se ha utilizado el siguiente problema:
Un fabricante de calzado desea mejorar la calidad de las suelas, las cuales se pueden hacer con uno de los cuatro tipos de cuero A, B, C y D disponibles en el mercado. Para ello, prueba los cueros con una máquina que hace pasar los zapatos por una superficie abrasiva; la suela de éstos se desgasta al pasarla por dicha superficie. Como criterio de desgaste se usa la pérdida de peso después de un número fijo de ciclos. Se prueban en orden aleatorio 24 zapatos, seis de cada tipo de cuero. Al hacer las pruebas en orden completamente al azar se evitan sesgos y las mediciones en un tipo de cuero resultan independientes de las demás. Los datos (en miligramos) sobre el desgaste de cada tipo de cuero se muestran en la tabla siguiente:
| Tipo de cuero | 1 | 2 | 3 | 4 | 5 | 6 |
|---|---|---|---|---|---|---|
| A | 264 | 260 | 258 | 241 | 262 | 255 |
| B | 208 | 220 | 216 | 200 | 213 | 206 |
| C | 220 | 263 | 219 | 225 | 230 | 228 |
| D | 217 | 226 | 215 | 227 | 220 | 222 |
El procedimiento para solucionarlo en R, es el siguiente:
# Realizamos análisis de varianza
cuero <- c("A", "A","A","A","A","A", "B", "B","B","B","B","B","C","C","C","C","C","C","D","D","D","D","D","D")
desgaste <- c(264,260,258,241,262,255,208,220,216,200,213,206,220,263,219,225,230,228,217,226,215,227,220,222)
cuero <- as.factor(cuero)
modelo <- lm(desgaste ~ cuero)
anova <- aov(modelo)
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## cuero 3 7019 2339.8 22.75 1.18e-06 ***
## Residuals 20 2056 102.8
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova$fitted.values## 1 2 3 4 5 6 7 8
## 256.6667 256.6667 256.6667 256.6667 256.6667 256.6667 210.5000 210.5000
## 9 10 11 12 13 14 15 16
## 210.5000 210.5000 210.5000 210.5000 230.8333 230.8333 230.8333 230.8333
## 17 18 19 20 21 22 23 24
## 230.8333 230.8333 221.1667 221.1667 221.1667 221.1667 221.1667 221.1667# Realizamos prueba LSD
library(agricolae)
pruebatukey <- TukeyHSD(anova, conf.level = 0.95)
pruebatukey## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = modelo)
##
## $cuero
## diff lwr upr p adj
## B-A -46.166667 -62.552998 -29.780336 0.0000008
## C-A -25.833333 -42.219664 -9.447002 0.0014117
## D-A -35.500000 -51.886331 -19.113669 0.0000349
## C-B 20.333333 3.947002 36.719664 0.0118160
## D-B 10.666667 -5.719664 27.052998 0.2926431
## D-C -9.666667 -26.052998 6.719664 0.3742863Prueba de Duncan
Las hipótesis a plantear corresponden a:
\[H_o: \mu_i = \mu_j \\ H_1 : u_i \neq u_j \\ \forall~i \neq j ,\\ i,j = 1, 2,...,a \]
Para \(a\) tratamientos se tienen en total \(\frac{a(a-1)}{2}\) pares de medias, o comparaciones pareadas a realizar.
2. Calculo de error estándar \(S_{\bar{y_{i.}}}\)
El error estándar \(S_{\bar{y_{i.}}}\) se calcula como sigue:
\[S_{\bar{y_{i.}}} = \sqrt{\frac{MS_{error}}{n}}\]
La anterior ecuación aplica para tamaños de muestras iguales, esto es, un diseño de experimentos de un factor con iguales réplicas por cada tratamiento.
Para tamaños de muestras distintos se sustituye \(n\) por la media armónica \(n_h\) que se calcula de la siguiente manera:
\[n_h=\frac{a}{\sum_{i=1}^a \frac{1}{n_i}}\]
3. Cálculo de rangos de significancia mínimos \(R_p\)
Los rangos de significancia mínimos \(R_p\) se calculan como sigue:
\[R_p=r_{\alpha}(p,f)S_{\bar{y_{i.}}}~,~p=2,3,...,a\\ Donde,~r_\alpha (p,f): rangos~significativos~de~tabla~Duncan \]
4. Organización de estimadores de medias de tratamientos y comparación de las diferencias entre estimadores de las medias organizadas y los \(R_p\).
Las estimaciones de las medias puntuales \(\bar{y_{i.}}\) se organizan de manera ascendente.
Posterior a la organización de las estimaciones de las medias puntuales, se calculan las diferencias entre medias, empezando con la más grande contra la menor, la cual se compararía con el rango mínimo de significación \(Ra\). Después se calcula la diferencia de la mayor y la segunda menor y se compara con el rango mínimo de significación \(R_a-1\). Estas comparaciones se continúan hasta que todas las medias se han comparado con la media mayor. Por último, se calcula la diferencia entre la segunda media mayor y la menor y se compara con el rango mínimo de significación \(R_a-1\). Este proceso se continúa hasta que se han considerado las diferencias entre todos los \(\frac{a(a - 1)}{2}\) pares de medias posibles.
Si una diferencia observada es mayor al rango de significancia mínima correspondiente se concluye que el par de medias en cuestión son diferentes, esto es:
\[Si,~\bar{y_{i.}}-\bar{y_{j.}}>R_p \rightarrow ~ Rechazo~ H_o\]
Ejemplo en \(RStudio\)
Vamos a tomar como ejemplo el caso de la fibra sintética fabricada con distintos porcentajes de algodón que se presume afecta la resistencia a la tensión:
| Peso porcentual algodón | 1 | 2 | 3 | 4 | 5 | Total | Promedio |
|---|---|---|---|---|---|---|---|
| 15% | 7 | 7 | 15 | 11 | 9 | 49 | 9.80 |
| 20% | 12 | 17 | 12 | 18 | 18 | 77 | 15.40 |
| 25% | 14 | 18 | 18 | 19 | 19 | 88 | 17.60 |
| 30% | 19 | 25 | 22 | 19 | 23 | 108 | 21.60 |
| 35% | 7 | 10 | 11 | 15 | 11 | 54 | 10.80 |
| 376 | 15.04 |
Para el ejemplo realizamos ANOVA
# Leemos los datos desde un archivo de excel creando un dataframe llamado "Datos_base"
library(readxl)
Datos_base <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2024__1/Diseno_de_Experimentos_11951/Clases/7_Pruebas_Rangos_Multiples/Datos_base_1.xlsx")
#Definimos el factor "Porcentaje de algodón" y planteamos el modelos de datos
Datos_base$Porc_Algodon<-as.factor(Datos_base$Porc_Algodon)
modelo <- lm(Datos_base$Resistencia~Datos_base$Porc_Algodon)
# Realizamos ANOVA
anova <- aov(modelo)
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## Datos_base$Porc_Algodon 4 475.8 118.94 14.76 9.13e-06 ***
## Residuals 20 161.2 8.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Prueba Duncan
library (agricolae)
prueba_duncan <- duncan.test(anova,"Datos_base$Porc_Algodon",alpha=0.05,group=TRUE)
prueba_duncan## $statistics
## MSerror Df Mean CV
## 8.06 20 15.04 18.87642
##
## $parameters
## test name.t ntr alpha
## Duncan Datos_base$Porc_Algodon 5 0.05
##
## $duncan
## Table CriticalRange
## 2 2.949998 3.745452
## 3 3.096506 3.931466
## 4 3.189616 4.049682
## 5 3.254648 4.132249
##
## $means
## Datos_base$Resistencia std r se Min Max Q25 Q50 Q75
## Tratamiento_15 9.8 3.346640 5 1.269646 7 15 7 9 11
## Tratamiento_20 15.4 3.130495 5 1.269646 12 18 12 17 18
## Tratamiento_25 17.6 2.073644 5 1.269646 14 19 18 18 19
## Tratamiento_30 21.6 2.607681 5 1.269646 19 25 19 22 23
## Tratamiento_35 10.8 2.863564 5 1.269646 7 15 10 11 11
##
## $comparison
## NULL
##
## $groups
## Datos_base$Resistencia groups
## Tratamiento_30 21.6 a
## Tratamiento_25 17.6 b
## Tratamiento_20 15.4 b
## Tratamiento_35 10.8 c
## Tratamiento_15 9.8 c
##
## attr(,"class")
## [1] "group"Los tratamientos pertenecientes al mismo grupo suponen medias estadísticamente iguales, los tratamientos en grupos diferentes suponen medias estadísticamente diferentes