Diseño factorial fraccionado
Ing. Anahy Hernández Trinidad
6/10/2021
1 Diseño factorial fraccionado
1.1 Ejemplo
8-6. R-D Snee (“Experimentación con un número grande de variables”, en Experiments in Industry: Design, Analysis and Interpretation of Results, de R.D Snee, L.B. Hare y J.B. Trout, editores, ASQC) describe un experimento en el que se usó un diseño \(2^5-^1\) con I=ABCDE para investigar los efectos de cinco factores sobre el color de un producto químico. Los factores son A=solvente/reactivo, B=catalizador/reactivo, C=temperatura, D=pureza del reactivo y E=pH del reactivo. Los resultados obtenidos fueron los siguientes:(Montgomery, 2004)
| e=-0.63 | d=6.79 |
| a=2.51 | ade=5.47 |
| b=-2.68 | bde=3.45 |
| abe=1.66 | adb=5.68 |
| c=2.06 | cde=5.22 |
| ace=1.22 | acd=4.38 |
| bce=-2.09 | bcd=4.30 |
| abc=1.93 | abcde=4.05 |
a) Construir una gráfica de probabilidad normal de los efectos ¿Qué efectos parecen estar activos?
b) Calcular los residuales. Construir una gráfica de probabilidad normal de los residuales y graficar los residuales contra los valores ajustados. Comentar las gráficas.
c) Si algunos de los factores son insignificantes, plegar el diseño \(2^5-^1\) a un diseño factorial completo en los factores activos. Comentar el diseño resultante e interpretar los resultados.
1.2 Solución del ejercicio
a) Construir una gráfica de probabilidad normal de los efectos ¿Qué efectos parecen estar activos?
library(printr)
library(FrF2)
datos=read.table("dataset.txt", header = TRUE)
str(datos)## 'data.frame': 16 obs. of 6 variables:
## $ solvente_reactivo : int -1 1 -1 1 -1 1 -1 1 -1 1 ...
## $ catalizador_reactivo: int -1 -1 1 1 -1 -1 1 1 -1 -1 ...
## $ temperatura : int -1 -1 -1 -1 1 1 1 1 -1 -1 ...
## $ pureza_reactivo : int -1 -1 -1 -1 -1 -1 -1 -1 -1 1 ...
## $ ph_reactivo : int 1 -1 -1 1 -1 1 1 -1 1 1 ...
## $ Resultados : num -0.63 2.51 -2.68 1.66 2.06 1.22 -2.09 1.93 6.79 5.47 ...View(datos)
attach(datos)A=factor(`solvente_reactivo`)
B=factor(`catalizador_reactivo`)
C=factor(`temperatura`)
D=factor(`pureza_reactivo`)
E=factor(`ph_reactivo`)
head(datos, n=16L)| solvente_reactivo | catalizador_reactivo | temperatura | pureza_reactivo | ph_reactivo | Resultados |
|---|---|---|---|---|---|
| -1 | -1 | -1 | -1 | 1 | -0.63 |
| 1 | -1 | -1 | -1 | -1 | 2.51 |
| -1 | 1 | -1 | -1 | -1 | -2.68 |
| 1 | 1 | -1 | -1 | 1 | 1.66 |
| -1 | -1 | 1 | -1 | -1 | 2.06 |
| 1 | -1 | 1 | -1 | 1 | 1.22 |
| -1 | 1 | 1 | -1 | 1 | -2.09 |
| 1 | 1 | 1 | -1 | -1 | 1.93 |
| -1 | -1 | -1 | -1 | 1 | 6.79 |
| 1 | -1 | -1 | 1 | 1 | 5.47 |
| -1 | 1 | -1 | 1 | 1 | 3.45 |
| 1 | 1 | -1 | 1 | -1 | 5.68 |
| -1 | -1 | 1 | 1 | 1 | 5.22 |
| 1 | -1 | 1 | 1 | -1 | 4.38 |
| -1 | 1 | 1 | 1 | -1 | 4.30 |
| 1 | 1 | 1 | 1 | 1 | 4.05 |
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 = Resultados)
halfnormal(experimento_respuesta, xlab = "Efectos activos")
En la gráfica se puede observar que los efectos activos en lo que respecta al diseño de experimentos es unicamente el factor D. Y para poder corroborar los datos que arroja la tabla se aplica en el gráfico de Daniel, el código es el siguiente:
DanielPlot(experimento_respuesta, main= "Gráfico de Daniel para Resultados")
En la gráfica de Daniel puede confirmarse que los datos que se obtuvieron en al análisis anterior corresponden al mismo, por lo que se concluye que los efectos activos, o en este caso, el efecto activo es el factor D.
En el siguiente inciso se pueden revisar los efectos principales para poder verificar cuales son los efectos de cinco factores sobre el color de un producto químico utilizando el siguiente comando:
b) Calcular los residuales. Construir una gráfica de probabilidad normal de los residuales y graficar los residuales contra los valores ajustados. Comentar las gráficas.
Prueba de Normalidad de Residuos de Shapiro-Wilk
\[H_{0}:R \, {\in} N\, ({\mu=0}, {\sigma^2}=Constante)\]
\[H_{1}:R \, {\notin}\, N\, ({\mu=0}, {\sigma^2}<Constante)\]
modelo= aov(Resultados~(A*B*C*D*E))
summary(modelo)## Df Sum Sq Mean Sq F value Pr(>F)
## A 1 6.86 6.86 0.249 0.705
## B 1 7.18 7.18 0.261 0.699
## C 1 0.09 0.09 0.003 0.964
## D 1 49.95 49.95 1.815 0.407
## E 1 0.01 0.01 0.001 0.986
## A:B 1 12.55 12.55 0.456 0.622
## A:C 1 0.28 0.28 0.010 0.936
## B:C 1 1.96 1.96 0.071 0.834
## A:D 1 3.91 3.91 0.142 0.771
## B:D 1 1.36 1.36 0.050 0.861
## C:D 1 0.51 0.51 0.018 0.914
## A:E 1 0.00 0.00 0.000 0.996
## B:E 1 0.01 0.01 0.000 0.990
## C:E 1 2.49 2.49 0.090 0.814
## Residuals 1 27.53 27.53Con los resultados arrojados por la tabla ANOVA, y contando con el valor de nivel de significancia de \({\propto}=0.05\), en la que se evaluan 5 factores a los cuales les afecta el color del producto químico, se puede observa que el \(Valor_p\) para todos los efectos produccidos de los diferentes factores no son significativos.
normalidad=shapiro.test(resid(modelo))
print(normalidad)##
## Shapiro-Wilk normality test
##
## data: resid(modelo)
## W = 0.51123, p-value = 2.63e-06Mediante la prueba de Shapiro-Wilk realizada para observar los residuos, utilizado el \(Valor_p < 0.05\) dado que el resultado es menor al valor antes mencionado, se dice que se rechaza la hipótesis nula y se puede deir que los residuos no tienen normalidad.
qqnorm(resid(modelo), main= "Gráfica de Probabilidad para los Residuales del Modelo", xlab="Cuantiles Teoricos", ylab = "Cuantiles de muestra")
qqline(resid(modelo))
La prueba de Shapiro-Wilk confirma que mediante la gráfica de probabilidad para los residuales del modelo muestra con cierta evidencia que hay una falta de normalidad en lo que respecta a los residuos, por lo que se puede decir que existen datos atípicos en el modelo propuesto, los cuales provocan que no sea idóneo para explicar la interacción que exise entre del color del quimico respescto a los factores que se evaluan.
Prueba Igualdad de Varianzas de Bartlett
\[H_{0}:\sigma^{2}_{i}=\sigma^{2}_{j}=Constante\] \[H_{1}:\sigma^{2}_{i} \neq \sigma^{2}_{j}\neq Constante\]
homocedasticidad=bartlett.test(resid(modelo),A,B,C,D,E,data=experimento_resp)
print(homocedasticidad)##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo) and A
## Bartlett's K-squared = 457.1, df = 1, p-value < 2.2e-16En base a los resultados obtenidos por el análisis de homocedasticidad, estos son homocedásticos, por lo que se puede decir que el modelo de regresión no es el adecuado.
c) Si algunos de los factores son insignificantes, plegar el diseño \(2^5-^1\) a un diseño factorial completo en los factores activos. Comentar el diseño resultante e interpretar los resultados.
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.29375 |
De la gráfica anteriormente hecha, se puede deducir que los efectos de cinco factores que generan un increment en el color de un producto químico es el factor D y el A, por otro lado hay una disminucion en los factores B, C y E. En el caso del factor A, se da un aumento del 1.31, en el caso del factor B se da un cambio de nivel máximo a mínimo, el cual provoca una disminución de 1.34, para el caso del C, este tiene una disminución de 0.1475, en el D se da un cambio del nivel de mínimo a máximo, el cual provoca un aumento del 4.42, mientras que en el E hay una disminución del 0.8275.
grafica_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 es de observarse en la gráfica de las interraciones existen algunas que se consideran fuertes, entre ellas estan los factores de A:D, B:D, C:D y D:E. Y para coprobar si estas son significativas se realiza un ANOVA del experiment, pero solo tomando unicamente las interacciones dobles.
modelo_a_d=lm(Resultados~(A*D),data = datos)
summary(modelo_a_d)##
## Call:
## lm.default(formula = Resultados ~ (A * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.3700 -0.8521 -0.0967 0.7062 6.1000
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.6900 1.0351 0.667 0.5176
## A1 1.1400 1.5526 0.734 0.4769
## D1 3.6333 1.6903 2.150 0.0527 .
## A1:D1 -0.5683 2.3527 -0.242 0.8132
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.314 on 12 degrees of freedom
## Multiple R-squared: 0.4395, Adjusted R-squared: 0.2994
## F-statistic: 3.137 on 3 and 12 DF, p-value: 0.06541anova_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 1.281 0.2797
## D 1 43.23 43.23 8.070 0.0149 *
## A:D 1 0.31 0.31 0.058 0.8132
## Residuals 12 64.28 5.36
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modelo_b_d=lm(Resultados~(B*D),data = datos)
summary(modelo_b_d)##
## Call:
## lm.default(formula = Resultados ~ (B * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.0200 -0.9825 -0.1950 0.6625 4.4000
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.3900 0.9217 2.593 0.0235 *
## B1 -2.6850 1.3826 -1.942 0.0760 .
## D1 2.6333 1.5052 1.749 0.1057
## B1:D1 2.0317 2.0952 0.970 0.3513
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.061 on 12 degrees of freedom
## Multiple R-squared: 0.5555, Adjusted R-squared: 0.4444
## F-statistic: 4.999 on 3 and 12 DF, p-value: 0.01778anova_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 1.691 0.21792
## D 1 52.53 52.53 12.366 0.00425 **
## B:D 1 3.99 3.99 0.940 0.35134
## Residuals 12 50.98 4.25
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modelo_c_d=lm(Resultados~(C*D),data = datos)
summary(modelo_c_d)##
## Call:
## lm.default(formula = Resultados ~ (C * D), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2100 -0.6823 0.2850 0.8550 5.2600
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.5300 1.0507 1.456 0.1710
## C1 -0.7500 1.5760 -0.476 0.6427
## D1 3.3367 1.7157 1.945 0.0756 .
## C1:D1 0.3708 2.3882 0.155 0.8792
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.349 on 12 degrees of freedom
## Multiple R-squared: 0.4225, Adjusted R-squared: 0.2781
## F-statistic: 2.926 on 3 and 12 DF, p-value: 0.07714anova_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.016 0.902
## D 1 48.23 48.23 8.739 0.012 *
## C:D 1 0.13 0.13 0.024 0.879
## Residuals 12 66.23 5.52
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modelo_d_e=lm(Resultados~(D*E),data = datos)
summary(modelo_d_e)##
## Call:
## lm.default(formula = Resultados ~ (D * E), data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.6350 -0.6475 0.0500 0.9356 5.4000
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.9550 1.1833 0.807 0.4353
## D1 3.8317 1.8075 2.120 0.0556 .
## E1 0.4350 1.5876 0.274 0.7887
## D1:E1 -0.6742 2.4057 -0.280 0.7841
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.367 on 12 degrees of freedom
## Multiple R-squared: 0.414, Adjusted R-squared: 0.2674
## F-statistic: 2.825 on 3 and 12 DF, p-value: 0.0836anova_d_e=aov(modelo_d_e)
summary(anova_d_e)## Df Sum Sq Mean Sq F value Pr(>F)
## D 1 46.96 46.96 8.384 0.0134 *
## E 1 0.08 0.08 0.014 0.9076
## D:E 1 0.44 0.44 0.079 0.7841
## Residuals 12 67.21 5.60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Con el análisis anterior se puede concluir que de los 5 factores que intervienen de manera significativa en el color del producto de un químico, considerando el \(Valor_p\) de la tabla anova de cada corrida, es el factor D correspondiente a la pureza del reactivo, por lo que se recomieda mantenerlo en un nivel alto.