Bac_Temp_BdBsal
Carly Muletz Wolz
3/31/2017
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
setwd("/Users/Carly/Documents/LipsLab/Gabe_project/")
#setwd("/Volumes/Carly/Documents/LipsLab/Gabe_project/")
library(ggplot2)## Warning: package 'ggplot2' was built under R version 3.2.4library(lme4)## Loading required package: Matrix## Warning: package 'Matrix' was built under R version 3.2.4library(gplots)##
## Attaching package: 'gplots'## The following object is masked from 'package:stats':
##
## lowesslibrary(lsmeans)## Loading required package: estimabilitylibrary(car)
library(plyr)
library(pvclust)
Bd_scores <- read.csv("Final_plate_reads/Inhibition_scores_final.csv")
Bd_scores_samples <- Bd_scores[Bd_scores$Bacteria %in% c("B1", "B2", "B3", "B4", "B5"),]
str(Bd_scores_samples)## 'data.frame': 867 obs. of 16 variables:
## $ OD : num 0.00025 -0.00075 -0.00075 -0.00075 0.00225 0.00625 0.00825 -0.00075 -0.00175 -0.00075 ...
## $ Isolate : Factor w/ 6 levels "Bsal","JEL404",..: 5 5 5 5 5 5 5 5 5 5 ...
## $ Genotype : Factor w/ 6 levels "BdBrazil-JEL649",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ Bacteria : Factor w/ 8 levels "B1","B2","B3",..: 4 4 4 4 5 5 5 2 2 2 ...
## $ Temp : int 12 12 12 12 12 12 12 12 12 12 ...
## $ Incubator : int 4 4 4 4 4 4 4 4 4 4 ...
## $ Plate : int 1 1 1 1 1 1 1 1 1 1 ...
## $ Column : int 1 1 1 1 1 1 1 2 2 2 ...
## $ Row : Factor w/ 8 levels "A","B","C","D",..: 1 2 3 4 5 6 8 2 3 4 ...
## $ Day : int 1 1 1 1 1 1 1 1 1 1 ...
## $ T_OD : num 0.00025 -0.00075 -0.00075 -0.00075 0.00225 ...
## $ Well : Factor w/ 96 levels "A1","A10","A11",..: 1 13 25 37 49 61 85 17 29 41 ...
## $ Inhibit_slope : num 0.000331 0.000274 0.000275 0.000064 0.000322 ...
## $ final_inhibition : num 0.112 0.265 0.263 0.829 0.137 ...
## $ R_squared : num 0.69 0.653 0.868 0.658 0.959 ...
## $ final_inhibition_percent: num 11.2 26.5 26.3 82.9 13.7 ...## Temp must be a factor in our analyses
Bd_scores_samples$Temp <- as.factor(Bd_scores_samples$Temp)
### This will be for making the heatmap later on, see code below for explanation
averages <- read.csv("stats_code/Final_code_from_Gabe/heatmap_averages_final.csv", header = T, row.names = 1)Basic data exploration
ggplot(Bd_scores_samples, aes(x = Bacteria, y = final_inhibition_percent)) + geom_boxplot(aes(color = Isolate)) + facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) + scale_x_discrete(labels = c("Ac", "Ps1", "Ps2", "Ch", "St"))
BACTERIA #1, Acinetobacter
B1_Ac <- Bd_scores_samples[Bd_scores_samples$Bacteria == "B1",]
## ~normal
hist(B1_Ac$final_inhibition)
qqp((B1_Ac$final_inhibition), "norm")
## This is what we are analyzing
ggplot(B1_Ac, aes(x = Isolate, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## make sure temp is a factor, after subsetting
class(B1_Ac$Temp)## [1] "factor"## MODEL
## See if temperature treatment or Bd isolate affects Inhibition, while correcting for noise from Incubator or Plate
m1_Ac <- lmer(final_inhibition_percent ~ Temp*Genotype + (1 | Incubator/Plate), data = B1_Ac, REML = F)
summary(m1_Ac)## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula:
## final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Data: B1_Ac
##
## AIC BIC logLik deviance df.resid
## 1740.0 1787.3 -855.0 1710.0 158
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1839 -0.4808 0.0308 0.6239 2.2883
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plate:Incubator (Intercept) 428.6 20.70
## Incubator (Intercept) 0.0 0.00
## Residual 978.1 31.27
## Number of obs: 173, groups: Plate:Incubator, 16; Incubator, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) -34.1871 14.6610 -2.332
## Temp18 -0.9948 19.6863 -0.051
## GenotypeBdGPL1-JEL404 -2.5005 19.5763 -0.128
## GenotypeBdGPL1-JEL647 -17.9798 12.9968 -1.383
## GenotypeBdGPL2-JEL423 -32.4674 19.9422 -1.628
## GenotypeBdGPL2-SRS810 -45.9142 12.9968 -3.533
## GenotypeBsal 79.8033 20.1454 3.961
## Temp18:GenotypeBdGPL1-JEL404 -39.9040 27.0936 -1.473
## Temp18:GenotypeBdGPL1-JEL647 -9.4755 17.1899 -0.551
## Temp18:GenotypeBdGPL2-JEL423 35.1410 27.1768 1.293
## Temp18:GenotypeBdGPL2-SRS810 21.7433 17.3200 1.255
## Temp18:GenotypeBsal 3.7862 27.5077 0.138
##
## Correlation of Fixed Effects:
## (Intr) Temp18 GBGPL1-JEL4 GBGPL1-JEL6 GBGPL2-J GBGPL2-S
## Temp18 -0.745
## GBGPL1-JEL4 -0.749 0.558
## GBGPL1-JEL6 -0.566 0.421 0.424
## GBGPL2-JEL4 -0.735 0.548 0.825 0.416
## GBGPL2-SRS8 -0.566 0.421 0.424 0.638 0.416
## GenotypeBsl -0.728 0.542 0.817 0.412 0.801 0.412
## T18:GBGPL1-JEL4 0.541 -0.727 -0.723 -0.306 -0.596 -0.306
## T18:GBGPL1-JEL6 0.428 -0.512 -0.320 -0.756 -0.314 -0.482
## T18:GBGPL2-J 0.539 -0.724 -0.605 -0.305 -0.734 -0.305
## T18:GBGPL2-S 0.424 -0.508 -0.318 -0.479 -0.312 -0.750
## Tmp18:GntyB 0.533 -0.716 -0.598 -0.301 -0.586 -0.301
## GntypB T18:GBGPL1-JEL4 T18:GBGPL1-JEL6 T18:GBGPL2-J
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4 -0.590
## T18:GBGPL1-JEL6 -0.311 0.372
## T18:GBGPL2-J -0.587 0.817 0.371
## T18:GBGPL2-S -0.309 0.369 0.582 0.368
## Tmp18:GntyB -0.732 0.807 0.366 0.804
## T18:GBGPL2-S
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4
## T18:GBGPL1-JEL6
## T18:GBGPL2-J
## T18:GBGPL2-S
## Tmp18:GntyB 0.364## not skewed
plot(m1_Ac)
## helpful site for evaluating which way to test you model
## http://glmm.wikidot.com/faq
## Easiest to use AND intrepret
Anova(m1_Ac, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition_percent
## Chisq Df Pr(>Chisq)
## Temp 0.0037 1 0.9513856
## Genotype 202.7352 5 < 2.2e-16 ***
## Temp:Genotype 25.2654 5 0.0001238 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## model reduction also indicates interaction
m1_Ac_r <- lmer(final_inhibition_percent ~ Temp+Genotype + (1 | Incubator/Plate), data = B1_Ac, REML = F)
## The p-value is less than 0.05 so this suggests that there is a significant interaction
anova(m1_Ac, m1_Ac_r)## Data: B1_Ac
## Models:
## m1_Ac_r: final_inhibition_percent ~ Temp + Genotype + (1 | Incubator/Plate)
## m1_Ac: final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## m1_Ac_r 10 1753.4 1785.0 -866.71 1733.4
## m1_Ac 15 1740.0 1787.3 -855.02 1710.0 23.376 5 0.000286 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## use lsmeans to look at what the interaction
p1_Ac <- lsmeans(m1_Ac, pairwise~Temp|Genotype, adjust = "tukey")
pairs(p1_Ac)## Genotype = BdBrazil-JEL649:
## contrast estimate SE df t.ratio p.value
## 12 - 18 0.9947926 22.00544 26.16 0.045 0.9643
##
## Genotype = BdGPL1-JEL404:
## contrast estimate SE df t.ratio p.value
## 12 - 18 40.8987575 20.96398 21.12 1.951 0.0645
##
## Genotype = BdGPL1-JEL647:
## contrast estimate SE df t.ratio p.value
## 12 - 18 10.4703107 20.70731 20.06 0.506 0.6186
##
## Genotype = BdGPL2-JEL423:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -34.1462570 21.07348 21.67 -1.620 0.1196
##
## Genotype = BdGPL2-SRS810:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -20.7484839 20.81783 20.56 -0.997 0.3305
##
## Genotype = Bsal:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -2.7914538 21.53128 23.83 -0.130 0.8979p1_Ac_2 <- lsmeans(m1_Ac, pairwise~Genotype|Temp, adjust = "tukey")
pairs(p1_Ac_2)## Temp = 12:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 2.500487 21.90091 37.20 0.114
## BdBrazil-JEL649 - BdGPL1-JEL647 17.979820 13.40267 170.21 1.342
## BdBrazil-JEL649 - BdGPL2-JEL423 32.467374 22.25032 39.84 1.459
## BdBrazil-JEL649 - BdGPL2-SRS810 45.914186 13.40267 170.21 3.426
## BdBrazil-JEL649 - Bsal -79.803307 22.44480 41.23 -3.556
## BdGPL1-JEL404 - BdGPL1-JEL647 15.479334 20.70731 29.52 0.748
## BdGPL1-JEL404 - BdGPL2-JEL423 29.966888 12.00306 165.51 2.497
## BdGPL1-JEL404 - BdGPL2-SRS810 43.413699 20.70731 29.52 2.097
## BdGPL1-JEL404 - Bsal -82.303794 12.36455 166.47 -6.656
## BdGPL1-JEL647 - BdGPL2-JEL423 14.487554 21.07348 31.97 0.687
## BdGPL1-JEL647 - BdGPL2-SRS810 27.934366 11.34789 165.31 2.462
## BdGPL1-JEL647 - Bsal -97.783127 21.28146 33.31 -4.595
## BdGPL2-JEL423 - BdGPL2-SRS810 13.446812 21.07348 31.97 0.638
## BdGPL2-JEL423 - Bsal -112.270681 13.01327 166.96 -8.627
## BdGPL2-SRS810 - Bsal -125.717493 21.28146 33.31 -5.907
## p.value
## 1.0000
## 0.7614
## 0.6913
## 0.0098
## 0.0116
## 0.9740
## 0.1311
## 0.3164
## <.0001
## 0.9821
## 0.1418
## 0.0008
## 0.9872
## <.0001
## <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 42.404452 21.07188 32.31 2.012
## BdBrazil-JEL649 - BdGPL1-JEL647 27.455338 11.54833 165.45 2.377
## BdBrazil-JEL649 - BdGPL2-JEL423 -2.673675 20.81783 30.64 -0.128
## BdBrazil-JEL649 - BdGPL2-SRS810 24.170910 11.75489 165.66 2.056
## BdBrazil-JEL649 - Bsal -83.589554 21.07443 32.33 -3.966
## BdGPL1-JEL404 - BdGPL1-JEL647 -14.949113 20.96398 31.58 -0.713
## BdGPL1-JEL404 - BdGPL2-JEL423 -45.078127 11.80975 165.99 -3.817
## BdGPL1-JEL404 - BdGPL2-SRS810 -18.233542 21.07443 32.33 -0.865
## BdGPL1-JEL404 - Bsal -125.994005 12.30347 167.08 -10.241
## BdGPL1-JEL647 - BdGPL2-JEL423 -30.129014 20.70731 29.91 -1.455
## BdGPL1-JEL647 - BdGPL2-SRS810 -3.284429 11.54833 165.45 -0.284
## BdGPL1-JEL647 - Bsal -111.044892 20.96398 31.58 -5.297
## BdGPL2-JEL423 - BdGPL2-SRS810 26.844585 20.81783 30.64 1.289
## BdGPL2-JEL423 - Bsal -80.915878 11.80975 165.99 -6.852
## BdGPL2-SRS810 - Bsal -107.760463 21.07188 32.31 -5.114
## p.value
## 0.3577
## 0.1703
## 1.0000
## 0.3156
## 0.0047
## 0.9789
## 0.0026
## 0.9521
## <.0001
## 0.6942
## 0.9997
## 0.0001
## 0.7883
## <.0001
## 0.0002
##
## P value adjustment: tukey method for comparing a family of 6 estimates#### IF you don't have an interaction, but have an isolate effect
#p_x <- lsmeans(m_x, "Isolate", method = "pairwise", adjust = "tukey")
#pairs(p)
d <- lsmip(m1_Ac, Genotype ~ Temp, xlab="Temperature", ylab="Percent inhibition",cex = 3, ylim = c(-90,70), lwd = 3)
BACTERIA #2, Pseudomonas sp. 1
B2_Ps1 <- Bd_scores_samples[Bd_scores_samples$Bacteria == "B2",]
## normal enough
hist(B2_Ps1$final_inhibition)
qqp((B2_Ps1$final_inhibition), "norm")
## This is what we are analyzing
ggplot(B2_Ps1, aes(x = Isolate, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## make sure temp is a factor, after subsetting
class(B2_Ps1$Temp)## [1] "factor"## MODEL
m2_Ps1 <- lmer(final_inhibition_percent ~ Temp*Genotype + (1 | Incubator/Plate), data = B2_Ps1, REML = F)
summary(m2_Ps1)## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula:
## final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Data: B2_Ps1
##
## AIC BIC logLik deviance df.resid
## 1679.0 1726.9 -824.5 1649.0 166
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5482 -0.5161 0.0030 0.5251 1.9707
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plate:Incubator (Intercept) 128.7 11.34
## Incubator (Intercept) 0.0 0.00
## Residual 467.6 21.62
## Number of obs: 181, groups: Plate:Incubator, 16; Incubator, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 122.997 9.082 13.543
## Temp18 -21.491 12.566 -1.710
## GenotypeBdGPL1-JEL404 -86.267 11.995 -7.192
## GenotypeBdGPL1-JEL647 -81.155 8.919 -9.099
## GenotypeBdGPL2-JEL423 -112.445 11.995 -9.374
## GenotypeBdGPL2-SRS810 -93.602 8.919 -10.495
## GenotypeBsal -58.323 11.995 -4.862
## Temp18:GenotypeBdGPL1-JEL404 -14.958 16.754 -0.893
## Temp18:GenotypeBdGPL1-JEL647 -15.237 12.329 -1.236
## Temp18:GenotypeBdGPL2-JEL423 35.444 16.754 2.116
## Temp18:GenotypeBdGPL2-SRS810 -58.850 12.329 -4.773
## Temp18:GenotypeBsal 4.305 16.754 0.257
##
## Correlation of Fixed Effects:
## (Intr) Temp18 GBGPL1-JEL4 GBGPL1-JEL6 GBGPL2-J GBGPL2-S
## Temp18 -0.723
## GBGPL1-JEL4 -0.757 0.547
## GBGPL1-JEL6 -0.621 0.449 0.470
## GBGPL2-JEL4 -0.757 0.547 0.797 0.470
## GBGPL2-SRS8 -0.621 0.449 0.470 0.633 0.470
## GenotypeBsl -0.757 0.547 0.797 0.470 0.797 0.470
## T18:GBGPL1-JEL4 0.542 -0.750 -0.716 -0.337 -0.571 -0.337
## T18:GBGPL1-JEL6 0.449 -0.604 -0.340 -0.723 -0.340 -0.458
## T18:GBGPL2-J 0.542 -0.750 -0.571 -0.337 -0.716 -0.337
## T18:GBGPL2-S 0.449 -0.604 -0.340 -0.458 -0.340 -0.723
## Tmp18:GntyB 0.542 -0.750 -0.571 -0.337 -0.571 -0.337
## GntypB T18:GBGPL1-JEL4 T18:GBGPL1-JEL6 T18:GBGPL2-J
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4 -0.571
## T18:GBGPL1-JEL6 -0.340 0.453
## T18:GBGPL2-J -0.571 0.792 0.453
## T18:GBGPL2-S -0.340 0.453 0.615 0.453
## Tmp18:GntyB -0.716 0.792 0.453 0.792
## T18:GBGPL2-S
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4
## T18:GBGPL1-JEL6
## T18:GBGPL2-J
## T18:GBGPL2-S
## Tmp18:GntyB 0.453# want random distribution of data pts
plot(m2_Ps1)
Anova(m2_Ps1, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition_percent
## Chisq Df Pr(>Chisq)
## Temp 22.813 1 1.785e-06 ***
## Genotype 479.078 5 < 2.2e-16 ***
## Temp:Genotype 57.389 5 4.203e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
m2_Ps1_r <- lmer(final_inhibition_percent ~ Temp+Genotype + (1 | Incubator/Plate), data = B2_Ps1, REML = F)
anova(m2_Ps1, m2_Ps1_r)## Data: B2_Ps1
## Models:
## m2_Ps1_r: final_inhibition_percent ~ Temp + Genotype + (1 | Incubator/Plate)
## m2_Ps1: final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## m2_Ps1_r 10 1718.5 1750.5 -849.25 1698.5
## m2_Ps1 15 1679.0 1726.9 -824.48 1649.0 49.528 5 1.731e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#p-value less than 0.05, suggesting there is a significant interaction
## use lsmeans to look at what the interaction
p2_Ps1 <- lsmeans(m2_Ps1, pairwise~Temp|Genotype, adjust = "tukey")
pairs(p2_Ps1)## Genotype = BdBrazil-JEL649:
## contrast estimate SE df t.ratio p.value
## 12 - 18 21.49059 13.88681 38.79 1.548 0.1299
##
## Genotype = BdGPL1-JEL404:
## contrast estimate SE df t.ratio p.value
## 12 - 18 36.44900 12.43315 25.36 2.932 0.0070
##
## Genotype = BdGPL1-JEL647:
## contrast estimate SE df t.ratio p.value
## 12 - 18 36.72744 12.43315 24.66 2.954 0.0068
##
## Genotype = BdGPL2-JEL423:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -13.95330 12.43315 25.36 -1.122 0.2723
##
## Genotype = BdGPL2-SRS810:
## contrast estimate SE df t.ratio p.value
## 12 - 18 80.34032 12.43315 24.66 6.462 <.0001
##
## Genotype = Bsal:
## contrast estimate SE df t.ratio p.value
## 12 - 18 17.18586 12.43315 25.36 1.382 0.1789p2_Ps1_2 <- lsmeans(m2_Ps1, pairwise~Genotype|Temp, adjust = "tukey")
pairs(p2_Ps1_2)## Temp = 12:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 86.266568 13.340777 48.55 6.466
## BdBrazil-JEL649 - BdGPL1-JEL647 81.154537 9.206304 179.16 8.815
## BdBrazil-JEL649 - BdGPL2-JEL423 112.444561 13.340777 48.55 8.429
## BdBrazil-JEL649 - BdGPL2-SRS810 93.602177 9.206304 179.16 10.167
## BdBrazil-JEL649 - Bsal 58.322635 13.340777 48.55 4.372
## BdGPL1-JEL404 - BdGPL1-JEL647 -5.112031 12.433151 37.06 -0.411
## BdGPL1-JEL404 - BdGPL2-JEL423 26.177993 7.833449 173.02 3.342
## BdGPL1-JEL404 - BdGPL2-SRS810 7.335609 12.433151 37.06 0.590
## BdGPL1-JEL404 - Bsal -27.943933 7.833449 173.02 -3.567
## BdGPL1-JEL647 - BdGPL2-JEL423 31.290024 12.433151 37.06 2.517
## BdGPL1-JEL647 - BdGPL2-SRS810 12.447640 7.833449 173.02 1.589
## BdGPL1-JEL647 - Bsal -22.831902 12.433151 37.06 -1.836
## BdGPL2-JEL423 - BdGPL2-SRS810 -18.842385 12.433151 37.06 -1.515
## BdGPL2-JEL423 - Bsal -54.121926 7.833449 173.02 -6.909
## BdGPL2-SRS810 - Bsal -35.279541 12.433151 37.06 -2.838
## p.value
## <.0001
## <.0001
## <.0001
## <.0001
## 0.0009
## 0.9984
## 0.0129
## 0.9911
## 0.0061
## 0.1454
## 0.6069
## 0.4560
## 0.6568
## <.0001
## 0.0734
##
## Temp = 18:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 101.224978 13.017315 45.21 7.776
## BdBrazil-JEL649 - BdGPL1-JEL647 96.391392 8.730991 174.17 11.040
## BdBrazil-JEL649 - BdGPL2-JEL423 77.000668 13.017315 45.21 5.915
## BdBrazil-JEL649 - BdGPL2-SRS810 152.451909 8.730991 174.17 17.461
## BdBrazil-JEL649 - Bsal 54.017907 13.017315 45.21 4.150
## BdGPL1-JEL404 - BdGPL1-JEL647 -4.833586 12.433151 37.21 -0.389
## BdGPL1-JEL404 - BdGPL2-JEL423 -24.224309 7.833449 173.02 -3.092
## BdGPL1-JEL404 - BdGPL2-SRS810 51.226932 12.433151 37.21 4.120
## BdGPL1-JEL404 - Bsal -47.207071 7.833449 173.02 -6.026
## BdGPL1-JEL647 - BdGPL2-JEL423 -19.390724 12.433151 37.21 -1.560
## BdGPL1-JEL647 - BdGPL2-SRS810 56.060517 7.833449 173.02 7.157
## BdGPL1-JEL647 - Bsal -42.373485 12.433151 37.21 -3.408
## BdGPL2-JEL423 - BdGPL2-SRS810 75.451241 12.433151 37.21 6.069
## BdGPL2-JEL423 - Bsal -22.982761 7.833449 173.02 -2.934
## BdGPL2-SRS810 - Bsal -98.434002 12.433151 37.21 -7.917
## p.value
## <.0001
## <.0001
## <.0001
## <.0001
## 0.0019
## 0.9988
## 0.0276
## 0.0026
## <.0001
## 0.6292
## <.0001
## 0.0183
## <.0001
## 0.0434
## <.0001
##
## P value adjustment: tukey method for comparing a family of 6 estimateslsmip(m2_Ps1, Genotype ~ Temp, xlab="Temperature", ylab="Percent inhibition",cex = 3, ylim = c(-120, 130), lwd = 3)
BACTERIA #4, Chryseobacterium sp.
B4_Ch <- Bd_scores_samples[Bd_scores_samples$Bacteria == "B4",]
## normal
hist(B4_Ch$final_inhibition)
qqp((B4_Ch$final_inhibition), "norm")
## This is what we are analyzing
ggplot(B4_Ch, aes(x = Isolate, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## make sure temp is a factor, after subsetting
class(B4_Ch$Temp)## [1] "factor"## MODEL
m4_Ch <- lmer(final_inhibition_percent ~ Temp*Genotype + (1 | Incubator/Plate), data = B4_Ch, REML = F)
summary(m4_Ch)## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula:
## final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Data: B4_Ch
##
## AIC BIC logLik deviance df.resid
## 1784.2 1832.7 -877.1 1754.2 172
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5174 -0.5491 0.0328 0.5997 2.3110
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plate:Incubator (Intercept) 452.8 21.28
## Incubator (Intercept) 0.0 0.00
## Residual 568.8 23.85
## Number of obs: 187, groups: Plate:Incubator, 16; Incubator, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 48.030 13.050 3.681
## Temp18 -5.306 17.862 -0.297
## GenotypeBdGPL1-JEL404 -43.875 17.862 -2.456
## GenotypeBdGPL1-JEL647 -80.530 9.625 -8.367
## GenotypeBdGPL2-JEL423 -135.351 17.862 -7.578
## GenotypeBdGPL2-SRS810 -132.131 9.625 -13.728
## GenotypeBsal 23.269 17.862 1.303
## Temp18:GenotypeBdGPL1-JEL404 -54.916 24.831 -2.212
## Temp18:GenotypeBdGPL1-JEL647 42.712 12.796 3.338
## Temp18:GenotypeBdGPL2-JEL423 56.828 24.831 2.289
## Temp18:GenotypeBdGPL2-SRS810 -19.728 12.796 -1.542
## Temp18:GenotypeBsal -22.792 24.831 -0.918
##
## Correlation of Fixed Effects:
## (Intr) Temp18 GBGPL1-JEL4 GBGPL1-JEL6 GBGPL2-J GBGPL2-S
## Temp18 -0.731
## GBGPL1-JEL4 -0.731 0.534
## GBGPL1-JEL6 -0.454 0.332 0.332
## GBGPL2-JEL4 -0.731 0.534 0.889 0.332
## GBGPL2-SRS8 -0.454 0.332 0.332 0.616 0.332
## GenotypeBsl -0.731 0.534 0.889 0.332 0.889 0.332
## T18:GBGPL1-JEL4 0.526 -0.719 -0.719 -0.239 -0.639 -0.239
## T18:GBGPL1-JEL6 0.342 -0.405 -0.250 -0.752 -0.250 -0.463
## T18:GBGPL2-J 0.526 -0.719 -0.639 -0.239 -0.719 -0.239
## T18:GBGPL2-S 0.342 -0.405 -0.250 -0.463 -0.250 -0.752
## Tmp18:GntyB 0.526 -0.719 -0.639 -0.239 -0.639 -0.239
## GntypB T18:GBGPL1-JEL4 T18:GBGPL1-JEL6 T18:GBGPL2-J
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4 -0.639
## T18:GBGPL1-JEL6 -0.250 0.292
## T18:GBGPL2-J -0.639 0.885 0.292
## T18:GBGPL2-S -0.250 0.292 0.566 0.292
## Tmp18:GntyB -0.719 0.885 0.292 0.885
## T18:GBGPL2-S
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4
## T18:GBGPL1-JEL6
## T18:GBGPL2-J
## T18:GBGPL2-S
## Tmp18:GntyB 0.292plot(m4_Ch)
Anova(m4_Ch, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition_percent
## Chisq Df Pr(>Chisq)
## Temp 0.1733 1 0.6772
## Genotype 932.8795 5 <2e-16 ***
## Temp:Genotype 122.0062 5 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
m4_Ch_r <- lmer(final_inhibition_percent ~ Temp+Genotype + (1 | Incubator/Plate), data = B4_Ch, REML = F)
anova(m4_Ch, m4_Ch_r)## Data: B4_Ch
## Models:
## m4_Ch_r: final_inhibition_percent ~ Temp + Genotype + (1 | Incubator/Plate)
## m4_Ch: final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## m4_Ch_r 10 1866.4 1898.7 -923.19 1846.4
## m4_Ch 15 1784.2 1832.7 -877.10 1754.2 92.194 5 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## use lsmeans to look at what the interaction
p4_Ch <- lsmeans(m4_Ch, pairwise~Temp|Genotype, adjust = "tukey")
pairs(p4_Ch)## Genotype = BdBrazil-JEL649:
## contrast estimate SE df t.ratio p.value
## 12 - 18 5.306101 20.17808 19.14 0.263 0.7954
##
## Genotype = BdGPL1-JEL404:
## contrast estimate SE df t.ratio p.value
## 12 - 18 60.222052 19.59554 16.82 3.073 0.0070
##
## Genotype = BdGPL1-JEL647:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -37.405649 19.59554 16.72 -1.909 0.0736
##
## Genotype = BdGPL2-JEL423:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -51.521665 19.59554 16.82 -2.629 0.0177
##
## Genotype = BdGPL2-SRS810:
## contrast estimate SE df t.ratio p.value
## 12 - 18 25.034062 19.59554 16.72 1.278 0.2189
##
## Genotype = Bsal:
## contrast estimate SE df t.ratio p.value
## 12 - 18 28.098247 19.59554 16.82 1.434 0.1699p4_Ch_2 <- lsmeans(m4_Ch, pairwise~Genotype|Temp, adjust = "tukey")
pairs(p4_Ch_2)## Temp = 12:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 43.8745490 20.178081 28.44 2.174
## BdBrazil-JEL649 - BdGPL1-JEL647 80.5301049 9.886986 182.40 8.145
## BdBrazil-JEL649 - BdGPL2-JEL423 135.3510636 20.178081 28.44 6.708
## BdBrazil-JEL649 - BdGPL2-SRS810 132.1306039 9.886986 182.40 13.364
## BdBrazil-JEL649 - Bsal -23.2685031 20.178081 28.44 -1.153
## BdGPL1-JEL404 - BdGPL1-JEL647 36.6555559 19.595536 24.94 1.871
## BdGPL1-JEL404 - BdGPL2-JEL423 91.4765146 8.636121 179.35 10.592
## BdGPL1-JEL404 - BdGPL2-SRS810 88.2560549 19.595536 24.94 4.504
## BdGPL1-JEL404 - Bsal -67.1430521 8.636121 179.35 -7.775
## BdGPL1-JEL647 - BdGPL2-JEL423 54.8209587 19.595536 24.94 2.798
## BdGPL1-JEL647 - BdGPL2-SRS810 51.6004990 8.636121 179.35 5.975
## BdGPL1-JEL647 - Bsal -103.7986080 19.595536 24.94 -5.297
## BdGPL2-JEL423 - BdGPL2-SRS810 -3.2204598 19.595536 24.94 -0.164
## BdGPL2-JEL423 - Bsal -158.6195667 8.636121 179.35 -18.367
## BdGPL2-SRS810 - Bsal -155.3991069 19.595536 24.94 -7.930
## p.value
## 0.2804
## <.0001
## <.0001
## <.0001
## 0.8547
## 0.4421
## <.0001
## 0.0017
## <.0001
## 0.0911
## <.0001
## 0.0002
## 1.0000
## <.0001
## <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 98.7904995 19.595536 25.08 5.041
## BdBrazil-JEL649 - BdGPL1-JEL647 37.8183541 8.636121 179.35 4.379
## BdBrazil-JEL649 - BdGPL2-JEL423 78.5232973 19.595536 25.08 4.007
## BdBrazil-JEL649 - BdGPL2-SRS810 151.8585649 8.636121 179.35 17.584
## BdBrazil-JEL649 - Bsal -0.4763574 19.595536 25.08 -0.024
## BdGPL1-JEL404 - BdGPL1-JEL647 -60.9721454 19.595536 25.08 -3.112
## BdGPL1-JEL404 - BdGPL2-JEL423 -20.2672022 8.636121 179.35 -2.347
## BdGPL1-JEL404 - BdGPL2-SRS810 53.0680654 19.595536 25.08 2.708
## BdGPL1-JEL404 - Bsal -99.2668569 8.636121 179.35 -11.494
## BdGPL1-JEL647 - BdGPL2-JEL423 40.7049432 19.595536 25.08 2.077
## BdGPL1-JEL647 - BdGPL2-SRS810 114.0402108 8.636121 179.35 13.205
## BdGPL1-JEL647 - Bsal -38.2947115 19.595536 25.08 -1.954
## BdGPL2-JEL423 - BdGPL2-SRS810 73.3352676 19.595536 25.08 3.742
## BdGPL2-JEL423 - Bsal -78.9996546 8.636121 179.35 -9.148
## BdGPL2-SRS810 - Bsal -152.3349222 19.595536 25.08 -7.774
## p.value
## 0.0004
## 0.0003
## 0.0057
## <.0001
## 1.0000
## 0.0468
## 0.1811
## 0.1088
## <.0001
## 0.3304
## <.0001
## 0.3948
## 0.0109
## <.0001
## <.0001
##
## P value adjustment: tukey method for comparing a family of 6 estimateslsmip(m4_Ch, Genotype ~ Temp, xlab="Temperature", ylab="percent inhibition",cex = 3, ylim = c(-120,130), lwd = 3)
BACTERIA #3, Pseudomonas sp. 2
B3_Ps2 <- Bd_scores_samples[Bd_scores_samples$Bacteria == "B3",]
## least normal, but all other bacteria were fine as is, so to keep consistent use data as is
hist((B3_Ps2$final_inhibition))
qqp((B3_Ps2$final_inhibition), "norm")
## This is what we are analyzing
ggplot(B3_Ps2, aes(x = Isolate, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## make sure temp is a factor, after subsetting
class(B3_Ps2$Temp)## [1] "factor"## MODEL
m3_Ps2 <- lmer(final_inhibition_percent ~ Temp*Genotype + (1 | Incubator/Plate), data = B3_Ps2, REML = F)
summary(m3_Ps2)## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula:
## final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Data: B3_Ps2
##
## AIC BIC logLik deviance df.resid
## 1217.9 1265.2 -594.0 1187.9 158
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.6487 -0.3921 0.0504 0.4697 3.1066
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plate:Incubator (Intercept) 29.75 5.455
## Incubator (Intercept) 0.00 0.000
## Residual 46.44 6.815
## Number of obs: 173, groups: Plate:Incubator, 16; Incubator, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 58.137 3.428 16.958
## Temp18 34.253 4.790 7.151
## GenotypeBdGPL1-JEL404 39.589 4.778 8.286
## GenotypeBdGPL1-JEL647 4.193 2.686 1.561
## GenotypeBdGPL2-JEL423 28.472 4.722 6.029
## GenotypeBdGPL2-SRS810 38.027 2.744 13.859
## GenotypeBsal 33.843 4.722 7.167
## Temp18:GenotypeBdGPL1-JEL404 -39.068 6.660 -5.866
## Temp18:GenotypeBdGPL1-JEL647 -13.316 3.802 -3.502
## Temp18:GenotypeBdGPL2-JEL423 -26.772 6.621 -4.044
## Temp18:GenotypeBdGPL2-SRS810 -34.259 3.882 -8.825
## Temp18:GenotypeBsal -33.917 6.621 -5.123
##
## Correlation of Fixed Effects:
## (Intr) Temp18 GBGPL1-JEL4 GBGPL1-JEL6 GBGPL2-J GBGPL2-S
## Temp18 -0.716
## GBGPL1-JEL4 -0.718 0.514
## GBGPL1-JEL6 -0.468 0.335 0.336
## GBGPL2-JEL4 -0.726 0.520 0.850 0.340
## GBGPL2-SRS8 -0.452 0.324 0.324 0.577 0.328
## GenotypeBsl -0.726 0.520 0.850 0.340 0.860 0.328
## T18:GBGPL1-JEL4 0.515 -0.719 -0.717 -0.241 -0.610 -0.233
## T18:GBGPL1-JEL6 0.331 -0.446 -0.238 -0.707 -0.240 -0.408
## T18:GBGPL2-J 0.518 -0.724 -0.606 -0.243 -0.713 -0.234
## T18:GBGPL2-S 0.320 -0.434 -0.229 -0.408 -0.232 -0.707
## Tmp18:GntyB 0.518 -0.724 -0.606 -0.243 -0.614 -0.234
## GntypB T18:GBGPL1-JEL4 T18:GBGPL1-JEL6 T18:GBGPL2-J
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4 -0.610
## T18:GBGPL1-JEL6 -0.240 0.320
## T18:GBGPL2-J -0.614 0.857 0.322
## T18:GBGPL2-S -0.232 0.312 0.555 0.314
## Tmp18:GntyB -0.713 0.857 0.322 0.863
## T18:GBGPL2-S
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4
## T18:GBGPL1-JEL6
## T18:GBGPL2-J
## T18:GBGPL2-S
## Tmp18:GntyB 0.314plot(m3_Ps2)
Anova(m3_Ps2, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition_percent
## Chisq Df Pr(>Chisq)
## Temp 10.436 1 0.001236 **
## Genotype 219.211 5 < 2.2e-16 ***
## Temp:Genotype 101.067 5 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
m3_Ps2_r <- lmer(final_inhibition_percent ~ Temp+Genotype + (1 | Incubator/Plate), data = B3_Ps2, REML = F)
anova(m3_Ps2, m3_Ps2_r)## Data: B3_Ps2
## Models:
## m3_Ps2_r: final_inhibition_percent ~ Temp + Genotype + (1 | Incubator/Plate)
## m3_Ps2: final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## m3_Ps2_r 10 1287.5 1319.1 -633.76 1267.5
## m3_Ps2 15 1217.9 1265.2 -593.97 1187.9 79.586 5 1.025e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## use lsmeans to look at what the interaction
p3_Ps2 <- lsmeans(m3_Ps2, pairwise~Temp|Genotype, adjust = "tukey")
pairs(p3_Ps2)## Genotype = BdBrazil-JEL649:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -34.253260568 5.387546 21.19 -6.358 <.0001
##
## Genotype = BdGPL1-JEL404:
## contrast estimate SE df t.ratio p.value
## 12 - 18 4.814246744 5.227136 18.89 0.921 0.3687
##
## Genotype = BdGPL1-JEL647:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -20.936959714 5.202597 18.20 -4.024 0.0008
##
## Genotype = BdGPL2-JEL423:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -7.481650655 5.173162 17.99 -1.446 0.1653
##
## Genotype = BdGPL2-SRS810:
## contrast estimate SE df t.ratio p.value
## 12 - 18 0.006083404 5.272429 19.39 0.001 0.9991
##
## Genotype = Bsal:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -0.336358372 5.173158 17.99 -0.065 0.9489p3_Ps2_2 <- lsmeans(m3_Ps2, pairwise~Genotype|Temp, adjust = "tukey")
pairs(p3_Ps2_2)## Temp = 12:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 -39.58938906 5.373739 31.09 -7.367
## BdBrazil-JEL649 - BdGPL1-JEL647 -4.19251339 2.771020 169.61 -1.513
## BdBrazil-JEL649 - BdGPL2-JEL423 -28.47171453 5.322487 29.78 -5.349
## BdBrazil-JEL649 - BdGPL2-SRS810 -38.02720338 2.827075 168.52 -13.451
## BdBrazil-JEL649 - Bsal -33.84278309 5.321423 29.76 -6.360
## BdGPL1-JEL404 - BdGPL1-JEL647 35.39687567 5.227136 27.69 6.772
## BdGPL1-JEL404 - BdGPL2-JEL423 11.11767453 2.674105 166.26 4.158
## BdGPL1-JEL404 - BdGPL2-SRS810 1.56218567 5.271022 28.79 0.296
## BdGPL1-JEL404 - Bsal 5.74660597 2.673729 166.26 2.149
## BdGPL1-JEL647 - BdGPL2-JEL423 -24.27920114 5.173162 26.40 -4.693
## BdGPL1-JEL647 - BdGPL2-SRS810 -33.83468999 2.564791 165.66 -13.192
## BdGPL1-JEL647 - Bsal -29.65026969 5.173158 26.40 -5.732
## BdGPL2-JEL423 - BdGPL2-SRS810 -9.55548885 5.217560 27.48 -1.831
## BdGPL2-JEL423 - Bsal -5.37106855 2.562023 165.72 -2.096
## BdGPL2-SRS810 - Bsal 4.18442030 5.217506 27.48 0.802
## p.value
## <.0001
## 0.6564
## 0.0001
## <.0001
## <.0001
## <.0001
## 0.0007
## 0.9997
## 0.2675
## 0.0009
## <.0001
## 0.0001
## 0.4633
## 0.2943
## 0.9647
##
## Temp = 18:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 -0.52188175 5.240629 28.03 -0.100
## BdBrazil-JEL649 - BdGPL1-JEL647 9.12378746 2.768863 168.43 3.295
## BdBrazil-JEL649 - BdGPL2-JEL423 -1.70010462 5.240629 28.03 -0.324
## BdBrazil-JEL649 - BdGPL2-SRS810 -3.76785941 2.826839 168.64 -1.333
## BdBrazil-JEL649 - Bsal 0.07411911 5.240629 28.03 0.014
## BdGPL1-JEL404 - BdGPL1-JEL647 9.64566921 5.202597 27.18 1.854
## BdGPL1-JEL404 - BdGPL2-JEL423 -1.17822287 2.473276 165.47 -0.476
## BdGPL1-JEL404 - BdGPL2-SRS810 -3.24597767 5.228523 27.82 -0.621
## BdGPL1-JEL404 - Bsal 0.59600086 2.473276 165.47 0.241
## BdGPL1-JEL647 - BdGPL2-JEL423 -10.82389208 5.202597 27.18 -2.080
## BdGPL1-JEL647 - BdGPL2-SRS810 -12.89164688 2.697811 165.58 -4.779
## BdGPL1-JEL647 - Bsal -9.04966835 5.202597 27.18 -1.739
## BdGPL2-JEL423 - BdGPL2-SRS810 -2.06775479 5.228523 27.82 -0.395
## BdGPL2-JEL423 - Bsal 1.77422373 2.473276 165.47 0.717
## BdGPL2-SRS810 - Bsal 3.84197852 5.228523 27.82 0.735
## p.value
## 1.0000
## 0.0150
## 0.9995
## 0.7663
## 1.0000
## 0.4502
## 0.9969
## 0.9885
## 0.9999
## 0.3266
## 0.0001
## 0.5191
## 0.9986
## 0.9796
## 0.9758
##
## P value adjustment: tukey method for comparing a family of 6 estimateslsmip(m3_Ps2, Genotype ~ Temp, xlab="Temperature", ylab="Percent inhibition", cex = 3, ylim = c(55,105), lwd = 3)
BACTERIA #5, Stenotrophomonas sp.
B5_St <- Bd_scores_samples[Bd_scores_samples$Bacteria == "B5",]
## normal enough
hist(B5_St$final_inhibition)
qqp((B5_St$final_inhibition), "norm")
## This is what we are analyzing
ggplot(B5_St, aes(x = Isolate, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## make sure temp is a factor, after subsetting
class(B5_St$Temp)## [1] "factor"## MODEL
m5_St <- lmer(final_inhibition_percent ~ Temp*Genotype + (1 | Incubator/Plate), data = B5_St, REML = F)
summary(m5_St)## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula:
## final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Data: B5_St
##
## AIC BIC logLik deviance df.resid
## 1169.3 1214.8 -569.7 1139.3 138
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.8560 -0.4115 0.0792 0.4232 2.4415
##
## Random effects:
## Groups Name Variance Std.Dev.
## Plate:Incubator (Intercept) 44.84 6.696
## Incubator (Intercept) 0.00 0.000
## Residual 83.18 9.120
## Number of obs: 153, groups: Plate:Incubator, 16; Incubator, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 73.2163 4.6665 15.690
## Temp18 27.0734 6.2092 4.360
## GenotypeBdGPL1-JEL404 27.0347 6.2093 4.354
## GenotypeBdGPL1-JEL647 0.5087 4.3324 0.117
## GenotypeBdGPL2-JEL423 15.3419 6.2848 2.441
## GenotypeBdGPL2-SRS810 24.6954 4.0450 6.105
## GenotypeBsal 19.1339 6.1793 3.096
## Temp18:GenotypeBdGPL1-JEL404 -30.7496 8.6182 -3.568
## Temp18:GenotypeBdGPL1-JEL647 -4.4296 5.6092 -0.790
## Temp18:GenotypeBdGPL2-JEL423 -19.7257 8.5713 -2.301
## Temp18:GenotypeBdGPL2-SRS810 -23.6801 5.3332 -4.440
## Temp18:GenotypeBsal -27.1122 8.5628 -3.166
##
## Correlation of Fixed Effects:
## (Intr) Temp18 GBGPL1-JEL4 GBGPL1-JEL6 GBGPL2-J GBGPL2-S
## Temp18 -0.752
## GBGPL1-JEL4 -0.752 0.565
## GBGPL1-JEL6 -0.487 0.366 0.366
## GBGPL2-JEL4 -0.743 0.558 0.846 0.361
## GBGPL2-SRS8 -0.555 0.417 0.417 0.560 0.412
## GenotypeBsl -0.755 0.568 0.860 0.368 0.849 0.419
## T18:GBGPL1-JEL4 0.541 -0.720 -0.720 -0.264 -0.610 -0.300
## T18:GBGPL1-JEL6 0.376 -0.444 -0.283 -0.772 -0.279 -0.432
## T18:GBGPL2-J 0.544 -0.724 -0.621 -0.265 -0.733 -0.302
## T18:GBGPL2-S 0.421 -0.484 -0.316 -0.425 -0.312 -0.758
## Tmp18:GntyB 0.545 -0.725 -0.620 -0.265 -0.613 -0.302
## GntypB T18:GBGPL1-JEL4 T18:GBGPL1-JEL6 T18:GBGPL2-J
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4 -0.619
## T18:GBGPL1-JEL6 -0.284 0.320
## T18:GBGPL2-J -0.623 0.827 0.321
## T18:GBGPL2-S -0.318 0.349 0.517 0.351
## Tmp18:GntyB -0.722 0.829 0.322 0.832
## T18:GBGPL2-S
## Temp18
## GBGPL1-JEL4
## GBGPL1-JEL6
## GBGPL2-JEL4
## GBGPL2-SRS8
## GenotypeBsl
## T18:GBGPL1-JEL4
## T18:GBGPL1-JEL6
## T18:GBGPL2-J
## T18:GBGPL2-S
## Tmp18:GntyB 0.351plot(m5_St)
Anova(m5_St, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition_percent
## Chisq Df Pr(>Chisq)
## Temp 5.5081 1 0.01893 *
## Genotype 40.5992 5 1.131e-07 ***
## Temp:Genotype 31.8633 5 6.323e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
m5_St_r <- lmer(final_inhibition_percent ~ Temp+Genotype + (1 | Incubator/Plate), data = B5_St, REML = F)
anova(m5_St, m5_St_r)## Data: B5_St
## Models:
## m5_St_r: final_inhibition_percent ~ Temp + Genotype + (1 | Incubator/Plate)
## m5_St: final_inhibition_percent ~ Temp * Genotype + (1 | Incubator/Plate)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## m5_St_r 10 1188.0 1218.3 -584.02 1168.0
## m5_St 15 1169.3 1214.8 -569.66 1139.3 28.708 5 2.646e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## use lsmeans to look at what the interaction
p5_St <- lsmeans(m5_St, pairwise~Temp|Genotype, adjust = "tukey")
pairs(p5_St)## Genotype = BdBrazil-JEL649:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -27.07341681 6.980869 24.59 -3.878 0.0007
##
## Genotype = BdGPL1-JEL404:
## contrast estimate SE df t.ratio p.value
## 12 - 18 3.67620852 6.741920 21.83 0.545 0.5911
##
## Genotype = BdGPL1-JEL647:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -22.64384258 7.009935 25.41 -3.230 0.0034
##
## Genotype = BdGPL2-JEL423:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -7.34768898 6.677615 20.88 -1.100 0.2837
##
## Genotype = BdGPL2-SRS810:
## contrast estimate SE df t.ratio p.value
## 12 - 18 -3.39331422 6.680693 20.45 -0.508 0.6169
##
## Genotype = Bsal:
## contrast estimate SE df t.ratio p.value
## 12 - 18 0.03882977 6.666744 20.72 0.006 0.9954p5_St_2 <- lsmeans(m5_St, pairwise~Genotype|Temp, adjust = "tukey")
pairs(p5_St_2)## Temp = 12:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 -27.0346801 6.980091 35.20 -3.873
## BdBrazil-JEL649 - BdGPL1-JEL647 -0.5087322 4.470214 146.91 -0.114
## BdBrazil-JEL649 - BdGPL2-JEL423 -15.3419148 7.051956 36.77 -2.176
## BdBrazil-JEL649 - BdGPL2-SRS810 -24.6954401 4.194408 151.92 -5.888
## BdBrazil-JEL649 - Bsal -19.1338961 6.952405 34.59 -2.752
## BdGPL1-JEL404 - BdGPL1-JEL647 26.5259480 6.896368 34.22 3.846
## BdGPL1-JEL404 - BdGPL2-JEL423 11.6927653 3.568094 145.59 3.277
## BdGPL1-JEL404 - BdGPL2-SRS810 2.3392401 6.602318 28.29 0.354
## BdGPL1-JEL404 - Bsal 7.9007840 3.381520 145.53 2.336
## BdGPL1-JEL647 - BdGPL2-JEL423 -14.8331826 6.969172 35.80 -2.128
## BdGPL1-JEL647 - BdGPL2-SRS810 -24.1867079 4.065040 147.30 -5.950
## BdGPL1-JEL647 - Bsal -18.6251639 6.868050 33.61 -2.712
## BdGPL2-JEL423 - BdGPL2-SRS810 -9.3535253 6.678377 29.79 -1.401
## BdGPL2-JEL423 - Bsal -3.7919813 3.527918 145.93 -1.075
## BdGPL2-SRS810 - Bsal 5.5615440 6.572542 27.70 0.846
## p.value
## 0.0055
## 1.0000
## 0.2734
## <.0001
## 0.0904
## 0.0061
## 0.0162
## 0.9992
## 0.1864
## 0.2963
## <.0001
## 0.0992
## 0.7264
## 0.8906
## 0.9559
##
## Temp = 18:
## contrast estimate SE df t.ratio
## BdBrazil-JEL649 - BdGPL1-JEL404 3.7149452 6.741702 31.36 0.551
## BdBrazil-JEL649 - BdGPL1-JEL647 3.9208421 3.674307 146.41 1.067
## BdBrazil-JEL649 - BdGPL2-JEL423 4.3838130 6.601360 28.52 0.664
## BdBrazil-JEL649 - BdGPL2-SRS810 -1.0153375 3.584545 146.12 -0.283
## BdBrazil-JEL649 - Bsal 7.9783505 6.696433 30.38 1.191
## BdGPL1-JEL404 - BdGPL1-JEL647 0.2058969 6.857425 33.82 0.030
## BdGPL1-JEL404 - BdGPL2-JEL423 0.6688678 3.795111 145.63 0.176
## BdGPL1-JEL404 - BdGPL2-SRS810 -4.7302827 6.817752 32.95 -0.694
## BdGPL1-JEL404 - Bsal 4.2634053 3.927533 145.52 1.086
## BdGPL1-JEL647 - BdGPL2-JEL423 0.4629710 6.719404 30.90 0.069
## BdGPL1-JEL647 - BdGPL2-SRS810 -4.9361796 3.781514 145.87 -1.305
## BdGPL1-JEL647 - Bsal 4.0575084 6.812825 32.81 0.596
## BdGPL2-JEL423 - BdGPL2-SRS810 -5.3991505 6.679437 30.07 -0.808
## BdGPL2-JEL423 - Bsal 3.5945375 3.708102 145.74 0.969
## BdGPL2-SRS810 - Bsal 8.9936880 6.773431 31.96 1.328
## p.value
## 0.9934
## 0.8936
## 0.9845
## 0.9997
## 0.8375
## 1.0000
## 1.0000
## 0.9814
## 0.8865
## 1.0000
## 0.7817
## 0.9906
## 0.9638
## 0.9269
## 0.7677
##
## P value adjustment: tukey method for comparing a family of 6 estimateslsmip(m5_St, Genotype ~ Temp, xlab="Temperature", ylab="Percent inhibition", cex = 3, lwd = 3, ylim = c(55,105))
CLUSTER ANALYSIS
### generate table of average inhibition values
## Table of values
Inh_avg <- ddply(Bd_scores_samples, .(Bacteria, Temp, Genotype), summarize, N=length(Bacteria), mean=mean(final_inhibition_percent))
Inh_avg <- Inh_avg[,c(3,5)]
## Probably better way to do this, but this works
df_Inh <- cbind(Inh_avg[1:6,1:2], Inh_avg[7:12,2], Inh_avg[13:18,2], Inh_avg[19:24,2], Inh_avg[25:30,2], Inh_avg[31:36,2], Inh_avg[37:42,2], Inh_avg[43:48,2], Inh_avg[49:54,2], Inh_avg[55:60,2])
rownames(df_Inh) <- df_Inh[,1]
df_Inh <- df_Inh[,-1]
dm_Inh <- as.matrix(df_Inh)
## Ward Hierarchical Clustering
#### make a distance matrix using euclidean distances
dist <- dist(dm_Inh, method = "euclidean")
clust <- hclust(dist, method="ward.D2")
plot(clust, hang = -1, las = 1, cex =0.8)
groups_Inh <- cutree(clust , k=2)
##### cut tree into 2 clusters
#### draw dendogram with red borders around the 2 clusters
rect.hclust(clust, k=2, border="red")
## Use pvclust to get p-values
dm_Inh_t <- t(dm_Inh)
fit_Inh <- pvclust(dm_Inh_t, method.hclust="ward.D2",method.dist="euclidean")## Bootstrap (r = 0.5)... Done.
## Bootstrap (r = 0.6)... Done.
## Bootstrap (r = 0.7)... Done.
## Bootstrap (r = 0.8)... Done.
## Bootstrap (r = 0.9)... Done.
## Bootstrap (r = 1.0)... Done.
## Bootstrap (r = 1.1)... Done.
## Bootstrap (r = 1.2)... Done.
## Bootstrap (r = 1.3)... Done.
## Bootstrap (r = 1.4)... Done.### Can read about this: http://www.sigmath.es.osaka-u.ac.jp/shimo-lab/prog/pvclust/
plot(fit_Inh, las = 1, cex = 0.8, hang = -1)
## Want to get rid of info on
plot(fit_Inh, las = 1, cex = 0.8, hang = -1, print.num = F, col.pv = c("red", "white")) # dendogram with p values
# add rectangles around groups highly supported by the data
pvrect(fit_Inh, alpha=.95, cex = 0.5)
## Also use k-means clustering
### Determine number of clusters
### 'Scree plot' indicates 2 clusters
wss_Inh <- (nrow(dm_Inh)-1)*sum(apply(dm_Inh,2,var))
for (i in 2:5) wss_Inh[i] <- sum(kmeans(dm_Inh, centers=i)$withinss)
plot(1:5, wss_Inh, type="b", xlab="Number of Clusters",
ylab="Within groups sum of squares")
## Make 2 clusters, based on scree plot
## p-values also support 2 clusters
cl2_Inh <- kmeans(dm_Inh,2)
# append cluster assignment
dm_Inh <- data.frame(dm_Inh, cl2_Inh$cluster)
## Same groups k-means as ward's D
cbind(rownames(dm_Inh), dm_Inh[,11])## [,1] [,2]
## [1,] "BdBrazil-JEL649" "1"
## [2,] "BdGPL1-JEL404" "2"
## [3,] "BdGPL1-JEL647" "2"
## [4,] "BdGPL2-JEL423" "2"
## [5,] "BdGPL2-SRS810" "2"
## [6,] "Bsal" "1"groups_Inh## BdBrazil-JEL649 BdGPL1-JEL404 BdGPL1-JEL647 BdGPL2-JEL423
## 1 2 2 2
## BdGPL2-SRS810 Bsal
## 2 1Heatmap figure for clustering analysis
## I could probably make this file somehow with code, but just created and editted it in excel
### Used the df_Inh and table of values "S" to edit and make a table (averages) in Excel
## Table of values
S <- ddply(Bd_scores_samples, .(Bacteria, Temp, Isolate), summarize, N=length(Bacteria), mean=mean(final_inhibition), sd = sd(final_inhibition), min = min(final_inhibition), max = max(final_inhibition), se = sd/sqrt(N))
S## Bacteria Temp Isolate N mean sd min
## 1 B1 12 Bsal 12 0.462266512 0.333258081 -0.23213605
## 2 B1 12 JEL404 16 -0.366876144 0.107116751 -0.56862794
## 3 B1 12 JEL423 13 -0.668362275 0.261885544 -1.06654635
## 4 B1 12 JEL647 16 -0.521669479 0.340568903 -1.17003623
## 5 B1 12 JEL649 10 -0.479361894 0.550638999 -1.35421622
## 6 B1 12 SRS810 16 -0.801013138 0.360875036 -1.36117345
## 7 B1 18 Bsal 14 0.529327061 0.246559881 -0.01209231
## 8 B1 18 JEL404 14 -0.752547075 0.717179267 -2.01291316
## 9 B1 18 JEL423 16 -0.325082451 0.204891349 -0.75530290
## 10 B1 18 JEL647 16 -0.626372586 0.407569816 -1.53035781
## 11 B1 18 JEL649 15 -0.355133576 0.342895277 -1.19436694
## 12 B1 18 SRS810 15 -0.605081357 0.482501425 -1.70840624
## 13 B2 12 Bsal 16 0.646748471 0.125131112 0.38718365
## 14 B2 12 JEL404 16 0.367309143 0.120612213 0.12419636
## 15 B2 12 JEL423 16 0.105529211 0.195713116 -0.20862548
## 16 B2 12 JEL647 16 0.418429455 0.123111136 0.14331911
## 17 B2 12 JEL649 10 1.191009851 0.041153825 1.11590367
## 18 B2 12 SRS810 16 0.293953057 0.486844443 -0.66310713
## 19 B2 18 Bsal 16 0.474889867 0.312026860 -0.02211940
## 20 B2 18 JEL404 16 0.002819161 0.385951277 -0.47408996
## 21 B2 18 JEL423 16 0.245062255 0.086862895 0.05687141
## 22 B2 18 JEL647 16 0.051155016 0.153830309 -0.23627132
## 23 B2 18 JEL649 11 1.037419765 0.049472270 0.97055741
## 24 B2 18 SRS810 16 -0.509450158 0.328739178 -1.00142614
## 25 B3 12 Bsal 15 0.919445769 0.025617280 0.86777329
## 26 B3 12 JEL404 13 0.978617096 0.009286975 0.95793354
## 27 B3 12 JEL423 15 0.865513642 0.052984939 0.74468994
## 28 B3 12 JEL647 16 0.623290613 0.141726077 0.33507788
## 29 B3 12 JEL649 12 0.570738918 0.265303802 0.05154906
## 30 B3 12 SRS810 14 0.950102629 0.033471155 0.90136180
## 31 B3 18 Bsal 16 0.923156893 0.055646342 0.81565317
## 32 B3 18 JEL404 16 0.929116902 0.047454505 0.85050552
## 33 B3 18 JEL423 16 0.940899130 0.025810534 0.88620144
## 34 B3 18 JEL647 14 0.831066193 0.070483296 0.69024560
## 35 B3 18 JEL649 13 0.920243164 0.042464258 0.86829670
## 36 B3 18 SRS810 13 0.962728940 0.042562426 0.92415790
## 37 B4 12 Bsal 16 0.712984175 0.187608188 0.28516271
## 38 B4 12 JEL404 16 0.041553654 0.343663270 -0.76309651
## 39 B4 12 JEL423 16 -0.873211492 0.393015093 -1.32837746
## 40 B4 12 JEL647 16 -0.325001905 0.239342030 -0.92564734
## 41 B4 12 JEL649 11 0.383283348 0.400558052 -0.24850070
## 42 B4 12 SRS810 16 -0.841006894 0.369221599 -1.54488753
## 43 B4 18 Bsal 16 0.432001705 0.399427634 -0.42206164
## 44 B4 18 JEL404 16 -0.560666864 0.409472024 -1.19111684
## 45 B4 18 JEL423 16 -0.357994842 0.223907817 -0.73201822
## 46 B4 18 JEL647 16 0.049054590 0.208708051 -0.37847433
## 47 B4 18 JEL649 16 0.427238131 0.283139331 -0.17033747
## 48 B4 18 SRS810 16 -1.091347518 0.393867201 -1.62173057
## 49 B5 12 Bsal 16 0.923501937 0.039615091 0.85767646
## 50 B5 12 JEL404 15 1.001616384 0.026403157 0.95196880
## 51 B5 12 JEL423 13 0.880035717 0.120928142 0.70467851
## 52 B5 12 JEL647 9 0.713813213 0.255678466 0.22542953
## 53 B5 12 JEL649 9 0.702203834 0.308306419 0.13705715
## 54 B5 12 SRS810 14 0.978031761 0.058439103 0.88425887
## 55 B5 18 Bsal 12 0.913668974 0.061687373 0.77299584
## 56 B5 18 JEL404 11 0.957219207 0.012622955 0.93999577
## 57 B5 18 JEL423 14 0.955011300 0.069456411 0.84644825
## 58 B5 18 JEL647 12 0.961006689 0.126092236 0.65752800
## 59 B5 18 JEL649 15 1.001971857 0.053779384 0.91837609
## 60 B5 18 SRS810 13 1.004896155 0.072449711 0.89731858
## max se
## 1 0.807493330 0.096203321
## 2 -0.132265486 0.026779188
## 3 -0.175488304 0.072633981
## 4 -0.022992630 0.085142226
## 5 0.607682173 0.174127340
## 6 -0.199670967 0.090218759
## 7 0.801475376 0.065895900
## 8 -0.085138434 0.191674221
## 9 0.000923426 0.051222837
## 10 0.101671921 0.101892454
## 11 0.259958012 0.088535180
## 12 -0.065145222 0.124581332
## 13 0.870763802 0.031282778
## 14 0.560654690 0.030153053
## 15 0.380960075 0.048928279
## 16 0.603960318 0.030777784
## 17 1.264446563 0.013013982
## 18 0.841587957 0.121711111
## 19 0.841600789 0.078006715
## 20 0.494790089 0.096487819
## 21 0.381427620 0.021715724
## 22 0.262395991 0.038457577
## 23 1.100353217 0.014916451
## 24 0.062483109 0.082184795
## 25 0.967895038 0.006614353
## 26 0.987429313 0.002575743
## 27 0.940550720 0.013680652
## 28 0.875528175 0.035431519
## 29 0.922672193 0.076586611
## 30 1.011637535 0.008945542
## 31 0.986067143 0.013911585
## 32 0.983286110 0.011863626
## 33 0.975639647 0.006452634
## 34 0.930075931 0.018837453
## 35 0.982588079 0.011777466
## 36 1.094017473 0.011804693
## 37 0.898369341 0.046902047
## 38 0.491454649 0.085915818
## 39 -0.192884153 0.098253773
## 40 0.016657359 0.059835508
## 41 0.914123197 0.120772797
## 42 -0.391517801 0.092305400
## 43 0.816951006 0.099856909
## 44 -0.079061500 0.102368006
## 45 -0.031962182 0.055976954
## 46 0.273325811 0.052177013
## 47 0.709067762 0.070784833
## 48 -0.485187864 0.098466800
## 49 0.985893178 0.009903773
## 50 1.029244822 0.006817266
## 51 1.092662550 0.033539432
## 52 1.092179253 0.085226155
## 53 1.091868005 0.102768806
## 54 1.052775754 0.015618507
## 55 0.985934964 0.017807611
## 56 0.979093981 0.003805964
## 57 1.030265622 0.018563007
## 58 1.075188263 0.036399693
## 59 1.069673808 0.013885777
## 60 1.125437531 0.020093934#write.csv(S, "stats_code//inhition_values_Bac_Bd_temp_meanmaxmin.csv")
## We brought this in in the first lines of code, I can't knit the .Rmd file if I bring it in here
## averages <- read.csv("stats_code/Final_code_from_Gabe/heatmap_averages_final.csv", header = T, row.names = 1)
# creates a own color palette from green to red
my_palette <- colorRampPalette(c("green", "yellow", "orange", "red"))(n = 399)
# (optional) defines the color breaks manually for a "skewed" color transition
col_breaks = c(seq(-110,0,length=100), # for green
seq(0.1,30, length = 100), # for yellow
seq(30.1,70,length=100), # for orange
seq(70.1,110,length=100)) # for red
avg_matrix <- data.matrix(averages)
## Error here when knitting, have to run code instead
distance <- dist(avg_matrix, method = "euclidean")
cluster <- hclust(distance, method = "ward.D2")
## Knitting the .Rmd gets upset at the margins on this figure, so will just have to run without making this figure
## This is the final figure for the heatmap in the publication
## heatmap.2(avg_matrix, Colv = F, Rowv = as.dendrogram(cluster), dendrogram = "row", margins = c(8,13), col = my_palette, breaks = col_breaks, trace = "none", density.info="none", lhei = c(1,4), cexRow = 1, cexCol = 1, key.xlab = "Percent inhibition", key.title = "", keysize = 1.6)
## Put in au scores from the pvclust analysis above
#text(locator(1),"98 *", cex = 0.9)
#text(locator(1),"95 *", cex =0.9)
#text(locator(1),"81", cex = 0.7)
#text(locator(1),"64", cex = 0.7)bacteria facultative???
## Have to compare the inhibition score of PC vs inhibition score of the sample
#Bd_scores <- read.csv("Final_plate_reads/Inhibition_scores_final.csv")
Bd_scores2 <- Bd_scores[Bd_scores$Bacteria %in% c("B1", "B2", "B3", "B4", "B5", "PC", "ND"),]
Bd_scores2$Temp <- as.factor(Bd_scores2$Temp)
ggplot(Bd_scores2, aes(x = Bacteria, y = final_inhibition_percent)) + geom_boxplot(aes(color = Genotype)) + facet_grid(.~Temp) + ylab("Percent inhibition") + scale_x_discrete(labels=c("Ac", "Ps1", "Ps2", "Chy", "Sten", "NDPC", "PC"))
ggplot(Bd_scores2, aes(x = Bacteria, y = final_inhibition_percent)) + geom_boxplot(aes(color = Temp)) + facet_grid(.~Genotype) + ylab("Percent inhibition") 
Bd_scores3 <- Bd_scores2[Bd_scores2$Temp == "12",]
ggplot(Bd_scores3, aes(x = Genotype, y = final_inhibition_percent)) + geom_boxplot(aes(color = Bacteria)) + facet_grid(.~Genotype, scales = "free_x") + ylab("Percent inhibition\n") + xlab("\nGenotype") + scale_color_discrete(labels=c("Ac", "Ps1", "Ps2", "Chy", "Sten", "NDPC", "PC"), guide = guide_legend(title = "Well type")) + theme(strip.text.x = element_blank())
Bd_scores4 <- Bd_scores2[Bd_scores2$Temp == "18",]
ggplot(Bd_scores4, aes(x = Genotype, y = final_inhibition_percent)) + geom_boxplot(aes(color = Bacteria)) + facet_grid(.~Genotype, scales = "free_x") + ylab("Percent inhibition\n") + xlab("\nGenotype") + scale_color_discrete(labels=c("Ac", "Ps1", "Ps2", "Chy", "Sten", "NDPC", "PC"), guide = guide_legend(title = "Well type")) + theme(strip.text.x = element_blank())
## Test each pathogen genotype
Bsal <- Bd_scores2[Bd_scores2$Isolate == "Bsal",]
J649 <- Bd_scores2[Bd_scores2$Isolate == "JEL649",]
J647 <- Bd_scores2[Bd_scores2$Isolate == "JEL647",]
SRS <- Bd_scores2[Bd_scores2$Isolate == "SRS810",]
J423 <- Bd_scores2[Bd_scores2$Isolate == "JEL423",]
J404 <- Bd_scores2[Bd_scores2$Isolate == "JEL404",]
## Bacteria factor is bacteria and pathogen wells -- bacteria 1-5, NDPC and PC, looking for bacteria effect. If bacteria is significant, it may be that one of the bacteria are different from the PC
## Bsal ##
m1_Bsal <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = Bsal, REML = F)
Anova(m1_Bsal, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 704.4531 6 < 2e-16 ***
## Temp 0.7613 1 0.38292
## Bacteria:Temp 13.7379 6 0.03271 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_Bsal <- lsmeans(m1_Bsal, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_Bsal)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.192373500 0.09395458 201.60 -2.048 0.3883
## B1 - B3 -0.459601923 0.09517695 201.35 -4.829 0.0001
## B1 - B4 -0.258609205 0.09395458 201.60 -2.752 0.0909
## B1 - B5 -0.469126966 0.09395458 201.60 -4.993 <.0001
## B1 - ND 0.462266512 0.10003997 200.72 4.621 0.0001
## B1 - PC 0.965993683 0.10003997 200.72 9.656 <.0001
## B2 - B3 -0.267228423 0.08812534 200.86 -3.032 0.0430
## B2 - B4 -0.066235705 0.08663716 200.72 -0.765 0.9880
## B2 - B5 -0.276753466 0.08663716 200.72 -3.194 0.0267
## B2 - ND 0.654640012 0.09395458 201.60 6.968 <.0001
## B2 - PC 1.158367183 0.09395458 201.60 12.329 <.0001
## B3 - B4 0.200992718 0.08812534 200.86 2.281 0.2585
## B3 - B5 -0.009525043 0.08812534 200.86 -0.108 1.0000
## B3 - ND 0.921868435 0.09517695 201.35 9.686 <.0001
## B3 - PC 1.425595606 0.09517695 201.35 14.978 <.0001
## B4 - B5 -0.210517761 0.08663716 200.72 -2.430 0.1917
## B4 - ND 0.720875717 0.09395458 201.60 7.673 <.0001
## B4 - PC 1.224602888 0.09395458 201.60 13.034 <.0001
## B5 - ND 0.931393479 0.09395458 201.60 9.913 <.0001
## B5 - PC 1.435120649 0.09395458 201.60 15.275 <.0001
## ND - PC 0.503727171 0.10003997 200.72 5.035 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 0.067360804 0.08994143 201.38 0.749 0.9892
## B1 - B3 -0.380906222 0.08994143 201.38 -4.235 0.0007
## B1 - B4 0.110248966 0.08994143 201.38 1.226 0.8834
## B1 - B5 -0.393841077 0.09668113 201.35 -4.074 0.0013
## B1 - ND 0.519827896 0.09668113 201.35 5.377 <.0001
## B1 - PC 1.047686855 0.09901266 201.27 10.581 <.0001
## B2 - B3 -0.448267026 0.08663716 200.72 -5.174 <.0001
## B2 - B4 0.042888162 0.08663716 200.72 0.495 0.9989
## B2 - B5 -0.461201881 0.09396112 201.63 -4.908 <.0001
## B2 - ND 0.452467093 0.09396112 201.63 4.815 0.0001
## B2 - PC 0.980326051 0.09627735 201.40 10.182 <.0001
## B3 - B4 0.491155189 0.08663716 200.72 5.669 <.0001
## B3 - B5 -0.012934855 0.09396112 201.63 -0.138 1.0000
## B3 - ND 0.900734119 0.09396112 201.63 9.586 <.0001
## B3 - PC 1.428593077 0.09627735 201.40 14.838 <.0001
## B4 - B5 -0.504090044 0.09396112 201.63 -5.365 <.0001
## B4 - ND 0.409578930 0.09396112 201.63 4.359 0.0004
## B4 - PC 0.937437888 0.09627735 201.40 9.737 <.0001
## B5 - ND 0.913668974 0.10003997 200.72 9.133 <.0001
## B5 - PC 1.441527932 0.10234954 200.84 14.084 <.0001
## ND - PC 0.527858958 0.10234954 200.84 5.157 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## NDPC and PC are different than all the bacteria, with greater inhibition than NDPC and PC
ggplot(Bsal, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## J649 ##
m1_J649 <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = J649, REML = F)
Anova(m1_J649, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 1435.7418 6 < 2.2e-16 ***
## Temp 7.9096 1 0.004917 **
## Bacteria:Temp 72.5828 6 1.206e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_J649 <- lsmeans(m1_J649, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_J649)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -1.61858761 0.1632539 185.08 -9.915 <.0001
## B1 - B3 -1.01557805 0.1557541 184.45 -6.520 <.0001
## B1 - B4 -0.84381465 0.1585614 184.01 -5.322 <.0001
## B1 - B5 -1.18242105 0.1666949 183.95 -7.093 <.0001
## B1 - ND -0.44483914 0.1557541 184.45 -2.856 0.0701
## B1 - PC 2.24740066 0.1557541 184.45 14.429 <.0001
## B2 - B3 0.60300955 0.1553649 183.99 3.881 0.0027
## B2 - B4 0.77477296 0.1588734 184.38 4.877 <.0001
## B2 - B5 0.43616656 0.1675374 184.89 2.603 0.1309
## B2 - ND 1.17374847 0.1553649 183.99 7.555 <.0001
## B2 - PC 3.86598827 0.1553649 183.99 24.883 <.0001
## B3 - B4 0.17176341 0.1514474 183.96 1.134 0.9168
## B3 - B5 -0.16684300 0.1603231 184.35 -1.041 0.9438
## B3 - ND 0.57073892 0.1480108 183.83 3.856 0.0030
## B3 - PC 3.26297871 0.1480108 183.83 22.046 <.0001
## B4 - B5 -0.33860640 0.1631276 184.02 -2.076 0.3715
## B4 - ND 0.39897551 0.1514474 183.96 2.634 0.1218
## B4 - PC 3.09121530 0.1514474 183.96 20.411 <.0001
## B5 - ND 0.73758191 0.1603231 184.35 4.601 0.0002
## B5 - PC 3.42982171 0.1603231 184.35 21.393 <.0001
## ND - PC 2.69223979 0.1480108 183.83 18.189 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -1.37687879 0.1444970 184.54 -9.529 <.0001
## B1 - B3 -1.27680197 0.1385165 185.24 -9.218 <.0001
## B1 - B4 -0.77775436 0.1303745 183.93 -5.966 <.0001
## B1 - B5 -1.34658658 0.1325864 184.09 -10.156 <.0001
## B1 - ND -0.35051623 0.1303745 183.93 -2.689 0.1070
## B1 - PC 1.07122033 0.1303745 183.93 8.216 <.0001
## B2 - B3 0.10007681 0.1508450 186.37 0.663 0.9944
## B2 - B4 0.59912443 0.1424729 184.42 4.205 0.0008
## B2 - B5 0.03029221 0.1442735 184.26 0.210 1.0000
## B2 - ND 1.02636256 0.1424729 184.42 7.204 <.0001
## B2 - PC 2.44809911 0.1424729 184.42 17.183 <.0001
## B3 - B4 0.49904761 0.1362715 184.98 3.662 0.0059
## B3 - B5 -0.06978460 0.1384851 185.21 -0.504 0.9988
## B3 - ND 0.92628574 0.1362715 184.98 6.797 <.0001
## B3 - PC 2.34802230 0.1362715 184.98 17.230 <.0001
## B4 - B5 -0.56883222 0.1303762 183.93 -4.363 0.0004
## B4 - ND 0.42723813 0.1281811 183.83 3.333 0.0176
## B4 - PC 1.84897468 0.1281811 183.83 14.425 <.0001
## B5 - ND 0.99607035 0.1303762 183.93 7.640 <.0001
## B5 - PC 2.41780690 0.1303762 183.93 18.545 <.0001
## ND - PC 1.42173655 0.1281811 183.83 11.092 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## greater inhibition than PC
ggplot(J649, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## J647 ##
m1_J647 <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = J647, REML = F)
Anova(m1_J647, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 2067.5726 6 <2e-16 ***
## Temp 0.5762 1 0.4478
## Bacteria:Temp 107.6639 6 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_J647 <- lsmeans(m1_J647, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_J647)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.940098935 0.07903079 215.73 -11.895 <.0001
## B1 - B3 -1.144960092 0.07903079 215.73 -14.488 <.0001
## B1 - B4 -0.196667575 0.07903079 215.73 -2.488 0.1687
## B1 - B5 -1.260355949 0.09350261 216.55 -13.479 <.0001
## B1 - ND -0.521669479 0.07903079 215.73 -6.601 <.0001
## B1 - PC 0.514000479 0.07903079 215.73 6.504 <.0001
## B2 - B3 -0.204861157 0.07903079 215.73 -2.592 0.1335
## B2 - B4 0.743431360 0.07903079 215.73 9.407 <.0001
## B2 - B5 -0.320257014 0.09350261 216.55 -3.425 0.0128
## B2 - ND 0.418429455 0.07903079 215.73 5.295 <.0001
## B2 - PC 1.454099413 0.07903079 215.73 18.399 <.0001
## B3 - B4 0.948292517 0.07903079 215.73 11.999 <.0001
## B3 - B5 -0.115395857 0.09350261 216.55 -1.234 0.8801
## B3 - ND 0.623290613 0.07903079 215.73 7.887 <.0001
## B3 - PC 1.658960571 0.07903079 215.73 20.991 <.0001
## B4 - B5 -1.063688375 0.09350261 216.55 -11.376 <.0001
## B4 - ND -0.325001905 0.07903079 215.73 -4.112 0.0011
## B4 - PC 0.710668053 0.07903079 215.73 8.992 <.0001
## B5 - ND 0.738686470 0.09350261 216.55 7.900 <.0001
## B5 - PC 1.774356428 0.09350261 216.55 18.977 <.0001
## ND - PC 1.035669958 0.07903079 215.73 13.105 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.677527602 0.07903079 215.73 -8.573 <.0001
## B1 - B3 -1.460447705 0.08199765 216.21 -17.811 <.0001
## B1 - B4 -0.675427176 0.07903079 215.73 -8.546 <.0001
## B1 - B5 -1.583832895 0.08552149 216.13 -18.520 <.0001
## B1 - ND -0.626372586 0.07903079 215.73 -7.926 <.0001
## B1 - PC 0.946416899 0.07903079 215.73 11.975 <.0001
## B2 - B3 -0.782920103 0.08199765 216.21 -9.548 <.0001
## B2 - B4 0.002100426 0.07903079 215.73 0.027 1.0000
## B2 - B5 -0.906305293 0.08552149 216.13 -10.597 <.0001
## B2 - ND 0.051155016 0.07903079 215.73 0.647 0.9951
## B2 - PC 1.623944501 0.07903079 215.73 20.548 <.0001
## B3 - B4 0.785020528 0.08199765 216.21 9.574 <.0001
## B3 - B5 -0.123385190 0.08853384 217.12 -1.394 0.8048
## B3 - ND 0.834075118 0.08199765 216.21 10.172 <.0001
## B3 - PC 2.406864603 0.08199765 216.21 29.353 <.0001
## B4 - B5 -0.908405719 0.08552149 216.13 -10.622 <.0001
## B4 - ND 0.049054590 0.07903079 215.73 0.621 0.9961
## B4 - PC 1.621844075 0.07903079 215.73 20.522 <.0001
## B5 - ND 0.957460309 0.08552149 216.13 11.196 <.0001
## B5 - PC 2.530249794 0.08552149 216.13 29.586 <.0001
## ND - PC 1.572789485 0.07903079 215.73 19.901 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## greater inhibition than PC
ggplot(J647, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## SRS ##
m1_SRS <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = SRS, REML = F)
Anova(m1_SRS, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 1660.4596 6 < 2.2e-16 ***
## Temp 0.8623 1 0.3531
## Bacteria:Temp 42.4055 6 1.529e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_SRS <- lsmeans(m1_SRS, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_SRS)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -1.09496619 0.1235889 217.75 -8.860 <.0001
## B1 - B3 -1.77950866 0.1280173 217.79 -13.901 <.0001
## B1 - B4 0.03999376 0.1235889 217.75 0.324 0.9999
## B1 - B5 -1.76938152 0.1280224 217.79 -13.821 <.0001
## B1 - ND -0.80101314 0.1235889 217.75 -6.481 <.0001
## B1 - PC 0.99870196 0.1235889 217.75 8.081 <.0001
## B2 - B3 -0.68454247 0.1280173 217.79 -5.347 <.0001
## B2 - B4 1.13495995 0.1235889 217.75 9.183 <.0001
## B2 - B5 -0.67441532 0.1280224 217.79 -5.268 <.0001
## B2 - ND 0.29395306 0.1235889 217.75 2.378 0.2128
## B2 - PC 2.09366815 0.1235889 217.75 16.941 <.0001
## B3 - B4 1.81950242 0.1280173 217.79 14.213 <.0001
## B3 - B5 0.01012715 0.1323024 217.83 0.077 1.0000
## B3 - ND 0.97849553 0.1280173 217.79 7.643 <.0001
## B3 - PC 2.77821062 0.1280173 217.79 21.702 <.0001
## B4 - B5 -1.80937527 0.1280224 217.79 -14.133 <.0001
## B4 - ND -0.84100689 0.1235889 217.75 -6.805 <.0001
## B4 - PC 0.95870820 0.1235889 217.75 7.757 <.0001
## B5 - ND 0.96836838 0.1280224 217.79 7.564 <.0001
## B5 - PC 2.76808347 0.1280224 217.79 21.622 <.0001
## ND - PC 1.79971509 0.1235889 217.75 14.562 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.10316266 0.1256993 217.79 -0.821 0.9826
## B1 - B3 -1.53319756 0.1325820 217.81 -11.564 <.0001
## B1 - B4 0.47873470 0.1256993 217.79 3.809 0.0034
## B1 - B5 -1.62674505 0.1327700 217.90 -12.252 <.0001
## B1 - ND -0.61261281 0.1256993 217.79 -4.874 <.0001
## B1 - PC 1.37199786 0.1256993 217.79 10.915 <.0001
## B2 - B3 -1.43003491 0.1308308 217.90 -10.930 <.0001
## B2 - B4 0.58189736 0.1235889 217.75 4.708 0.0001
## B2 - B5 -1.52358239 0.1308271 217.89 -11.646 <.0001
## B2 - ND -0.50945016 0.1235889 217.75 -4.122 0.0010
## B2 - PC 1.47516052 0.1235889 217.75 11.936 <.0001
## B3 - B4 2.01193227 0.1308308 217.90 15.378 <.0001
## B3 - B5 -0.09354749 0.1377330 218.04 -0.679 0.9936
## B3 - ND 0.92058475 0.1308308 217.90 7.036 <.0001
## B3 - PC 2.90519543 0.1308308 217.90 22.206 <.0001
## B4 - B5 -2.10547975 0.1308271 217.89 -16.094 <.0001
## B4 - ND -1.09134752 0.1235889 217.75 -8.830 <.0001
## B4 - PC 0.89326316 0.1235889 217.75 7.228 <.0001
## B5 - ND 1.01413224 0.1308271 217.89 7.752 <.0001
## B5 - PC 2.99874291 0.1308271 217.89 22.921 <.0001
## ND - PC 1.98461068 0.1235889 217.75 16.058 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## greater inhibition than PC
ggplot(SRS, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## J423 ##
m1_J423 <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = J423, REML = F)
Anova(m1_J423, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 1678.502 6 < 2.2e-16 ***
## Temp 30.675 1 3.05e-08 ***
## Bacteria:Temp 96.704 6 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_J423 <- lsmeans(m1_J423, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_J423)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.77928818 0.1117205 219.86 -6.975 <.0001
## B1 - B3 -1.53422192 0.1133810 219.94 -13.532 <.0001
## B1 - B4 0.19945253 0.1117205 219.86 1.785 0.5594
## B1 - B5 -1.55876034 0.1178139 221.65 -13.231 <.0001
## B1 - ND -0.67375896 0.1117205 219.86 -6.031 <.0001
## B1 - PC 1.46460041 0.1117205 219.86 13.110 <.0001
## B2 - B3 -0.75493374 0.1075199 219.81 -7.021 <.0001
## B2 - B4 0.97874070 0.1057146 219.56 9.258 <.0001
## B2 - B5 -0.77947217 0.1118993 220.70 -6.966 <.0001
## B2 - ND 0.10552921 0.1057146 219.56 0.998 0.9540
## B2 - PC 2.24388858 0.1057146 219.56 21.226 <.0001
## B3 - B4 1.73367445 0.1075199 219.81 16.124 <.0001
## B3 - B5 -0.02453842 0.1135786 220.66 -0.216 1.0000
## B3 - ND 0.86046295 0.1075199 219.81 8.003 <.0001
## B3 - PC 2.99882232 0.1075199 219.81 27.891 <.0001
## B4 - B5 -1.75821287 0.1118993 220.70 -15.712 <.0001
## B4 - ND -0.87321149 0.1057146 219.56 -8.260 <.0001
## B4 - PC 1.26514788 0.1057146 219.56 11.968 <.0001
## B5 - ND 0.88500138 0.1118993 220.70 7.909 <.0001
## B5 - PC 3.02336075 0.1118993 220.70 27.019 <.0001
## ND - PC 2.13835937 0.1057146 219.56 20.228 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.57014471 0.1057146 219.56 -5.393 <.0001
## B1 - B3 -1.26598158 0.1057146 219.56 -11.975 <.0001
## B1 - B4 0.03291239 0.1057146 219.56 0.311 0.9999
## B1 - B5 -1.28045626 0.1095037 219.94 -11.693 <.0001
## B1 - ND -0.32508245 0.1057146 219.56 -3.075 0.0377
## B1 - PC 0.63324648 0.1057146 219.56 5.990 <.0001
## B2 - B3 -0.69583688 0.1057146 219.56 -6.582 <.0001
## B2 - B4 0.60305710 0.1057146 219.56 5.705 <.0001
## B2 - B5 -0.71031156 0.1095037 219.94 -6.487 <.0001
## B2 - ND 0.24506226 0.1057146 219.56 2.318 0.2401
## B2 - PC 1.20339119 0.1057146 219.56 11.383 <.0001
## B3 - B4 1.29889397 0.1057146 219.56 12.287 <.0001
## B3 - B5 -0.01447468 0.1095037 219.94 -0.132 1.0000
## B3 - ND 0.94089913 0.1057146 219.56 8.900 <.0001
## B3 - PC 1.89922806 0.1057146 219.56 17.966 <.0001
## B4 - B5 -1.31336865 0.1095037 219.94 -11.994 <.0001
## B4 - ND -0.35799484 0.1057146 219.56 -3.386 0.0145
## B4 - PC 0.60033409 0.1057146 219.56 5.679 <.0001
## B5 - ND 0.95537381 0.1095037 219.94 8.725 <.0001
## B5 - PC 1.91370274 0.1095037 219.94 17.476 <.0001
## ND - PC 0.95832893 0.1057146 219.56 9.065 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## greater inhibition than PC
ggplot(J423, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 
## JEL404 ##
m1_404 <- lmer(final_inhibition ~ Bacteria*Temp + (1 | Incubator/Plate), data = J404, REML = F)
Anova(m1_404, type = "II", test = "Chisq")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: final_inhibition
## Chisq Df Pr(>Chisq)
## Bacteria 1094.9804 6 < 2.2e-16 ***
## Temp 2.5614 1 0.1095019
## Bacteria:Temp 23.9712 6 0.0005287 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## this indicates a significant interaction
## use lsmeans to look at what the interaction
p1_404 <- lsmeans(m1_404, pairwise~Bacteria|Temp, adjust = "tukey")
pairs(p1_404)## Temp = 12:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.734185288 0.1181327 217.73 -6.215 <.0001
## B1 - B3 -1.331994000 0.1250397 217.91 -10.653 <.0001
## B1 - B4 -0.408429799 0.1181327 217.73 -3.457 0.0115
## B1 - B5 -1.369318255 0.1201444 217.77 -11.397 <.0001
## B1 - ND -0.366876144 0.1181327 217.73 -3.106 0.0345
## B1 - PC 0.888914748 0.1181327 217.73 7.525 <.0001
## B2 - B3 -0.597808713 0.1250397 217.91 -4.781 0.0001
## B2 - B4 0.325755489 0.1181327 217.73 2.758 0.0894
## B2 - B5 -0.635132968 0.1201444 217.77 -5.286 <.0001
## B2 - ND 0.367309143 0.1181327 217.73 3.109 0.0342
## B2 - PC 1.623100036 0.1181327 217.73 13.740 <.0001
## B3 - B4 0.923564202 0.1250397 217.91 7.386 <.0001
## B3 - B5 -0.037324255 0.1270688 218.02 -0.294 0.9999
## B3 - ND 0.965117856 0.1250397 217.91 7.718 <.0001
## B3 - PC 2.220908748 0.1250397 217.91 17.762 <.0001
## B4 - B5 -0.960888457 0.1201444 217.77 -7.998 <.0001
## B4 - ND 0.041553654 0.1181327 217.73 0.352 0.9998
## B4 - PC 1.297344547 0.1181327 217.73 10.982 <.0001
## B5 - ND 1.002442111 0.1201444 217.77 8.344 <.0001
## B5 - PC 2.258233003 0.1201444 217.77 18.796 <.0001
## ND - PC 1.255790892 0.1181327 217.73 10.630 <.0001
##
## Temp = 18:
## contrast estimate SE df t.ratio p.value
## B1 - B2 -0.793639488 0.1225721 217.93 -6.475 <.0001
## B1 - B3 -1.719937229 0.1225721 217.93 -14.032 <.0001
## B1 - B4 -0.230153463 0.1225721 217.93 -1.878 0.4973
## B1 - B5 -1.676827510 0.1349328 217.90 -12.427 <.0001
## B1 - ND -0.790820327 0.1225721 217.93 -6.452 <.0001
## B1 - PC 0.531583143 0.1225721 217.93 4.337 0.0004
## B2 - B3 -0.926297741 0.1181327 217.73 -7.841 <.0001
## B2 - B4 0.563486025 0.1181327 217.73 4.770 0.0001
## B2 - B5 -0.883188022 0.1312392 217.95 -6.730 <.0001
## B2 - ND 0.002819161 0.1181327 217.73 0.024 1.0000
## B2 - PC 1.325222631 0.1181327 217.73 11.218 <.0001
## B3 - B4 1.489783766 0.1181327 217.73 12.611 <.0001
## B3 - B5 0.043109719 0.1312392 217.95 0.328 0.9999
## B3 - ND 0.929116902 0.1181327 217.73 7.865 <.0001
## B3 - PC 2.251520372 0.1181327 217.73 19.059 <.0001
## B4 - B5 -1.446674046 0.1312392 217.95 -11.023 <.0001
## B4 - ND -0.560666864 0.1181327 217.73 -4.746 0.0001
## B4 - PC 0.761736606 0.1181327 217.73 6.448 <.0001
## B5 - ND 0.886007182 0.1312392 217.95 6.751 <.0001
## B5 - PC 2.208410653 0.1312392 217.95 16.827 <.0001
## ND - PC 1.322403470 0.1181327 217.73 11.194 <.0001
##
## P value adjustment: tukey method for comparing a family of 7 estimates## greater inhibition than PC
ggplot(J404, aes(x = Bacteria, y = final_inhibition)) + geom_boxplot(aes(color = Isolate)) +facet_grid(.~Temp) + ylab("% inhibition") + theme(text = element_text(size = 18)) 