library(pid)
## Loading required package: ggplot2
## Loading required package: png
## Loading required package: FrF2
## Loading required package: DoE.base
## Loading required package: grid
## Loading required package: 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
library(readxl)
library(ggplot2)

library(multcomp)
## Loading required package: mvtnorm
## Loading required package: survival
## Loading required package: TH.data
## Loading required package: MASS
## 
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
## 
##     geyser
library(multcompView)
library(sandwich)

library(FrF2)
# Función auxiliar para hallar el error estándar de los resdiuales
# de un modelo "aov"
options(width=120)
errsdred<-function(modelo){
  aovm<-anova(modelo)
  return(sqrt(aovm$"Mean Sq"[length(aovm$"Mean Sq")]))
}

Diseño fraccionado (Primera etapa)

Experimento de la tasa de filtración.

Un producto químico se produce en una cámara de presión. Se llevó a cabo un experimento factorial fraccionado para estudiar los factores que influencian en la tasa de filtración del producto obtenido. Los cuatro factores son: temperatura (A), presión (B), concentración de formaldehido (C) y velocidad de agitación (D). Cada uno de los factores se tomó a dos niveles.

Se realizó un experimento \(2^{4-1}_{IV}\).

Efectos dobles y sus alias

library(pid)
tradeoff(runs = 8, factors = 4)
## With 8 experiments, and 4 factors:
##   Resolution: IV
##   Generator: D=ABC
##   Aliasing (related ONLY to main effects and 2-factor interactions):
##       Main effects are not aliased with 2-factor interactions.
##       AB=CD
##       AC=BD
##       AD=BC

Trabla de diseños fraccionados

library(pid)
tradeOffTable()

La resolución se refiere al grado de alias que tienen las interacciones con los factores principales y entre si.

Así por ejemplo:

Resolución III. Ningún efecto principal se alia con otro efecto principal, pero los efectos principales se alian con interacciones dobles y estas dobles entre si.

Resolución IV. Los efectos principales no se alian con otros principales con con interacciones dobles. Pero las interacciones dobles pueden aliarse entre ellas.

Resolución V. Los efectos principales no se alian con efectos principales no con interacciones dobles. Y ninguna interacción doble se alía entre ellas.

Lectura de la base de datos.

library(readxl)
tasaFilt <- read_excel("tasaFiltracion.xlsx")
tasaFilt
## # A tibble: 8 x 5
##       A     B     C     D    TF
##   <dbl> <dbl> <dbl> <dbl> <dbl>
## 1   -1.   -1.   -1.   -1.   45.
## 2    1.   -1.   -1.    1.  100.
## 3   -1.    1.   -1.    1.   45.
## 4    1.    1.   -1.   -1.   65.
## 5   -1.   -1.    1.    1.   75.
## 6    1.   -1.    1.   -1.   60.
## 7   -1.    1.    1.   -1.   80.
## 8    1.    1.    1.    1.   96.
tasaFilt$A <- factor(tasaFilt$A, labels = c("b-T", "a-T"))
tasaFilt$B <- factor(tasaFilt$B, labels = c("b-P", "a-P"))
tasaFilt$C <- factor(tasaFilt$C, labels = c("b-H2C=O", "a-H2C=O"))
tasaFilt$D <- factor(tasaFilt$D, labels = c("b-V", "a-V"))
tasaFilt
## # A tibble: 8 x 5
##   A     B     C       D        TF
##   <fct> <fct> <fct>   <fct> <dbl>
## 1 b-T   b-P   b-H2C=O b-V     45.
## 2 a-T   b-P   b-H2C=O a-V    100.
## 3 b-T   a-P   b-H2C=O a-V     45.
## 4 a-T   a-P   b-H2C=O b-V     65.
## 5 b-T   b-P   a-H2C=O a-V     75.
## 6 a-T   b-P   a-H2C=O b-V     60.
## 7 b-T   a-P   a-H2C=O b-V     80.
## 8 a-T   a-P   a-H2C=O a-V     96.

Gráfica exploratoria

library(ggplot2)
g1 <- ggplot(tasaFilt, aes(A, TF, col = B))
g1 + geom_point(size = 5) + facet_grid(C ~ D)

Gráfica por factor principal.

ggplot(tasaFilt, aes(A, TF)) +
  geom_point() +
  stat_summary(fun.y = mean, geom = "line", lwd = 2, aes(group = 1))

ggplot(tasaFilt, aes(B, TF)) +
  geom_point() +
  stat_summary(fun.y = mean, geom = "line", lwd = 2, aes(group = 1))

ggplot(tasaFilt, aes(C, TF)) +
  geom_point() +
  stat_summary(fun.y = mean, geom = "line", lwd = 2, aes(group = 1))

ggplot(tasaFilt, aes(D, TF)) +
  geom_point() +
  stat_summary(fun.y = mean, geom = "line", lwd = 2, aes(group = 1))

Modelización

\[ y_{ijklm} = \mu + \alpha_i + \beta_j + \gamma_k + \delta_l + \\ (\alpha\beta)_{ij} + (\alpha\gamma)_{ik} + (\alpha\delta)_{il} + \\ \epsilon_{ijklm} \]

Donde: \[ \epsilon_{ijklm} \sim \mathcal{N}(0, \sigma^2) \quad i.i.d. \]

Además, se hacen los siguientes supuestos con respecto a las interacciones dobles:

\[ (\alpha\beta)_{ij} = (\gamma\delta)_{kl} \\ (\alpha\gamma)_{ik} = (\beta\delta)_{jl} \\ (\alpha\delta)_{il} = (\beta\gamma)_{jk} \] Todas las demás interacciones de orden superior se suponen no significactivas.

modelo1 <- aov(TF ~ A + B + C + D + A:B + A:C + A:D, data = tasaFilt)
summary(modelo1)
##             Df Sum Sq Mean Sq
## A            1  722.0   722.0
## B            1    4.5     4.5
## C            1  392.0   392.0
## D            1  544.5   544.5
## A:B          1    2.0     2.0
## A:C          1  684.5   684.5
## A:D          1  722.0   722.0

Determinación de interacciones significativas

library(pid)
paretoPlot(modelo1)

efectos <- unlist(model.tables(modelo1))
efectos <- efectos[substring(names(efectos),1,6) == "tables"]
names(efectos)<-gsub("tables.", "", names(efectos))

qqR <- qqnorm(efectos, xlim = c(-3, 3), ylim = c(-15,15), pch = 19, cex = 0.5)
qqline(efectos)
with(qqR,text(x, y, labels = names(efectos), cex = 0.7, pos = 3))

modelo2 <- aov(TF ~ A + B + C + D  + A:C + A:D, data = tasaFilt)
summary(modelo2)
##             Df Sum Sq Mean Sq F value Pr(>F)  
## A            1  722.0   722.0  361.00 0.0335 *
## B            1    4.5     4.5    2.25 0.3743  
## C            1  392.0   392.0  196.00 0.0454 *
## D            1  544.5   544.5  272.25 0.0385 *
## A:C          1  684.5   684.5  342.25 0.0344 *
## A:D          1  722.0   722.0  361.00 0.0335 *
## Residuals    1    2.0     2.0                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Modelización definitiva

modelo3 <- aov(TF ~ A  + C + D  + A:C + A:D, data = tasaFilt)
summary(modelo3)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## A            1  722.0   722.0   222.2 0.00447 **
## C            1  392.0   392.0   120.6 0.00819 **
## D            1  544.5   544.5   167.5 0.00592 **
## A:C          1  684.5   684.5   210.6 0.00471 **
## A:D          1  722.0   722.0   222.2 0.00447 **
## Residuals    2    6.5     3.3                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Diagnósticos

# Normalidad
residualesEstand <- rstandard(modelo3)
shapiro.test(residualesEstand)
## 
##  Shapiro-Wilk normality test
## 
## data:  residualesEstand
## W = 0.90039, p-value = 0.2913
# Igualdad de varianza
bartlett.test(TF ~ A, data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by A
## Bartlett's K-squared = 0.021022, df = 1, p-value = 0.8847
bartlett.test(TF ~ B, data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by B
## Bartlett's K-squared = 0.014884, df = 1, p-value = 0.9029
bartlett.test(TF ~ C, data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by C
## Bartlett's K-squared = 0.76337, df = 1, p-value = 0.3823
bartlett.test(TF ~ D, data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by D
## Bartlett's K-squared = 0.75823, df = 1, p-value = 0.3839
bartlett.test(TF ~ paste(A,C), data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by paste(A, C)
## Bartlett's K-squared = Inf, df = 3, p-value < 2.2e-16
bartlett.test(TF ~ paste(A,D), data = tasaFilt)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by paste(A, D)
## Bartlett's K-squared = 3.722, df = 3, p-value = 0.2931

Comparaciones de Tukey

library(multcomp)
library(multcompView)
library(sandwich)
tasaFilt$AC <- with(tasaFilt, interaction(A, C))

modelo3_1 <- aov(TF ~ -1 + AC, data = tasaFilt)
summary(modelo3_1)
##           Df Sum Sq Mean Sq F value  Pr(>F)   
## AC         4  41843   10461   32.87 0.00256 **
## Residuals  4   1273     318                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
compMultHetero <- glht(modelo3_1, linfct = mcp(AC = "Tukey"), vcov = sandwich)
summary(compMultHetero)
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: aov.default(formula = TF ~ -1 + AC, data = tasaFilt)
## 
## Linear Hypotheses:
##                                Estimate Std. Error t value Pr(>|t|)    
## a-T.b-H2C=O - b-T.b-H2C=O == 0   37.500     12.374   3.030    0.106    
## b-T.a-H2C=O - b-T.b-H2C=O == 0   32.500      1.768  18.385   <0.001 ***
## a-T.a-H2C=O - b-T.b-H2C=O == 0   33.000     12.728   2.593    0.161    
## b-T.a-H2C=O - a-T.b-H2C=O == 0   -5.000     12.500  -0.400    0.971    
## a-T.a-H2C=O - a-T.b-H2C=O == 0   -4.500     17.752  -0.253    0.992    
## a-T.a-H2C=O - b-T.a-H2C=O == 0    0.500     12.850   0.039    1.000    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
op <- par(no.readonly = TRUE) 
par(mar=op$mar+c(0,12,0,0))
plot(compMultHetero)

par(op)

Comportamiento del modelo

tasaFilt$TFp <- predict(modelo3)
g1 <- ggplot(tasaFilt, aes(A, TF, col = C))
g1 + geom_point() + facet_grid(. ~ D) +
  geom_line(aes(A, TFp, group = C)) +
  geom_point(aes(A, TFp, col = C), shape = 17)

Valor de los efectos

model.tables(modelo3, se = TRUE)
## Tables of effects
## 
##  A 
## A
##  b-T  a-T 
## -9.5  9.5 
## 
##  C 
## C
## b-H2C=O a-H2C=O 
##      -7       7 
## 
##  D 
## D
##   b-V   a-V 
## -8.25  8.25 
## 
##  A:C 
##      C
## A     b-H2C=O a-H2C=O
##   b-T -9.25    9.25  
##   a-T  9.25   -9.25  
## 
##  A:D 
##      D
## A     b-V  a-V 
##   b-T  9.5 -9.5
##   a-T -9.5  9.5
## 
## Standard errors of effects
##              A      C      D    A:C    A:D
##         0.9014 0.9014 0.9014 1.2748 1.2748
## replic.      4      4      4      2      2
model.tables(modelo3, se = TRUE, type = "mean")
## Tables of means
## Grand mean
##       
## 70.75 
## 
##  A 
## A
##   b-T   a-T 
## 61.25 80.25 
## 
##  C 
## C
## b-H2C=O a-H2C=O 
##   63.75   77.75 
## 
##  D 
## D
##  b-V  a-V 
## 62.5 79.0 
## 
##  A:C 
##      C
## A     b-H2C=O a-H2C=O
##   b-T 45.0    77.5   
##   a-T 82.5    78.0   
## 
##  A:D 
##      D
## A     b-V  a-V 
##   b-T 62.5 60.0
##   a-T 62.5 98.0
## 
## Standard errors for differences of means
##             A     C     D   A:C   A:D
##         1.275 1.275 1.275 1.803 1.803
## replic.     4     4     4     2     2

Estimación de los errores

errsdred(modelo3)
## [1] 1.802776

Diseño factorial (\(2^k\))

Realización de los datos complementarios

library(readxl)
tasaFilt2 <- read_excel("tasaFiltracion2.xlsx")
tasaFilt2
## # A tibble: 16 x 6
##        A     B     C     D    TF Etapa
##    <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
##  1   -1.   -1.   -1.   -1.   45.    1.
##  2    1.   -1.   -1.    1.  100.    1.
##  3   -1.    1.   -1.    1.   45.    1.
##  4    1.    1.   -1.   -1.   65.    1.
##  5   -1.   -1.    1.    1.   75.    1.
##  6    1.   -1.    1.   -1.   60.    1.
##  7   -1.    1.    1.   -1.   80.    1.
##  8    1.    1.    1.    1.   96.    1.
##  9   -1.   -1.   -1.    1.   43.    2.
## 10    1.   -1.   -1.   -1.   71.    2.
## 11   -1.    1.   -1.   -1.   48.    2.
## 12    1.    1.   -1.    1.  104.    2.
## 13   -1.   -1.    1.   -1.   68.    2.
## 14    1.   -1.    1.    1.   86.    2.
## 15   -1.    1.    1.    1.   70.    2.
## 16    1.    1.    1.   -1.   65.    2.
tasaFilt2$A <- factor(tasaFilt2$A, labels = c("b-T", "a-T"))
tasaFilt2$B <- factor(tasaFilt2$B, labels = c("b-P", "a-P"))
tasaFilt2$C <- factor(tasaFilt2$C, labels = c("b-H2C=O", "a-H2C=O"))
tasaFilt2$D <- factor(tasaFilt2$D, labels = c("b-V", "a-V"))
tasaFilt2$Etapa <- factor(tasaFilt2$Etapa, labels = c("e-1", "e-2"))
tasaFilt2
## # A tibble: 16 x 6
##    A     B     C       D        TF Etapa
##    <fct> <fct> <fct>   <fct> <dbl> <fct>
##  1 b-T   b-P   b-H2C=O b-V     45. e-1  
##  2 a-T   b-P   b-H2C=O a-V    100. e-1  
##  3 b-T   a-P   b-H2C=O a-V     45. e-1  
##  4 a-T   a-P   b-H2C=O b-V     65. e-1  
##  5 b-T   b-P   a-H2C=O a-V     75. e-1  
##  6 a-T   b-P   a-H2C=O b-V     60. e-1  
##  7 b-T   a-P   a-H2C=O b-V     80. e-1  
##  8 a-T   a-P   a-H2C=O a-V     96. e-1  
##  9 b-T   b-P   b-H2C=O a-V     43. e-2  
## 10 a-T   b-P   b-H2C=O b-V     71. e-2  
## 11 b-T   a-P   b-H2C=O b-V     48. e-2  
## 12 a-T   a-P   b-H2C=O a-V    104. e-2  
## 13 b-T   b-P   a-H2C=O b-V     68. e-2  
## 14 a-T   b-P   a-H2C=O a-V     86. e-2  
## 15 b-T   a-P   a-H2C=O a-V     70. e-2  
## 16 a-T   a-P   a-H2C=O b-V     65. e-2

Gráfica exploratoria 2

library(ggplot2)
g1 <- ggplot(tasaFilt2, aes(A, TF, col = B))
g1 + geom_point(size = 5) + facet_grid(C ~ D)

Modelización

\[ y_{ijklmo} = \mu + \alpha_i + \beta_j + \gamma_k + \delta_l + \kappa_{o}\\ (\alpha\beta)_{ij} + (\alpha\gamma)_{ik} + (\alpha\delta)_{il} + \\ \epsilon_{ijklmo} \]

Donde: \[ \epsilon_{ijklmo} \sim \mathcal{N}(0, \sigma^2) \quad i.i.d. \]

modelo4 <- aov(TF ~ A + B + C + D + A*B*C*D + Etapa, data = tasaFilt2)
summary(modelo4)
##             Df Sum Sq Mean Sq
## A            1 1870.6  1870.6
## B            1   39.1    39.1
## C            1  390.1   390.1
## D            1  855.6   855.6
## Etapa        1    7.6     7.6
## A:B          1    0.1     0.1
## A:C          1 1314.1  1314.1
## B:C          1   22.6    22.6
## A:D          1 1105.6  1105.6
## B:D          1    0.6     0.6
## C:D          1    5.1     5.1
## A:B:C        1   14.1    14.1
## A:B:D        1   68.1    68.1
## A:C:D        1   10.6    10.6
## B:C:D        1   27.6    27.6

Determinación de interacciones significativas

library(pid)
paretoPlot(modelo4)

efectos <- unlist(model.tables(modelo4))
efectos <- efectos[substring(names(efectos),1,6)=="tables"]
names(efectos)<-gsub("tables.","",names(efectos))
qqR<-qqnorm(efectos, xlim = c(-3, 3), ylim = c(-15,15), pch = 19, cex = 0.5)
qqline(efectos)
with(qqR,text(x,y,labels=names(efectos),cex=0.7,pos=3))

modelo5 <- aov(TF ~ A + B + C + D  + A*B*C*D + Etapa - A:B:C:D,
               data = tasaFilt2)
summary(modelo5)
##             Df Sum Sq Mean Sq
## A            1 1870.6  1870.6
## B            1   39.1    39.1
## C            1  390.1   390.1
## D            1  855.6   855.6
## Etapa        1    7.6     7.6
## A:B          1    0.1     0.1
## A:C          1 1314.1  1314.1
## B:C          1   22.6    22.6
## A:D          1 1105.6  1105.6
## B:D          1    0.6     0.6
## C:D          1    5.1     5.1
## A:B:C        1   14.1    14.1
## A:B:D        1   68.1    68.1
## A:C:D        1   10.6    10.6
## B:C:D        1   27.6    27.6
modelo6 <- aov(TF ~ A + B + C + D  + A*B*C*D + Etapa - A:B:C:D - A:C:D,
               data = tasaFilt2)
summary(modelo6)
##             Df Sum Sq Mean Sq F value Pr(>F)  
## A            1 1870.6  1870.6 177.095 0.0477 *
## B            1   39.1    39.1   3.698 0.3053  
## C            1  390.1   390.1  36.929 0.1038  
## D            1  855.6   855.6  81.000 0.0704 .
## Etapa        1    7.6     7.6   0.716 0.5529  
## A:B          1    0.1     0.1   0.006 0.9511  
## A:C          1 1314.1  1314.1 124.408 0.0569 .
## B:C          1   22.6    22.6   2.136 0.3820  
## A:D          1 1105.6  1105.6 104.669 0.0620 .
## B:D          1    0.6     0.6   0.053 0.8556  
## C:D          1    5.1     5.1   0.479 0.6145  
## A:B:C        1   14.1    14.1   1.331 0.4546  
## A:B:D        1   68.1    68.1   6.444 0.2389  
## B:C:D        1   27.6    27.6   2.609 0.3529  
## Residuals    1   10.6    10.6                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Modelización definitiva 2

modelo6 <- aov(TF ~ A + B + C + D  + A*B*C*D + Etapa -
                A:B:C:D - A:C:D - A:B:C - B:C:D - A:B:D -
                A:B - B:D - C:D - B:C - Etapa - B,
               data = tasaFilt2)
summary(modelo6)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## A            1 1870.6  1870.6   95.86 1.93e-06 ***
## C            1  390.1   390.1   19.99   0.0012 ** 
## D            1  855.6   855.6   43.85 5.92e-05 ***
## A:C          1 1314.1  1314.1   67.34 9.41e-06 ***
## A:D          1 1105.6  1105.6   56.66 2.00e-05 ***
## Residuals   10  195.1    19.5                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Diagnósticos

Normalidad

# Normalidad
residualesStud <- rstudent(modelo6)
shapiro.test(residualesStud)
## 
##  Shapiro-Wilk normality test
## 
## data:  residualesStud
## W = 0.95297, p-value = 0.5381

Igualdad de varianza

# Igualdad de varianza
bartlett.test(TF ~ A, data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by A
## Bartlett's K-squared = 0.12214, df = 1, p-value = 0.7267
bartlett.test(TF ~ B, data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by B
## Bartlett's K-squared = 0.040578, df = 1, p-value = 0.8404
bartlett.test(TF ~ C, data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by C
## Bartlett's K-squared = 3.3143, df = 1, p-value = 0.06868
bartlett.test(TF ~ D, data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by D
## Bartlett's K-squared = 3.0781, df = 1, p-value = 0.07935
bartlett.test(TF ~ paste(A,C), data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by paste(A, C)
## Bartlett's K-squared = 11.332, df = 3, p-value = 0.01006
bartlett.test(TF ~ paste(A,D), data = tasaFilt2)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  TF by paste(A, D)
## Bartlett's K-squared = 5.0452, df = 3, p-value = 0.1685

Comparaciones de Tukey

library(multcomp)
library(multcompView)
library(sandwich)
tasaFilt2$AC <- with(tasaFilt, interaction(A, C))
modelo6_1 <- aov(TF ~ -1 + AC, data = tasaFilt2)
summary(modelo6_1)
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## AC         4  82115   20529   114.2 1.92e-09 ***
## Residuals 12   2156     180                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
compMultHetero <- glht(modelo6_1, linfct = mcp(AC = "Tukey"), vcov = sandwich)
summary(compMultHetero)
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: aov.default(formula = TF ~ -1 + AC, data = tasaFilt2)
## 
## Linear Hypotheses:
##                                Estimate Std. Error t value Pr(>|t|)    
## a-T.b-H2C=O - b-T.b-H2C=O == 0   39.750      8.641   4.600  0.00237 ** 
## b-T.a-H2C=O - b-T.b-H2C=O == 0   28.000      2.494  11.228  < 0.001 ***
## a-T.a-H2C=O - b-T.b-H2C=O == 0   31.500      7.448   4.229  0.00445 ** 
## b-T.a-H2C=O - a-T.b-H2C=O == 0  -11.750      8.905  -1.320  0.53487    
## a-T.a-H2C=O - a-T.b-H2C=O == 0   -8.250     11.338  -0.728  0.86880    
## a-T.a-H2C=O - b-T.a-H2C=O == 0    3.500      7.752   0.451  0.96313    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
op <- par(no.readonly = TRUE) 
par(mar = op$mar + c(0, 12, 0, 0))
plot(compMultHetero)

par(op)

Comportamiento del modelo

tasaFilt2$TFp <- predict(modelo6)
g1 <- ggplot(tasaFilt2, aes(A, TF, col = C))
g1 + geom_point() + facet_grid(. ~ D) +
  geom_line(aes(A, TFp, group = C)) +
  geom_point(aes(A, TFp, col = C), shape = 17)

Valor de los efectos

model.tables(modelo6, se = TRUE)
## Tables of effects
## 
##  A 
## A
##     b-T     a-T 
## -10.812  10.812 
## 
##  C 
## C
## b-H2C=O a-H2C=O 
##  -4.938   4.938 
## 
##  D 
## D
##    b-V    a-V 
## -7.313  7.313 
## 
##  A:C 
##      C
## A     b-H2C=O a-H2C=O
##   b-T -9.062   9.062 
##   a-T  9.062  -9.063 
## 
##  A:D 
##      D
## A     b-V    a-V   
##   b-T  8.312 -8.312
##   a-T -8.312  8.312
## 
## Standard errors of effects
##             A     C     D   A:C   A:D
##         1.562 1.562 1.562 2.209 2.209
## replic.     8     8     8     4     4
model.tables(modelo6, se = TRUE, type = "mean")
## Tables of means
## Grand mean
##         
## 70.0625 
## 
##  A 
## A
##   b-T   a-T 
## 59.25 80.88 
## 
##  C 
## C
## b-H2C=O a-H2C=O 
##   65.12   75.00 
## 
##  D 
## D
##   b-V   a-V 
## 62.75 77.38 
## 
##  A:C 
##      C
## A     b-H2C=O a-H2C=O
##   b-T 45.25   73.25  
##   a-T 85.00   76.75  
## 
##  A:D 
##      D
## A     b-V   a-V  
##   b-T 60.25 58.25
##   a-T 65.25 96.50
## 
## Standard errors for differences of means
##             A     C     D   A:C   A:D
##         2.209 2.209 2.209 3.123 3.123
## replic.     8     8     8     4     4

Estimación de los errores

errsdred(modelo6)
## [1] 4.417296

Generación de diseños.

Diseños fraccionados.

Diseño \(2^{4-1}_{IV}\).

library(FrF2)
FrF2(8, 4)
##    A  B  C  D
## 1  1 -1  1 -1
## 2  1  1 -1 -1
## 3 -1 -1  1  1
## 4 -1  1 -1  1
## 5 -1 -1 -1 -1
## 6  1 -1 -1  1
## 7  1  1  1  1
## 8 -1  1  1 -1
## class=design, type= FrF2

Diseño \(2^{8-2}_{V}\).

library(FrF2)
FrF2(64, 8)
##     A  B  C  D  E  F  G  H
## 1   1 -1  1 -1  1 -1  1  1
## 2   1  1 -1 -1 -1  1  1 -1
## 3  -1 -1 -1 -1 -1  1  1 -1
## 4  -1 -1 -1 -1  1 -1  1 -1
## 5   1  1 -1  1  1  1 -1  1
## 6   1 -1 -1 -1 -1  1 -1  1
## 7   1  1 -1 -1  1  1  1  1
## 8  -1  1  1 -1  1  1  1 -1
## 9  -1  1 -1 -1  1  1 -1 -1
## 10  1  1  1 -1  1  1 -1  1
## 11 -1  1  1 -1 -1 -1  1 -1
## 12  1 -1  1 -1  1  1  1 -1
## 13  1  1 -1 -1  1 -1  1 -1
## 14 -1  1 -1  1  1 -1  1  1
## 15  1 -1 -1  1 -1  1  1  1
## 16  1  1 -1  1 -1  1 -1 -1
## 17  1 -1  1 -1 -1 -1  1 -1
## 18 -1 -1  1  1 -1 -1  1  1
## 19  1  1  1 -1  1 -1 -1 -1
## 20 -1 -1 -1  1 -1  1 -1 -1
## 21  1 -1  1 -1 -1  1  1  1
## 22 -1 -1 -1  1  1 -1 -1 -1
## 23 -1  1  1  1  1 -1 -1  1
## 24 -1 -1 -1  1  1  1 -1  1
## 25  1  1  1  1 -1  1  1 -1
## 26 -1 -1  1 -1  1  1 -1  1
## 27  1 -1  1  1 -1 -1 -1 -1
## 28 -1  1  1  1 -1  1 -1  1
## 29 -1  1 -1 -1 -1  1 -1  1
## 30 -1  1 -1  1  1  1  1 -1
## 31  1 -1 -1  1  1 -1  1  1
## 32 -1  1 -1  1 -1  1  1  1
## 33  1  1  1 -1 -1  1 -1 -1
## 34  1 -1 -1 -1  1 -1 -1  1
## 35  1 -1 -1  1  1  1  1 -1
## 36 -1  1  1 -1  1 -1  1  1
## 37  1  1  1 -1 -1 -1 -1  1
## 38 -1  1  1  1 -1 -1 -1 -1
## 39 -1 -1  1 -1  1 -1 -1 -1
## 40 -1 -1  1  1 -1  1  1 -1
## 41  1  1  1  1  1  1  1  1
## 42 -1  1  1 -1 -1  1  1  1
## 43 -1 -1  1  1  1 -1  1 -1
## 44 -1 -1  1  1  1  1  1  1
## 45 -1  1 -1 -1 -1 -1 -1 -1
## 46  1  1 -1  1  1 -1 -1 -1
## 47  1 -1 -1 -1  1  1 -1 -1
## 48  1 -1  1  1  1 -1 -1  1
## 49  1  1 -1  1 -1 -1 -1  1
## 50  1 -1  1  1  1  1 -1 -1
## 51 -1 -1 -1  1 -1 -1 -1  1
## 52  1  1 -1 -1 -1 -1  1  1
## 53 -1  1 -1  1 -1 -1  1 -1
## 54  1 -1 -1 -1 -1 -1 -1 -1
## 55  1 -1 -1  1 -1 -1  1 -1
## 56 -1 -1 -1 -1 -1 -1  1  1
## 57 -1 -1 -1 -1  1  1  1  1
## 58 -1 -1  1 -1 -1  1 -1 -1
## 59 -1  1 -1 -1  1 -1 -1  1
## 60  1  1  1  1 -1 -1  1  1
## 61 -1 -1  1 -1 -1 -1 -1  1
## 62  1 -1  1  1 -1  1 -1  1
## 63 -1  1  1  1  1  1 -1 -1
## 64  1  1  1  1  1 -1  1 -1
## class=design, type= FrF2

Diseño de Plackett-Burman.

El diseño de Plackett-Burman es un diseño de dos niveles de resolución III. Recordemos que en este diseño se supone que las interacciones dobles no son significativas porque están aliadas con los efectos principales.

¿Cuándo usar el diseño Plackett-Burman?:

  • En selección o tamizado de variables.
  • Cuando las interacciones de cualquier orden se consideran no significativas.
  • Experimentos multifactorial de dos niveles cada factor.
  • Se usa cuando hay más de cuatro factores.
  • Se buscan los factores principales más relevantes.
  • El número de experimentos deberá ser de 12, 20, 24, 28, 36, … (Múltiplos de 4).
library(FrF2)
pb(12, 11)
##     A  B  C  D  E  F  G  H  J  K  L
## 1  -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1
## 2   1 -1  1  1 -1  1  1  1 -1 -1 -1
## 3  -1  1 -1  1  1 -1  1  1  1 -1 -1
## 4   1  1  1 -1 -1 -1  1 -1  1  1 -1
## 5   1 -1  1  1  1 -1 -1 -1  1 -1  1
## 6  -1  1  1  1 -1 -1 -1  1 -1  1  1
## 7  -1  1  1 -1  1  1  1 -1 -1 -1  1
## 8   1  1 -1 -1 -1  1 -1  1  1 -1  1
## 9  -1 -1 -1  1 -1  1  1 -1  1  1  1
## 10 -1 -1  1 -1  1  1 -1  1  1  1 -1
## 11  1  1 -1  1  1  1 -1 -1 -1  1 -1
## 12  1 -1 -1 -1  1 -1  1  1 -1  1  1
## class=design, type= pb