2^k Factorial Designs

STA6246

Author

Dr. Cohen

Catalyst of Chemical Process

  • An experimenter wanted to study the Yield of a chemical process

  • The engineer decides to test 2 reactant concentrations (15 vs 25) and 2 catalysts (1 vs 2).

  • 3 replicates were run

Data description

library(tidyverse)
library(gtsummary)
# 2^2 Factorial design
# 2^2 Data Yield in a chemical process

Y=c(28,36,18,31,25,32,19,30,27,32,23,29)
A=gl(n=2,k=1,length = 4,labels = c("L","H"))
B=gl(n=2,k=2,length = 4,labels = c("Lo","Hi"))

# block
bl=gl(n=3,k=4,length = 12,labels=c(1,2,3))
# data
mdata=data.frame(Y,A,B)

# ANOVA
# EA= 8.33
# EB = -5
model1 <- aov(Y ~ A+B, data = mdata)
summary(model1)
            Df Sum Sq Mean Sq F value   Pr(>F)    
A            1 208.33  208.33   47.27 7.27e-05 ***
B            1  75.00   75.00   17.02  0.00258 ** 
Residuals    9  39.67    4.41                     
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Computing Effects

Estimate the effects from the regression

The contrasts

getOption('contrasts')
        unordered           ordered 
"contr.treatment"      "contr.poly"

Orthogonal - effect coding - half-effects

options(contrasts=c("contr.sum","contr.poly"))
contrasts(A)
  [,1]
L    1
H   -1
m2=lm(Y ~ A*B, data = mdata)
coef(m2)
(Intercept)          A1          B1       A1:B1 
 27.5000000  -4.1666667   2.5000000   0.8333333 
# change reference level
co=coef(lm(Y ~ relevel(A, ref = "H")*relevel(B, ref = "Hi"),data = mdata))

# Effects
Eff=2*co[2:4]
names(Eff) = c("A","B","AB")
Eff
        A         B        AB 
 8.333333 -5.000000  1.666667

Reference coding - effects

It works only with no interactions

options (contrasts = c("contr.treatment", "contr.poly"))
contrasts(A)
  H
L 0
H 1
m3 = lm(Y ~ A+ B ,data=mdata)
#By default The intercept: overall mean - (EffectA+ EffectB)/2
# and the regression coefficients are the Effect Estimates

Eff2=coef(m3)[2:3]
names(Eff2) <- c("A","B") 
Eff2
        A         B 
 8.333333 -5.000000 
# Graphics
interaction.plot(x.factor = mdata$A, 
                 trace.factor = mdata$B, 
                 response = mdata$Y, 
                 fun = mean, 
                 type = "b", 
                 col=c("black" ,"red"),
                 fixed = T,
                 pch = c(19,17))

Etch Rate data

An etch process on a single wafer plasma etching tool. The design factors:

  • The gap between the electrodes (A)

  • The gas flow (B)

  • RF power applied to the cathode. (C)

  • Each factor is run at two levels

  • The design is replicated twice n=2.

Etch=c(550,669,633,642,1037,749,1075,729,604,650,601,635,1052,868,1063,860)
A=gl(n=2,k=1,length = 16,labels = c("L","H"))
B=gl(n=2,k=2,length = 16,labels = c("L","H"))
C=gl(n=2,k=4,length = 16,labels = c("L","H"))
df=tibble(Etch,A,B,C)


model1 <- aov(Etch ~ A*B*C, data = df)
summary(model1)
            Df Sum Sq Mean Sq F value   Pr(>F)    
A            1  41311   41311  18.339 0.002679 ** 
B            1    218     218   0.097 0.763911    
C            1 374850  374850 166.411 1.23e-06 ***
A:B          1   2475    2475   1.099 0.325168    
A:C          1  94403   94403  41.909 0.000193 ***
B:C          1     18      18   0.008 0.930849    
A:B:C        1    127     127   0.056 0.818586    
Residuals    8  18020    2253                     
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Get estimates
getOption("contrasts")
[1] "contr.treatment" "contr.poly"     
options(contrasts=c("contr.sum","contr.poly")) # provide half-effects
co=coef(lm(Etch ~ relevel(A, ref = "H")*relevel(B, ref = "H")*relevel(C, ref = "H")))
# Effects
Eff=2*co[2:8]
names(Eff) <- c("A","B","C","AB","AC","BC","ABC") 
Eff
       A        B        C       AB       AC       BC      ABC 
-101.625    7.375  306.125  -24.875 -153.625   -2.125    5.625 
qq=qqnorm(Eff,cex=1,ylim = c(-200,300),xlim=c(-3,1.5))
qqline(Eff)
text(qq$x,qq$y,names(Eff),pos=1,cex=1,col="red")

  • The factors A, C, and interactions AC are the most important terms.