Diseños factoriales
Omar Jhondely Méndez Juárez
9/8/2021
DISEÑO FACTORIAL \(2^{5-1}\)
El producto químico
Considere un experimiento \(2^{5-1}\) con I=ABCDE que fue utilidado para investigar los efectos de cinco factores sobre el el color de un producto químico. Los factores son A=solvente/reactante, B=catalizador/reactante, C=temperatura, D=pureza del reactante y E=acidez del reactante. Los resultados obtenidos son los siguientes:
combinaciones
a) Calcule los efectos y grafíquelos en Pareto y en papel normal ¿Cuáles parecen significativos?
b)Obtenga el mejor análisis de varianza.¿Con cuáles efectos se está contruyendo el error?
c)Represente gráficamente cada efecto significativo e interpretelo a detalle.
d) Determine el mejor tratamiento y la respuesta predicha por el modelo
e) Haga un analisis de residuos y comente los resultados
a)
Lo primero que se debe realizar es mandar a instruir las librerías que se van a utilizar para poder realizar nuestro diseño factorial, tenemos la librería “readx” que llama los datos de excel, y la librería “FrF2” que hará nuestro analísis estadístico.
library(readxl)
library(FrF2)## Loading required package: DoE.base## Loading required package: grid## Loading required package: conf.design## Registered S3 method overwritten by 'DoE.base':
## method from
## factorize.factor conf.design##
## Attaching package: 'DoE.base'## The following objects are masked from 'package:stats':
##
## aov, lm## The following object is masked from 'package:graphics':
##
## plot.design## The following object is masked from 'package:base':
##
## lengthsdatos=read_excel(path = "dataset.xlsx")
View(datos)
attach(datos)
str(datos)## tibble [16 x 6] (S3: tbl_df/tbl/data.frame)
## $ A : num [1:16] -1 1 -1 1 -1 1 -1 1 -1 1 ...
## $ B : num [1:16] -1 -1 1 1 -1 -1 1 1 -1 -1 ...
## $ C : num [1:16] -1 -1 -1 -1 1 1 1 1 -1 -1 ...
## $ D : num [1:16] -1 -1 -1 -1 -1 -1 -1 -1 1 1 ...
## $ E : num [1:16] 1 -1 -1 1 -1 1 1 -1 -1 1 ...
## $ color: num [1:16] -0.63 2.51 -2.68 1.66 2.06 1.22 -2.09 1.93 6.79 5.47 ...head(datos,n=16L)## # A tibble: 16 x 6
## A B C D E color
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 -1 -1 -1 -1 1 -0.63
## 2 1 -1 -1 -1 -1 2.51
## 3 -1 1 -1 -1 -1 -2.68
## 4 1 1 -1 -1 1 1.66
## 5 -1 -1 1 -1 -1 2.06
## 6 1 -1 1 -1 1 1.22
## 7 -1 1 1 -1 1 -2.09
## 8 1 1 1 -1 -1 1.93
## 9 -1 -1 -1 1 -1 6.79
## 10 1 -1 -1 1 1 5.47
## 11 -1 1 -1 1 1 3.45
## 12 1 1 -1 1 -1 5.68
## 13 -1 -1 1 1 1 5.22
## 14 1 -1 1 1 -1 4.38
## 15 -1 1 1 1 -1 4.3
## 16 1 1 1 1 1 4.05se requerirá obtener la gráfica de Pareto y la de Daniel, estas gráficas darán evidencia de cuales son los efectos relevantes del experimento, con base en el Principio de Escasez. Para obtenerlas debemos correr el siguiente código:
library(FrF2)
experimento=FrF2(nruns = 16, nfactors = 5, factor.names = list(A=c(-1,1),B=c(-1,1),C=c(-1,1),D=c(-1,1),E=c(-1,1)),generators = "ABCD",replications = 1,randomize = FALSE)
experimento_respuesta=add.response(design = experimento,response = color)
halfnormal(experimento_respuesta, xlab = "Efectos activos")##
## Significant effects (alpha=0.05, Lenth method):## [1] D1
En esta gráfica se puede observar que los efectos activos en el diseño experimental es unicamente el efecto del factor D. Para confirmar la información que arroja la tabla, se verifica el gráfico de Daniel, mismo que se despliega utilizando el siguiente código:
DanielPlot(experimento_respuesta, main= "Gráfico de Daniel para el color de un producto químico")
En esta grafica puede confirmarse la información que se obtuvo del análisis anterior, por lo que se concluye que el unico efecto activo es el D.
b)
Procederemos a realizar la tabla ANOVA del experimento, pero considerando unidamente las interacciones dobles combinadas con el factor D, esto se realiza mediante el siguiente código:
#A:D#
modelo_a_d=lm(color~(A*D),data = datos)
summary(modelo_a_d)##
## Call:
## lm.default(formula = color ~ (A * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8450 -0.6913 -0.0350 0.6012 2.8950
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.7075 0.3410 7.939 4.07e-06 ***
## A 0.6550 0.3410 1.921 0.0789 .
## D 2.2100 0.3410 6.480 3.02e-05 ***
## A:D -0.6775 0.3410 -1.987 0.0703 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.364 on 12 degrees of freedom
## Multiple R-squared: 0.8053, Adjusted R-squared: 0.7566
## F-statistic: 16.54 on 3 and 12 DF, p-value: 0.0001459anova_a_d=aov(modelo_a_d)
summary(anova_a_d)## Df Sum Sq Mean Sq F value Pr(>F)
## A 1 6.86 6.86 3.689 0.0789 .
## D 1 78.15 78.15 41.991 3.02e-05 ***
## A:D 1 7.34 7.34 3.946 0.0703 .
## Residuals 12 22.33 1.86
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#B:D#
modelo_b_d=lm(color~(B*D),data = datos)
summary(modelo_b_d)##
## Call:
## lm.default(formula = color ~ (B * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.3850 -0.9613 -0.0700 1.2425 2.2250
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.7075 0.3894 6.952 1.53e-05 ***
## B -0.6700 0.3894 -1.720 0.111008
## D 2.2100 0.3894 5.675 0.000103 ***
## B:D 0.1225 0.3894 0.315 0.758497
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.558 on 12 degrees of freedom
## Multiple R-squared: 0.7461, Adjusted R-squared: 0.6826
## F-statistic: 11.75 on 3 and 12 DF, p-value: 0.0006947anova_b_d=aov(modelo_b_d)
summary(anova_b_d)## Df Sum Sq Mean Sq F value Pr(>F)
## B 1 7.18 7.18 2.960 0.111008
## D 1 78.15 78.15 32.205 0.000103 ***
## B:D 1 0.24 0.24 0.099 0.758497
## Residuals 12 29.12 2.43
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#C:D#
modelo_c_d=lm(color~(C*D),data = datos)
summary(modelo_c_d)##
## Call:
## lm.default(formula = color ~ (C * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.8950 -0.5394 0.2275 1.1825 2.2950
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.70750 0.42342 6.394 3.43e-05 ***
## C -0.07375 0.42342 -0.174 0.864631
## D 2.21000 0.42342 5.219 0.000215 ***
## C:D -0.35625 0.42342 -0.841 0.416599
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.694 on 12 degrees of freedom
## Multiple R-squared: 0.6999, Adjusted R-squared: 0.6248
## F-statistic: 9.327 on 3 and 12 DF, p-value: 0.001847anova_c_d=aov(modelo_c_d)
summary(anova_c_d)## Df Sum Sq Mean Sq F value Pr(>F)
## C 1 0.09 0.09 0.030 0.864631
## D 1 78.15 78.15 27.242 0.000215 ***
## C:D 1 2.03 2.03 0.708 0.416599
## Residuals 12 34.42 2.87
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#D:E#
modelo_d_e=lm(color~(D*E),data = datos)
summary(modelo_d_e)##
## Call:
## lm.default(formula = color ~ (D * E), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.6350 -0.9275 0.5325 1.1238 1.6200
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.70750 0.41939 6.456 3.13e-05 ***
## D 2.21000 0.41939 5.270 0.000198 ***
## E -0.41375 0.41939 -0.987 0.343341
## D:E 0.04375 0.41939 0.104 0.918641
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.678 on 12 degrees of freedom
## Multiple R-squared: 0.7055, Adjusted R-squared: 0.6319
## F-statistic: 9.584 on 3 and 12 DF, p-value: 0.001653anova_d_e=aov(modelo_d_e)
summary(anova_d_e)## Df Sum Sq Mean Sq F value Pr(>F)
## D 1 78.15 78.15 27.768 0.000198 ***
## E 1 2.74 2.74 0.973 0.343341
## D:E 1 0.03 0.03 0.011 0.918641
## Residuals 12 33.77 2.81
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1En nuestro análisis anterior podemos concluir que el factor que intervienen de manera significativa para maximizar el color de un producto Químico, considerando el valor p de la tabla anova de cada corrida, es: el facor D, el cual es la pureza del reactante, este factor principalmente influye significativamente en la coloración del producto químico.
c)
Para nuestro análisis gráfico, utlizaremos el siguiente coddígo:
efectos_principales=MEPlot(experimento_respuesta)
head(efectos_principales)## A B C D E
## - 2.0525 3.3775 2.78125 0.4975 3.12125
## + 3.3625 2.0375 2.63375 4.9175 2.29375grafica_interacciones=IAPlot(experimento_respuesta)
head(grafica_interacciones)## A:B A:C A:D A:E B:C B:D B:E C:D C:E D:E
## -:- 3.360 1.7325 -0.835 2.6175 3.5350 1.290 3.9350 0.2150 3.0750 0.9550
## +:- 3.395 3.8300 1.830 3.6250 2.0275 -0.295 2.3075 0.7800 3.1675 5.2875
## -:+ 0.745 2.3725 4.940 1.4875 3.2200 5.465 2.8200 5.3475 2.4875 0.0400
## +:+ 3.330 2.8950 4.895 3.1000 2.0475 4.370 1.7675 4.4875 2.1000 4.5475- Como era previsto, tanto en la gráfica de efectos individuales como la de efectos interactivos podemos observar que el factor D provoca un efecto significativo en el color de un producto Químico. Cabe destacar que el factor C y E son los que menos interacción y efecto significativo tienen.
D)
No cabe duda, que el factor D el cual es la pureza del reactante es el mejor tratimiento y la respuesta predicha por el modelo, que tiene un gran impacto en el color del producto Químico, esto debido a que en todas las pruebas estadisticas que se utilizaron en este modelo, fue el unico factor que dio respuestas diferentes a las demas, y se pudo observar en la grafíca de efectos individuales como influye este factor en el proceso.
E)
Para hacer nuestro analisis residual, debemos hacer 2 pruebas, el de normalidad y de homocedasticidad, para ello utlizaremos la siguiente codificación:
normalidad=shapiro.test(resid(modelo_a_d))
print(normalidad)##
## Shapiro-Wilk normality test
##
## data: resid(modelo_a_d)
## W = 0.95541, p-value = 0.5799normalidad2=shapiro.test(resid(modelo_b_d))
print(normalidad2)##
## Shapiro-Wilk normality test
##
## data: resid(modelo_b_d)
## W = 0.96133, p-value = 0.6861normalidad3=shapiro.test(resid(modelo_c_d))
print(normalidad3)##
## Shapiro-Wilk normality test
##
## data: resid(modelo_c_d)
## W = 0.92953, p-value = 0.2396normalidad4=shapiro.test(resid(modelo_d_e))
print(normalidad4)##
## Shapiro-Wilk normality test
##
## data: resid(modelo_d_e)
## W = 0.89309, p-value = 0.06235Homosedasticidad=bartlett.test(resid(modelo_a_d),A,D,data= datos)
print(Homosedasticidad)##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo_a_d) and A
## Bartlett's K-squared = 5.385, df = 1, p-value = 0.02031Homosedasticidad2=bartlett.test(resid(modelo_b_d),B,D,data= datos)
print(Homosedasticidad2)##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo_b_d) and B
## Bartlett's K-squared = 1.1263, df = 1, p-value = 0.2886Homosedasticidad3=bartlett.test(resid(modelo_c_d),C,D,data= datos)
print(Homosedasticidad3)##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo_c_d) and C
## Bartlett's K-squared = 0.58703, df = 1, p-value = 0.4436Homosedasticidad4=bartlett.test(resid(modelo_d_e),D,E,data= datos)
print(Homosedasticidad4)##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo_d_e) and D
## Bartlett's K-squared = 2.7658, df = 1, p-value = 0.0963Prueba de normalidad Shapiro-Wilks
- Observamos que nuestro valor p de nuestras 4 pruebas es mayor al nivel de significancia 0.05, por lo tanto no rechazamos la hipótesis nula en ningún caso, esto quiere decir que los residuos provienen de una población con distribución normal.
Prueba de Homocedasticidad
- Observamos que nuestro valores p solo varía en una combinación el cual es el AD,los demas modelos sus valores p es mucho mayor al nivel de significancia 0.05, por lo tanto, no se rechaza la hipótesis nula, esto quiere decir, que las muestras presentan varianzas iguales en estas 3 interacciones.