El rendimiento de un proceso químico se midió utilizando cinco lotes de materia prima, cinco concentraciones de ácido, cinco tiempos de procesamiento (A,B,C,D y E) y cinco concentraciones de catalizador (\(\alpha, \beta, \gamma, \delta, \epsilon\)). Se usó el cuadrado grecolatino siguiente. Analizar los datos de este experimento (utilizar \(\alpha =0.05\)) y sacar concluciones.
Nuestras hipótesis a considerar en este estudio seran las siguientes:
La hipótesis nula será que sus medias son iguales.
La hipótesis alternativa será que al menos una de ellas difiere.
Datos del estudio:
Lote_mat_pri<-factor(c(1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5))
Concen_d_aci<-factor(c(1,1,1,1,1,2,2,2,2,2,3,3,3,3,3,4,4,4,4,4,5,5,5,5,5))
tiem_proce<-factor(c("A","B","C","D","E",
"B","C","D","E","A",
"C","D","E","A","B",
"D","E","A","B","C",
"E","A","B","C","D"))
concen_d_cata<-factor(c("a","c","e","b","d",
"b","d","a","c","e",
"c","e","b","d","a",
"d","a","c","e","b",
"e","b","d","a","c"))
proce_quim<-c(26,18,20,15,10,
16,21,12,15,24,
19,18,16,22,17,
16,11,25,14,17,
13,21,13,17,14)
data = data.frame(Lote_mat_pri,Concen_d_aci,concen_d_cata,tiem_proce,proce_quim)
data
## Lote_mat_pri Concen_d_aci concen_d_cata tiem_proce proce_quim
## 1 1 1 a A 26
## 2 2 1 c B 18
## 3 3 1 e C 20
## 4 4 1 b D 15
## 5 5 1 d E 10
## 6 1 2 b B 16
## 7 2 2 d C 21
## 8 3 2 a D 12
## 9 4 2 c E 15
## 10 5 2 e A 24
## 11 1 3 c C 19
## 12 2 3 e D 18
## 13 3 3 b E 16
## 14 4 3 d A 22
## 15 5 3 a B 17
## 16 1 4 d D 16
## 17 2 4 a E 11
## 18 3 4 c A 25
## 19 4 4 e B 14
## 20 5 4 b C 17
## 21 1 5 e E 13
## 22 2 5 b A 21
## 23 3 5 d B 13
## 24 4 5 a C 17
## 25 5 5 c D 14
Mostraremos la tabla de ANOVA correspondiente a los datos anteriores
modelo<-lm(proce_quim~Lote_mat_pri+Concen_d_aci+concen_d_cata+tiem_proce)
anova<-aov(modelo)
summary(anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Lote_mat_pri 4 10.0 2.50 0.427 0.785447
## Concen_d_aci 4 24.4 6.10 1.043 0.442543
## concen_d_cata 4 12.0 3.00 0.513 0.728900
## tiem_proce 4 342.8 85.70 14.650 0.000941 ***
## Residuals 8 46.8 5.85
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Teniendo en cuenta los datos obtenidos de la tabla de ANOVA, podemos hacer las siguientes afirmaciones:
Teniendo en cuenta que el valor p (del lote de materia prima, concentracionde acido, concentracion de catalizador) obtenido en la tabla de ANOVA es mayor que el \(\alpha\) considerado en el problema, no hay suficiente evidencia para rechazar la hipótesis nula, por lo tanto consideraremos que sus medias son iguales.
Teniendo en cuenta que el valor p (del tiempo de procesamiento) obtenido en la tabla de ANOVA es menor que el \(\alpha\) considerado en el problema, entonces rechazamos la hipótesis nula y consideramos que al menos una de sus medias difiere.
En los tiempos de procesamiento nos dice que con respecto a la tabla de ANOVA al menos una de sus medias no son iguales por lo tanto analizamos las diferencias por medio del procedimiento LSD
library(agricolae)
lsd_modelo<-LSD.test(y=anova, trt = "tiem_proce", group = T,console = T)
##
## Study: anova ~ "tiem_proce"
##
## LSD t Test for proce_quim
##
## Mean Square Error: 5.85
##
## tiem_proce, means and individual ( 95 %) CI
##
## proce_quim std r LCL UCL Min Max
## A 23.6 2.073644 5 21.10568 26.09432 21 26
## B 15.6 2.073644 5 13.10568 18.09432 13 18
## C 18.8 1.788854 5 16.30568 21.29432 17 21
## D 15.0 2.236068 5 12.50568 17.49432 12 18
## E 13.0 2.549510 5 10.50568 15.49432 10 16
##
## Alpha: 0.05 ; DF Error: 8
## Critical Value of t: 2.306004
##
## least Significant Difference: 3.527508
##
## Treatments with the same letter are not significantly different.
##
## proce_quim groups
## A 23.6 a
## C 18.8 b
## B 15.6 bc
## D 15.0 c
## E 13.0 c
bar.group(x=lsd_modelo$groups, horiz = T,col="red",
xlab="Rendimiento del proceso químico",
ylab="Tiempo de procesamiento",
xlim=c(0,28),
main="Comparación de los Tiempos de procesamiento\n por medio del procedimiento LSD")
De la anterior grafica podemos hacer las siguientes concluciones
A continuación, analizaremos si existe una diferencia en tiempo de procesamiento por pares. Mostraremos una tabla de datos y un gráfico para ilustrar.
modelo2<-lm(proce_quim~tiem_proce)
anova2<-aov(modelo2)
summary(anova2)
## Df Sum Sq Mean Sq F value Pr(>F)
## tiem_proce 4 342.8 85.70 18.39 1.77e-06 ***
## Residuals 20 93.2 4.66
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
intervals = TukeyHSD(anova2)
intervals
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = modelo2)
##
## $tiem_proce
## diff lwr upr p adj
## B-A -8.0 -12.0854407 -3.9145593 0.0000875
## C-A -4.8 -8.8854407 -0.7145593 0.0164823
## D-A -8.6 -12.6854407 -4.5145593 0.0000339
## E-A -10.6 -14.6854407 -6.5145593 0.0000017
## C-B 3.2 -0.8854407 7.2854407 0.1723691
## D-B -0.6 -4.6854407 3.4854407 0.9916448
## E-B -2.6 -6.6854407 1.4854407 0.3469798
## D-C -3.8 -7.8854407 0.2854407 0.0760762
## E-C -5.8 -9.8854407 -1.7145593 0.0032224
## E-D -2.0 -6.0854407 2.0854407 0.5954156
plot(intervals)
De los datos obtenidos anteriormente, se puede concluir que el grupo A es el que difiere más de los demás grupos. Dependiendo del objetivo del estudio, se podría considerar que el grupo A es el que genera más ruido en el análisis y por lo tanto, debería ser retirado del estudio para poder realizar un análisis más preciso. Sin embargo, si se quisiera analizar el grupo con mayor tiempo de procesamiento, entonces la mejor opción sería mantener el grupo A en el estudio.
A partir de los residuos del modelo, comprobaremos si el modelo ANOVA es adecuado. Los tres supuestos que se deben cumplir son: independencia, homocedasticidad y normalidad.
Probaremos la independencia con la siguiente gráfica:
plot(anova$residuals,main = "Residuos vs Observaciones")
Como se puede observar al ver la gráfica, se puede apreciar que los datos están dispersos y no tienen una forma específica. Por lo tanto, podemos afirmar que son independientes.
Para probar la normalidad en el estudio, lo haremos a través de una gráfica y utilizando el test de Shapiro-Wilk
hist(anova$residuals, col="lightgreen", main = "Histograma de los Residuos",
freq = F, xlab="Residuos",ylab="Densidad")
lines(density(anova$residuals), col="red", lwd=3)
Como podemos observar de la figura anterior, el histograma intenta tomar forma de campana de Gauss, lo cual es una indicación de que hay normalidad en los datos.
Si utilizamos el test de Shapiro-Wilk, tendremos lo siguiente:
Nuestras hipótesis a considerar en la validación de la normalidad son los siguientes:
La hipótesis nula será que tiene una distribucion normal.
La hipótesis alternativa será que no sigue una distribucion normal.
shapiro.test(anova$res)
##
## Shapiro-Wilk normality test
##
## data: anova$res
## W = 0.94721, p-value = 0.2167
Con los resultados obtenidos del test, se puede observar que el valor p es de 0.2167, lo cual es mayor que el nivel de significancia (\(\alpha\)) establecido. En consecuencia, no hay suficiente evidencia para rechazar la hipótesis nula, por lo tanto, se puede afirmar que existe normalidad en los datos.
Para demostrar si las varianzas son iguales o difieren, utilizaremos el test de Bartlett, ya que estamos en una distribución que es normal.
Nuestras hipótesis a considerar en la validación de la homocedasticidad son los siguientes:
La hipótesis nula será que sus varianzas son iguales.
La hipótesis alternativa será que al menos una de ellas difiere.
bartlett.test(anova$residuals ~ Lote_mat_pri)
##
## Bartlett test of homogeneity of variances
##
## data: anova$residuals by Lote_mat_pri
## Bartlett's K-squared = 1.1088, df = 4, p-value = 0.8929
bartlett.test(anova$residuals ~ concen_d_cata)
##
## Bartlett test of homogeneity of variances
##
## data: anova$residuals by concen_d_cata
## Bartlett's K-squared = 4.7772, df = 4, p-value = 0.3109
bartlett.test(anova$residuals ~ Concen_d_aci)
##
## Bartlett test of homogeneity of variances
##
## data: anova$residuals by Concen_d_aci
## Bartlett's K-squared = 2.9543, df = 4, p-value = 0.5655
bartlett.test(anova$residuals ~ tiem_proce)
##
## Bartlett test of homogeneity of variances
##
## data: anova$residuals by tiem_proce
## Bartlett's K-squared = 1.423, df = 4, p-value = 0.8402
Teniendo en cuenta que el valor p obtenido en el test de bartlett entonces para el lote de materia prima, la concentración de ácido y la concentración de catalizador es mayor que el nivel de significancia \(\alpha\) establecido para el problema, no hay suficiente evidencia para rechazar la hipótesis nula. Por lo tanto, se considera que las varianzas de dichas variables son iguales.