Note: Revised with revised data
This may be somewhat of an R-specific answer.
By default, R uses Type I sums of squares. This is fine for balanced designs or when sequential sums of squares are desired, but people often want to use Type II or Type III sums of squares. As a matter of course I recommend using the lm function which fits a general linear model, and the Anova function in the car package for the anova table. This function defaults to Type II sums of squares, which is a good choice for a default. Another advantage of this approach is that it allows you to use the emmeans package for post-hoc comparisons.
In the case of your data, the effects are not completely crossed. That is, Enzymes AH and AL are measured only in Strain 1. This will create some problems in estimating effects.
The following code will run at: rdrr.io/snippets/, or in R.
Install packages and assemble data
if(!require(car)){install.packages("car")}
if(!require(FSA)){install.packages("FSA")}
if(!require(emmeans)){install.packages("emmeans")}
Butanol =c(11.7462,11.7904,11.9162,11.8732,11.8583,11.8677,11.8697,
11.7289,11.9296,11.7722,11.9813,11.9873,11.8058,11.8711,11.937,
11.8628,11.7786,11.7649,11.8459,11.9139,12.0537,12.1359,11.9949,
12.0752,11.9993,12.2802,12.2227,12.1274,12.1408,11.9896,12.1362,
12.265,12.1353,12.0812,12.511,12.1871,11.7881,12.1962,12.2482,12.189)
Strain = c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,
1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2)
Enzyme = c("AH","AH","AH","AH","AH","AA","AA","AA","AA","AA",
"AL","AL","AL","AL","AL","B","B","B","B","B","C","C","C","C","C",
"AA","AA","AA","AA","AA","B","B","B","B","B","C","C","C","C","C")
myData = data.frame(Butanol, Strain, Enzyme)
myData$Strain = factor(myData$Strain)
Summarize data. Note that Enzyme AH and AL are measured only in Strain 1.
library(FSA)
Summarize(Butanol ~ Strain + Enzyme, myData)
# Strain Enzyme n mean sd min Q1 median Q3 max
# 1 1 AA 5 11.83362 0.0812620 11.7289 11.7722 11.8677 11.8697 11.9296
# 2 2 AA 5 12.15214 0.1101567 11.9896 12.1274 12.1408 12.2227 12.2802
# 3 1 AH 5 11.83686 0.0679207 11.7462 11.7904 11.8583 11.8732 11.9162
# 4 1 AL 5 11.91650 0.0773750 11.8058 11.8711 11.9370 11.9813 11.9873
# 5 1 B 5 11.83322 0.0616360 11.7649 11.7786 11.8459 11.8628 11.9139
# 6 2 B 5 12.22574 0.1732074 12.0812 12.1353 12.1362 12.2650 12.5110
# 7 1 C 5 12.05180 0.0583478 11.9949 11.9993 12.0537 12.0752 12.1359
# 8 2 C 5 12.12172 0.1881807 11.7881 12.1871 12.1890 12.1962 12.2482
Fit a general linear model. Despite the fact that R returns an anova table, caution should be used since the design is not fully crossed.
model = lm(Butanol ~ Strain*Enzyme, myData)
library(car)
Anova(model)
# Anova Table (Type II tests)
#
# Response: Butanol
# Sum Sq Df F value Pr(>F)
# Strain 0.50825 1 39.8805 4.373e-07 ***
# Enzyme 0.06587 4 1.2921 0.293837
# Strain:Enzyme 0.14279 2 5.6021 0.008203 **
# Residuals 0.40782 32
Compare estimated marginal means. Note some cannot be estimated.
library(emmeans)
marginal = emmeans(model, ~ Strain:Enzyme)
pairs(marginal)
CLD(marginal, Letters=letters)
# Strain Enzyme emmean SE df lower.CL upper.CL .group
# 1 B 11.8 0.0505 32 11.7 11.9 a
# 1 AA 11.8 0.0505 32 11.7 11.9 a
# 1 AH 11.8 0.0505 32 11.7 11.9 a
# 1 AL 11.9 0.0505 32 11.8 12.0 ab
# 1 C 12.1 0.0505 32 11.9 12.2 abc
# 2 C 12.1 0.0505 32 12.0 12.2 bc
# 2 AA 12.2 0.0505 32 12.0 12.3 bc
# 2 B 12.2 0.0505 32 12.1 12.3 c
# 2 AH nonEst NA NA NA NA
# 2 AL nonEst NA NA NA NA
Also note that the joint_tests function warns of problems.
library(emmeans)
joint_tests(model)
# model term df1 df2 F.ratio p.value note
# Enzyme 2 32 1.757 0.1888 e
# Strain:Enzyme 2 32 5.602 0.0082 e
#
# e: df1 reduced due to non-estimability
Perhaps one idea is to look at Enzymes only within Strains.
library(emmeans)
marginal = emmeans(model, ~ Enzyme|Strain)
pairs(marginal)
CLD(marginal, Letters=letters)
# Strain = 1:
# Enzyme emmean SE df lower.CL upper.CL .group
# B 11.8 0.0505 32 11.7 11.9 a
# AA 11.8 0.0505 32 11.7 11.9 a
# AH 11.8 0.0505 32 11.7 11.9 a
# AL 11.9 0.0505 32 11.8 12.0 ab
# C 12.1 0.0505 32 11.9 12.2 b
#
# Strain = 2:
# Enzyme emmean SE df lower.CL upper.CL .group
# C 12.1 0.0505 32 12.0 12.2 a
# AA 12.2 0.0505 32 12.0 12.3 a
# B 12.2 0.0505 32 12.1 12.3 a
# AH nonEst NA NA NA NA
# AL nonEst NA NA NA NA
#
# Confidence level used: 0.95
# P value adjustment: tukey method for comparing a family of 5 estimates
# significance level used: alpha = 0.05
Check some model assumptions
hist(residuals(model), col="darkgray")
plot(predict(model), residuals(model))
Type III anova appears to not produce output in this case
model3 = lm(Butanol ~ Strain*Enzyme, myData,
contrasts=list(Strain=contr.sum, Enzyme=contr.sum))
library(car)
Anova(model3, type=3)
# Error in Anova.III.lm(mod, error, singular.ok = singular.ok, ...) :
# there are aliased coefficients in the model
Perhaps one approach would be to fit a model with only the interaction as an independent variable. In general, this isn't an ideal approach, but in this case it might be a reasonable approach.
myData$Int = interaction(myData$Strain, myData$Enzyme)
model2 = lm(Butanol ~ Int, myData)
Anova(model2)
# Anova Table (Type II tests)
#
# Response: Butanol
# Sum Sq Df F value Pr(>F)
# Int 0.90294 7 10.121 1.294e-06 ***
# Residuals 0.40782 32
marginal2 = emmeans(model2, ~ Int)
CLD(marginal2, Letters=letters)
# Int emmean SE df lower.CL upper.CL .group
# 1.B 11.8 0.0505 32 11.7 11.9 a
# 1.AA 11.8 0.0505 32 11.7 11.9 a
# 1.AH 11.8 0.0505 32 11.7 11.9 a
# 1.AL 11.9 0.0505 32 11.8 12.0 ab
# 1.C 12.1 0.0505 32 11.9 12.2 abc
# 2.C 12.1 0.0505 32 12.0 12.2 bc
# 2.AA 12.2 0.0505 32 12.0 12.3 c
# 2.B 12.2 0.0505 32 12.1 12.3 c
#
# Confidence level used: 0.95
# P value adjustment: tukey method for comparing a family of 8 estimates
# significance level used: alpha = 0.05
Check some model assumptions
hist(residuals(model2), col="darkgray")
plot(predict(model2), residuals(model2))
