Skip to content

Factorial Experiments

"Factorial experiment" is a general term used for experiments with multiple factors, designed so that all possible combinations of levels of factors are tested. In theory, such experiments could be quite large. The number of possible combinations of levels is given by the product of the number of levels in each factor. For example if

  • then the total number of combinations of levels is .

Note that this number grows quite quickly... Even if only factors are considered, each with just levels, the total number of combinations would be . Due to this inflation, we will mostly only consider factorial experiments with (or at most ) levels in each factor.

Observations are called replicates if they involve the same combination of levels/treatments of each factor.

Notation
  • Upper case letters for factors,
    e.g. , , when there are few factors
    and when there are many factors
  • Lower case letters for the number of levels/treatments in each factor
    e.g. or
  • Numbers for individual levels/treatments in factors
    e.g. , or etc.
    If there are only two levels, then we may use different notation like or else (see the section on factorial experiments)

  • total number of observations
  • number of replicate observations for each combination of levels.
        So .
  • will be indices for factors (unless we use for total number of factors)

  • Observations are written as follows (for a three factor experiment):
    • replicate observation with factor levels , , and
          i.e. first factor has level , second factor has level , etc.
  • Sample means are as follows (for experiment with factors , , ):
    • sample mean of all replicates with , ,
          This is the mean of the cell
    • sample mean of all observations with , and
          This is the mean of the line
    • sample mean of all observations with
          This is the mean of the plane
    • etc.

We mentally organize data into a -dimensional array (with one dimension for each factor) of lists (with one element for each replicate) of observations.

Interaction between factors

If the effects of changing two factors are not independent, then we say there is interaction between them. Factorial experiments are designed to look for interaction between factors.

The factorial model assumes that expected values for each combination of factor values can be written as the overall mean, plus a main (treatment) effect of changing individual factors, plus interaction effects of changing combinations of factors.

Our general model looks something like the following:

  • "Model:"
    • is the over-all mean
    • are the "main treatment effects" for the different levels of corresponding factors
    • are the "interaction effects" of the corresponding factors
    • is the error, distributed as .

If two factors have no interaction, then when we plot their sample means we should get a series of parallel shapes. If there is interaction, then slopes will be unequal.

### Example: two factors
A <- c(1,2,3)    # Levels of factor A
B <- c(1,2,3,4)  # Levels of factor B

a <- length(A)
b <- length(B)

## choose mean and treatment effects
mu <- runif(1, -2,2)   # mean
tA <- runif(a, -4,4)   # treatment A effects
tB <- runif(b, -4,4)   # treatment B effects

## combine mean and treatment effects for each level 
Ybar <- mu + 
        matrix(rep(tA,b), nrow=a, byrow=FALSE) +
        matrix(rep(tB,a), nrow=a, byrow=TRUE)

cat('means of experiment (no interaction)')
Ybar

## make factor indicators for each data point
fA <- matrix(rep(A,b), nrow=a, byrow=FALSE)
fB <- matrix(rep(B,a), nrow=a, byrow=TRUE)

## check that indicators are correct :)
cat('factor A indicator')
fA
cat('factor B indicator')
fB

## plot resulting means (no interaction)
plot(fA,Ybar,main="No interaction (parallel)")

## use lines to show B factor value
for (j in 1:b) {
	lines(fA[,j],Ybar[,j], lty=j)
}

## now let's add an interaction term!

tAB <- matrix(runif(a*b, -4,4), nrow=a)

Ybar <- Ybar + tAB

## plot resulting means (with interaction)
plot(fA,Ybar,main="With interaction (not parallel)")

## use lines to show B factor value
for (j in 1:b) {
	lines(fA[,j],Ybar[,j], lty=j)
}

Two Factor Experiments

Example: Observations for 2 factor experiment
FactorB
Level123
A1 (1,1) observations (1,2) observations (1,3) observations
2 (2,1) observations (2,2) observations (2,3) observations
3 (3,1) observations (3,2) observations (3,3) observations

where " observations" is the set of replicates

For example, the upper-right corner of observations for and is

  • The mean of each box is its replicate mean
  • The means of the rows are the means of observations: , ,
  • The means of the columns are the means of observations: , ,
  • The overall mean is
MeansMeans
Means

(This table should remind you of the marginal probability tables that we made for joint random variables!)

Model

  • where

Alternately, the could be absorbed into the yielding the model:

  • where

Point estimators

  • -- the overall mean.

  • -- the main treatment effect at
  • -- the main treatment effect at

  • -- the interaction effect at and .

Comparison to "single factor with blocking"

Note that a two factor experiment is very similar to a single factor experiment with blocking. There are two main differences, however.

  • single factor with blocking experiments had no replicates -- within each level, each observation had a unique block.
  • single factor with blocking experiments were assumed to have "level" and "block" effects independent -- there was no interaction effect

Sum of Squares

Setting variances equal on either side of the model equation leads to the "sum of square" terms.

  • From we get variance of Factor A (row) means:
  • From we get variance of Factor B (column) means:
  • From we get variance of box means not accounted for by Factor A and B (row and column) variances:
  • From we get pooled variance within boxes:
  • On the other side of the equation, gives the variance of all observations:

Sources of Variation

  SourcedfSum of SquareMean Square
Factor A
Factor B
Interaction
Error
Total (*)

(*) Not used.

Note:

F Tests

There are three different -tests associated to the two factor model.

  • Model:

There are two tests for the treatment effects.

  • H0: for all



  • H0: for all

And one test for interaction effect.

  • H0: for all
### Example: two factors (random experiment)
A <- c(1,2,3)    # Levels of factor A
B <- c(1,2,3,4)  # Levels of factor B

a <- length(A)
b <- length(B)

## choose random mean and treatment effects
mu  <- runif(1, -2,2)   # overall mean
tA  <- runif(a, -2,2)   # treatment A effects
tB  <- runif(b, -2,2)   # treatment B effects
tAB <- runif(a*b, -2,2) # interaction effects 

## combine mean and treatment effects for each level 
Ybar <- mu + 
        matrix(rep(tA,b), nrow=a, byrow=FALSE) +
        matrix(rep(tB,a), nrow=a, byrow=TRUE)  +
		matrix(tAB,       nrow=a)

## make factor indicators for each data point
fA <- matrix(rep(A,b), nrow=a, byrow=FALSE)
fB <- matrix(rep(B,a), nrow=a, byrow=TRUE)

## now add replicates for each level combination
n  <- 50

y  <- as.vector(rep(Ybar, n)) + rnorm(a*b*n, 0,5) # add error
fA <- factor(   rep(fA,   n))
fB <- factor(   rep(fB,   n))

data.AB <- data.frame(obs=y, A=fA, B=fB) 
head(data.AB,10)
Example: Testing for interactions by hand!
### Let's compute the source of variance table by hand!!!!
ybar_A  <- aggregate(obs ~ A,     data=data.AB, mean)
ybar_B  <- aggregate(obs ~ B,     data=data.AB, mean)
ybar_AB <- aggregate(obs ~ A + B, data=data.AB, mean)

## sum of squares A and B are just var of row and col means
SS_A  <- n*b * (a-1) * sd( ybar_A$obs )^2
SS_B  <- n*a * (b-1) * sd( ybar_B$obs )^2

## sum of square interaction is "unexplained" variance of box means
SS_AB <- n  * (a*b-1)* sd( ybar_AB$obs )^2 - SS_A - SS_B

## total sum of squares is just variance of data
SS_T  <- (a*b*n-1)   * sd( y )^2

## sum of square error is the "unexplained" part of total
SS_E  <- SS_T - SS_A - SS_B - SS_AB

## convert to mean square
MS_A  <-  SS_A / (a-1)
MS_B  <-  SS_B / (b-1)
MS_AB <- SS_AB / ((a-1)*(b-1))
MS_E  <-  SS_E / (a*b*(n-1))
MS_T  <-  SS_T / (a*b*n-1)

## F-statistics and p-values!
F_A  <-  MS_A / MS_E
F_B  <-  MS_B / MS_E
F_AB <- MS_AB / MS_E

p_A  <- 1 - pf(F_A, a-1, a*b*(n-1))
p_B  <- 1 - pf(F_B, b-1, a*b*(n-1))
p_AB <- 1 - pf(F_AB, (a-1)*(b-1), a*b*(n-1))

## pretty table!
data.frame(
	Source=c("Factor A","Factor B","Interaction", "Error", "Total"),
	df    =c(   a-1,       b-1,     (a-1)*(b-1), a*b*(n-1),a*b*n-1),
	SS    =c(  SS_A,      SS_B,        SS_AB,      SS_E,    SS_T  ),
	MS    =c(  MS_A,      MS_B,        MS_AB,      MS_E,    MS_T  ),
	F     =c(  F_A,        F_B,         F_AB,        '',     ''   ),
	p     =c(  p_A,        p_B,         p_AB,        '',     ''   )
)
Example: Testing for interactions using R commands
## This is very quick and easy to do in R...
##  interactions terms are indicated by multiplication

summary( aov(obs ~ A + B + A*B, data=data.AB) )
Example: Testing for interactions that don't exist

If interactions don't exist then we should be able to see this from the corresponding -value. In this case, we can remove the interaction effect and rerun the anova to look more carefully at the treatment effects of the other variables.

An interesting thing happens when we remove treatment or interaction effects from our model. The corresponding row of the "Sources of Variance" table is removed and added to the Error row, without changing any of the other rows! This is called model collapse.

## What if there are no interactions??
### Example: two factors (random experiment)
A <- c(1,2,3)    # Levels of factor A
B <- c(1,2,3,4)  # Levels of factor B

a <- length(A)
b <- length(B)

## choose random mean and treatment effects
mu  <- runif(1, -2,2)   # overall mean
tA  <- runif(a, -2,2)   # treatment A effects
tB  <- runif(b, -2,2)   # treatment B effects
tAB <- rep(0, a*b   )   # NO INTERACTION EFFECTS!!!!

## combine mean and treatment effects for each level 
Ybar <- mu + 
        matrix(rep(tA,b), nrow=a, byrow=FALSE) +
        matrix(rep(tB,a), nrow=a, byrow=TRUE)  +
		matrix(tAB,       nrow=a)

## make factor indicators for each data point
fA <- matrix(rep(A,b), nrow=a, byrow=FALSE)
fB <- matrix(rep(B,a), nrow=a, byrow=TRUE)

## now add replicates for each level combination
n  <- 50

y  <- as.vector(rep(Ybar, n)) + rnorm(a*b*n, 0,5) # add error
fA <- factor(   rep(fA,   n))
fB <- factor(   rep(fB,   n))

data.noAB <- data.frame(obs=y, A=fA, B=fB) 

cat('aov test with interactions')
summary( aov(obs ~ A + B + A*B, data=data.noAB) )

cat('aov test without interactions')
summary( aov(obs ~ A + B,       data=data.noAB) )

Two factor experiments with no replicates ("one observation per cell")

  • This is also called an "unreplicated experiment"

DO NOT DO THIS with only two factors!!!!!

In this case, the degrees of freedom of the error variance is . This is bad!

If you assume that there is no interaction effect, then the interaction row will collapse onto the error row, so that we get some degrees of freedom. This is equivalent to one factor with blocking:

Model

Sources of Variation (no replicates, no interaction)

  SourcedfSum of SquareMean Square
Factor A
Factor B
Error
Total (*)

(*) if then .

... but one of the main points of factorial experiments is to look for interaction effects. So this is stupid.

Summary:

Without replicates, you cannot check for interaction effects in two factor experiments!

Three Factor Experiments

... in this case, data is a 3d array of lists of replicates.

Model

As well as the three main treatment effects, there are now three pair-wise interaction effects as well as one triple interaction effect.

Note that if there are factors, then there will be interaction effects. This number grows very quickly. Once there are a lot of factors, we may be satisfied with only checking as far as pair-wise interaction effects (that way we only have interaction terms).

Point estimators

Point estimators of main effects and interactions are given by alternating sums:







Sum of Squares

Setting variances equal on either side of the model equation leads to the "sum of square" terms from the point estimators. For example:



Sources of Variation

  SourcedfSum of SquareMean Square
Factor A
Factor B
Factor C
Mix AB
Mix BC
Mix AC
Mix ABC
Error
Total (*)

(*) Not used.

Now there are seven -tests. Three for the main treatment effects ; three for the two-factor interactions ; and one for the three-factor interaction .

Example. Random three factor experiment

First I'll generate treatment effects and data.

### Example: three factors (random experiment)
A <- c(1,2,3)    # Levels of factor A
B <- c(1,2,3,4)  # Levels of factor B
C <- c(1,2,3,4)  # Levels of factor C

a <- length(A)
b <- length(B)
c <- length(C)

## choose random mean and treatment effects
mu  <- runif(1, -2,2)  # overall mean
tA  <- runif(a, -2,2)   # treatment A main effects
tB  <- runif(b, -2,2)   # treatment B main effects
tC  <- runif(c, -2,2)   # treatment C main effects
tAB <- runif(a*b, -2,2)  # AB interaction effects
tAC <- runif(a*c, -2,2)  # AC interaction effects
tBC <- runif(b*c, -2,2)  # BC interaction effects
tABC<- runif(a*b*c, -2,2) # ABC interaction effects

## combine mean and treatment effects 
##  use array(.., dim=c(..)) to arrange into 3d arrays
##  use aperm(.., c(..)) to reorder dimensions of arrays
Ybar <- mu + 
              array(tA, dim=c(a,b,c))            +
		aperm(array(tB, dim=c(b,c,a)), c(3,1,2)) +
		aperm(array(tC, dim=c(c,a,b)), c(2,3,1)) +
		      array(tAB, dim=c(a,b,c))            +
		aperm(array(tAC, dim=c(a,c,b)), c(1,3,2)) +
        aperm(array(tBC, dim=c(b,c,a)), c(3,1,2)) +
		      array(tABC, dim=c(a,b,c))

## compare to simpler versions missing some interactions
Ybar.noABC <- mu + 
              array(tA, dim=c(a,b,c))            +
		aperm(array(tB, dim=c(b,c,a)), c(3,1,2)) +
		aperm(array(tC, dim=c(c,a,b)), c(2,3,1)) +
		      array(tAB, dim=c(a,b,c))            +
		aperm(array(tAC, dim=c(a,c,b)), c(1,3,2)) +
        aperm(array(tBC, dim=c(b,c,a)), c(3,1,2))

Ybar.noBC <-  mu + 
              array(tA, dim=c(a,b,c))            +
		aperm(array(tB, dim=c(b,c,a)), c(3,1,2)) +
		aperm(array(tC, dim=c(c,a,b)), c(2,3,1)) +
		      array(tAB, dim=c(a,b,c))            +
		aperm(array(tAC, dim=c(a,c,b)), c(1,3,2)) 

Ybar.noB <- mu + 
              array(tA, dim=c(a,b,c))            +
		aperm(array(tC, dim=c(c,a,b)), c(2,3,1)) +
		aperm(array(tAC, dim=c(a,c,b)), c(1,3,2))

Ybar.onlyAC <- mu + 
              array(tA, dim=c(a,b,c))            +
		aperm(array(tB, dim=c(b,c,a)), c(3,1,2)) +
		aperm(array(tC, dim=c(c,a,b)), c(2,3,1)) +
		aperm(array(tAC, dim=c(a,c,b)), c(1,3,2))


## make factor indicators for each data point
fA <-       array(A, dim=c(a,b,c))
fB <- aperm(array(B, dim=c(b,c,a)), c(3,1,2))
fC <- aperm(array(C, dim=c(c,a,b)), c(2,3,1))

## now add replicates
n  <- 50

## random error for replicates
error <- rnorm(a*b*c*n, 0,5)

## observations
y  <- as.vector(rep(Ybar, n)) + error 

## also add replicates for reduced versions
y.noABC  <- as.vector(rep(Ybar.noABC, n))  + error 
y.noBC   <- as.vector(rep(Ybar.noBC, n))   + error 
y.noB    <- as.vector(rep(Ybar.noB, n))    + error 
y.onlyAC <- as.vector(rep(Ybar.onlyAC, n)) + error 

## now add replicates for factors
fA <- factor(rep(fA, n))
fB <- factor(rep(fB, n))
fC <- factor(rep(fC, n))

## combine everything into a big data frame!
data.ABC <- data.frame(
	               obs    = y,
				   noABC  = y.noABC,
				   noBC   = y.noBC,
				   noB    = y.noB,
				   onlyAC = y.onlyAC,
				   A=fA, B=fB, C=fC
)

Now let's look at some anova tables!