En una empresa lechera se han tenido problemas con la viscosidad de cierta bebida de chocolate. Se cree que con tres ingredientes que se agregan en pequeñas cantidades se puede resolver este problema, por lo que es necesario explorar la situación; para ello se corre un experimento \(2^{3}\) con 2 réplicas. A continuación se aprecián los resultados obtenidos:
resultados de la tabla
Lo primero que se debe realizar es mandar a instruir las librerías que se van a utlizar 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)
datos=read_excel(path = "dataset.xlsx")
View(datos)
datos=read.table("dataset.txt",header = TRUE)
str(datos)
## 'data.frame': 16 obs. of 4 variables:
## $ Viscosidad : num 13.3 14.7 14.6 14.3 16.9 15.5 17.4 18.9 13.9 14.4 ...
## $ ingredienteA: int -1 1 -1 1 -1 1 -1 1 -1 1 ...
## $ ingredienteB: int -1 -1 1 1 -1 -1 1 1 -1 -1 ...
## $ ingredienteC: int -1 -1 -1 -1 1 1 1 1 -1 -1 ...
attach(datos)
Para efectos de resolver adecuadamente el ejercicio, reescribimos la tabla, de manera que sea interpretada como un conjunto de variables, quedando de la siguiente manera:
| Ingrediente A | Ingradiente B | Ingrediente C | Viscosidad |
|---|---|---|---|
| -1 | -1 | -1 | 13.3,13.9 |
| +1 | -1 | -1 | 14.7,14.4 |
| -1 | +1 | -1 | 14.6,14.9 |
| +1 | +1 | -1 | 14.3,14.1 |
| -1 | -1 | +1 | 16.9,17.2 |
| +1 | -1 | +1 | 15.5,15.1 |
| -1 | +1 | +1 | 17.4,17.1 |
| +1 | +1 | +1 | 18.9,19.2 |
H0:A=0 El ingrediente A no afecta a la viscosidad de la bebida del chocolate .
H1:A≠0 El ingrediente A afecta a la viscosidad de la bebida del chocolate.
H0:B=0 El ingrediente B no afecta a la viscosidad de la bebida del chocolate .
H1:B≠0 El ingrediente B afecta a la viscosidad de la bebida del chocolate.
H0:C=0 El ingrediente C no afecta a la viscosidad de la bebida del chocolate .
H1:C≠0 El ingrediente C afecta a la viscosidad de la bebida del chocolate.
H0:AB=0 El ingrediente A y B no afecta a la viscosidad de la bebida del chocolate .
H1:AB≠0 El ingrediente A y B afecta a la viscosidad de la bebida del chocolate.
H0:AC=0 El ingrediente A y C no afecta a la viscosidad de la bebida del chocolate .
H1:AC≠0 El ingrediente A y C afecta a la viscosidad de la bebida del chocolate.
H0:BC=0 El ingrediente B y C no afecta a la viscosidad de la bebida del chocolate .
H1:BC≠0 El ingrediente B y C afecta a la viscosidad de la bebida del chocolate.
H0:ABC=0 El ingrediente A, B y C no afecta a la viscosidad de la bebida del chocolate .
H1:ABC≠0 El ingrediente A, B y C afecta a la viscosidad de la bebida del chocolate.*
El modelo estadístico para este diseño es:
Para determinar los efectos, utilizamos el siguiente código:
f_a=factor(`ingredienteA`)
f_b=factor(`ingredienteB`)
f_c=factor(`ingredienteC`)
f_viscosidad=factor(Viscosidad)
modelo=lm(Viscosidad~(f_a+f_b+f_c+f_a*f_b+f_a*f_c+f_b*f_c+f_a*f_b*f_c))
summary(modelo)
##
## Call:
## lm(formula = Viscosidad ~ (f_a + f_b + f_c + f_a * f_b + f_a *
## f_c + f_b * f_c + f_a * f_b * f_c))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.30 -0.15 0.00 0.15 0.30
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13.6000 0.1777 76.551 9.45e-13 ***
## f_a1 0.9500 0.2512 3.781 0.00538 **
## f_b1 1.1500 0.2512 4.577 0.00181 **
## f_c1 3.4500 0.2512 13.732 7.63e-07 ***
## f_a1:f_b1 -1.5000 0.3553 -4.222 0.00291 **
## f_a1:f_c1 -2.7000 0.3553 -7.599 6.31e-05 ***
## f_b1:f_c1 -0.9500 0.3553 -2.674 0.02820 *
## f_a1:f_b1:f_c1 5.0500 0.5025 10.050 8.18e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2512 on 8 degrees of freedom
## Multiple R-squared: 0.9898, Adjusted R-squared: 0.9809
## F-statistic: 110.8 on 7 and 8 DF, p-value: 2.493e-07
En esta parte podemos observar el r cuadrado ajustado donde podemos obtener la variabilidad de las respuestas es de 98.09% por los factores involucrados en el análisis, en este caso el ingrediente A, B y C, y el 1.99% esta explicados por otros factores no considerados.
Utilizaremos los siguientes comandos para obtener nuestra tabla Anova y hacer nuestra prueba de hipótesis y comprobar si hay un efecto significativo tanto en los factores como en la interacción
anova=aov(modelo)
summary(anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## f_a 1 0.05 0.05 0.802 0.396647
## f_b 1 5.64 5.64 89.356 1.29e-05 ***
## f_c 1 33.35 33.35 528.327 1.36e-08 ***
## f_a:f_b 1 1.05 1.05 16.644 0.003536 **
## f_a:f_c 1 0.03 0.03 0.485 0.505830
## f_b:f_c 1 2.48 2.48 39.297 0.000241 ***
## f_a:f_b:f_c 1 6.38 6.38 101.000 8.18e-06 ***
## Residuals 8 0.50 0.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Aplicando un ANOVA mejorado:
f_a2=factor(`ingredienteA`)
f_b2=factor(`ingredienteB`)
f_c2=factor(`ingredienteC`)
modelo2=lm(Viscosidad~(f_b2+f_c2+f_a2*f_b2+f_b2*f_c2+f_a2*f_b2*f_c2))
anova2=aov(modelo2)
summary(anova2)
## Df Sum Sq Mean Sq F value Pr(>F)
## f_b2 1 5.64 5.64 89.356 1.29e-05 ***
## f_c2 1 33.35 33.35 528.327 1.36e-08 ***
## f_a2 1 0.05 0.05 0.802 0.396647
## f_b2:f_a2 1 1.05 1.05 16.644 0.003536 **
## f_b2:f_c2 1 2.48 2.48 39.297 0.000241 ***
## f_c2:f_a2 1 0.03 0.03 0.485 0.505830
## f_b2:f_c2:f_a2 1 6.38 6.38 101.000 8.18e-06 ***
## Residuals 8 0.50 0.06
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
CONCLUSIÓN: Considerando un ANOVA MEJORADO vemos que los efectos antes mencionados siguen siendo significativos, por lo tanto no hay un tratamiento para minimizar dichos efectos.
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':
##
## lengths
experimento=FrF2(nruns=8,nfactors=3, factor.names= list(f_a=c(-1, 1),f_b=c(-1, 1),f_c=c(-1, 1)),replications=2,randomize=FALSE)
## creating full factorial with 8 runs ...
experimento_resp=add.response(desig=experimento, response=Viscosidad)
MEPlot(experimento_resp)
IAPlot(experimento_resp)
#———-Analisis residuales———-#
normalidad=shapiro.test(resid(modelo))
print(normalidad)
##
## Shapiro-Wilk normality test
##
## data: resid(modelo)
## W = 0.87343, p-value = 0.0307
Viscosidad=bartlett.test(resid(modelo),f_a,f_b,f_c,data=experimento_resp)
print(Viscosidad)
##
## Bartlett test of homogeneity of variances
##
## data: resid(modelo) and f_a
## Bartlett's K-squared = 0.41308, df = 1, p-value = 0.5204