Finalize code for first submission of manuscript

code/Rfunctions.R

@@ -13,7 +13,7 @@ loadPhyloseq <- function(){
     #Remove non-target taxa
     prune_taxa(taxa_names(.)!="Melampsora" & !grepl("OTU",taxa_names(.)),.) %>%
     #Remove two outlier samples that cause problems for model fitting
-    prune_samples(!(grepl("G6.T3.R4.TP1",sample_names(.))),.) %>%
+    #prune_samples(!(grepl("G6.T3.R4.TP1",sample_names(.))),.) #%>%
     #This sample has an order of magnitude greater Trichoderma which causes instability in the mvabund model
     prune_samples(!(grepl("G4.T2.R5.TP1",sample_names(.))),.)
 }
@@ -32,8 +32,6 @@ loadRust <- function(){
   return(dat)
     }
 
-
-
 #______________________________________________#
 ### apply bias correction to phyloseq object ###
 unbias <- function(phy.in,bias){
@@ -81,7 +79,7 @@ annotation_custom2 <- function (grob, xmin = -Inf, xmax = Inf, ymin = -Inf, ymax
 }
 
 #_____________________________#
-### alculate geometric mean ###
+### calculate geometric mean ###
 gm_mean <- function(x, na_rm = FALSE) {
   exp(mean(log(x), na.rm = na_rm))
 }

code/compile.R

@@ -85,9 +85,9 @@ phy %<>% prune_taxa(!contam.dat$contam,.)
 phy %<>% prune_samples(!grepl("NEG",sample_names(.)),.)
 
 #' ### Look at sequencing depth
-minDepth <- 4000
-data.frame(SeqDepth=sample_sums(phy), TP=factor(sample_data(phy)$Timepoint)) %>%
-  mutate(cutoff=SeqDepth>minDepth) %>%
+minDepth <- 5000
+data.frame(SeqDepth=sample_sums(phy), TP=factor(sample_data(phy)$Timepoint),Type=sample_data(phy)$Treatment ) %>%
+  mutate(cutoff=SeqDepth>minDepth,neg=ifelse(Type=="Negative","neg","pos"),TP=paste(TP,neg)) %>%
   ggplot(aes(x=TP, y=SeqDepth)) +
   geom_violin() +
   geom_point(aes(color=cutoff),position=position_jitter(width=0.1)) +

code/jsdModels.R

@@ -5,7 +5,7 @@ library(mvabund)
 source("code/Rfunctions.R")
 
 #set seed for reproducability
-set.seed(32453)
+set.seed(32576)
 
 # load phyloseq data
 (phy <- loadPhyloseq())
@@ -46,7 +46,7 @@ mv.full <- manyglm(mvDat ~ Genotype*Treatment,
 # Test with anova.manyglm 
 # Using unstructured correlation matrix and wald tests.  
 # Including univariate test with adjustment for multiple testing.  
-mv.anova <- anova(mv.full, nBoot=4999, p.uni="adjusted", cor.type="R", test="wald")
+mv.anova <- anova(mv.full, nBoot=4999, p.uni="adjusted", cor.type="shrink", test="wald")
 
 # Save results 
 mv.anova$table %>% write.csv("output/tabs/mv.genotype.csv")
@@ -63,12 +63,12 @@ mv.region <- manyglm(mvDat ~ Region*Treatment,
                      data=data.frame(sample_data(phy)))
 
 #' ## Check model assumptions
-plot(mv.region)
-meanvar.plot(mvDat~sample_data(phy)$Region)
+#plot(mv.region)
+#meanvar.plot(mvDat~sample_data(phy)$Region)
 
 #' ## Test with anova.manyglm 
 #+ cache=T, results='asis'
-mv.region.anova <- anova(mv.region, nBoot=4999, p.uni="adjusted", cor.type="R", test="wald")
+mv.region.anova <- anova(mv.region, nBoot=4999, p.uni="adjusted", cor.type="shrink", test="wald")
 
 mv.region.anova$table %>% write.csv("output/tabs/mv.region.csv")
 saveRDS(mv.region.anova, "output/rds/mv.region.rds")

code/mapFigS1.R

@@ -21,6 +21,6 @@ map <- make_bbox(points$Longitude,points$Latitud,f=0.35) %>%
 ggdraw()+
   draw_plot(map)+
   draw_image("data/500px-Populus_trichocarpa_range_map.svg.png",scale=1,width=0.32,x=0.52,y=0.24)
-ggsave("output/figs/Figs.S1.pdf",width=20,height=16,units="cm")
+ggsave("output/figs/Fig.S1.pdf",width=20,height=16,units="cm")
 
 

code/ordinations.R

@@ -44,7 +44,7 @@ dbrda.treatment.r2 <- (eigenvals(dbrda.treatment)/sum(eigenvals(dbrda.treatment)
     ggplot(aes(x=x,y=y,fill=Treatment)) +
     geom_errorbar(aes(ymin=y-sd2,ymax=y+sd2),width=0)+
     geom_errorbarh(aes(xmin=x-sd1,xmax=x+sd1),height=0)+
-    geom_point(size=5,shape=21)+
+    geom_point(size=4,shape=21)+
     annotate("text", x = -Inf, y = Inf, label = "(a)",hjust=-0.4, vjust=1.75, size=5)+
     scale_fill_manual(values=pal.treatment)+
     coord_fixed(eigenvals(dbrda.treatment)[2]/eigenvals(dbrda.treatment)[1])+
@@ -74,7 +74,7 @@ dbrda.genotype.r2 <- (eigenvals(dbrda.genotype)/sum(eigenvals(dbrda.genotype)))
     ggplot(aes(x=x,y=y,fill=Region)) +
     geom_errorbar(aes(ymin=y-sd2,ymax=y+sd2),width=0)+
     geom_errorbarh(aes(xmin=x-sd1,xmax=x+sd1),height=0)+
-    geom_point(size=5,shape=21)+
+    geom_point(size=4,shape=21)+
     annotate("text", x = -Inf, y = Inf, label = "(b)",hjust=-0.4, vjust=1.75, size=5)+
     scale_fill_manual(values=pal.region)+
     coord_fixed(eigenvals(dbrda.genotype)[2]/eigenvals(dbrda.genotype)[1])+

code/rustCor.R

@@ -82,7 +82,7 @@ dat$variable_lab %<>% factor(.,levels=unique(.)[c(grep("*",.,fixed=T,invert = T)
 # Merge correlation test results and format for markdown
 cor.dat <- dat %>% select(variable,variable_lab) %>% unique() %>%
   left_join(pearsons) %>%
-  mutate(stats=paste0("*p*=",pval,"; ","*r*^(2)=",cor))
+  mutate(stats=paste0("*p*=",pval,"; ","*r*=",cor))
 #make figure
 ggplot(dat, aes(x=value,y=resid))+
   geom_point()+

output/compiled/taxonomy.table.csv

@@ -53,6 +53,7 @@
 "Aureobasidium","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Dothideales","f__Aureobasidiaceae",NA,NA
 "OTU.103","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Filobasidiales","f__Filobasidiaceae","g__Filobasidium","s__wieringae"
 "OTU.79","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Tremellales","f__Bulleribasidiaceae","g__Vishniacozyma","s__victoriae"
+"OTU.88","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales",NA,NA,NA
 "OTU.255","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Cystofilobasidiales","f__Cystofilobasidiaceae","g__Cystofilobasidium","s__macerans"
 "OTU.140","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Cantharellales_fam_Incertae_sedis","g__Burgoa",NA
 "OTU.93","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Polyporales","f__Phanerochaetaceae","g__Hyphodermella","s__corrugata"
@@ -63,13 +64,12 @@
 "OTU.321","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Trechisporales","f__Hydnodontaceae","g__Brevicellicium","s__olivascens"
 "OTU.226","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Cantharellales_fam_Incertae_sedis","g__Sistotrema","s__brinkmannii"
 "OTU.182","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Cantharellales_fam_Incertae_sedis","g__Sistotrema","s__oblongisporum"
-"OTU.125","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales",NA,NA,NA
 "OTU.291","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Aspergillus",NA
 "OTU.181","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Ceratobasidiaceae",NA,NA
 "OTU.316","k__Fungi","p__Ascomycota","c__Leotiomycetes","o__Phacidiales","f__Phacidiaceae","g__Phacidium","s__pseudophacidioides"
 "OTU.284","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales","f__Cordycipitaceae","g__Simplicillium","s__cylindrosporum"
-"OTU.102","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Polyporales","f__Steccherinaceae","g__Steccherinum",NA
-"OTU.308","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Teratosphaeriaceae","g__Catenulostroma",NA
+"OTU.63","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Pleosporales","f__Phaeosphaeriaceae","g__Phaeosphaeria","s__tofieldiae"
+"OTU.308","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Teratosphaeriaceae",NA,NA
 "OTU.188","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Aspergillus","s__penicillioides"
 "OTU.264","k__Fungi","p__Basidiomycota","c__Exobasidiomycetes","o__Golubeviales","f__Golubeviaceae","g__Golubevia","s__pallescens"
 "OTU.292","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Agaricales",NA,NA,NA
@@ -78,7 +78,7 @@
 "OTU.281","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Botryobasidiaceae","g__Botryobasidium","s__conspersum"
 "OTU.174","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Trechisporales","f__Hydnodontaceae","g__Subulicystidium","s__perlongisporum"
 "OTU.156","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales",NA,NA,NA
-"OTU.63","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Pleosporales","f__Phaeosphaeriaceae","g__Phaeosphaeria","s__tofieldiae"
+"OTU.102","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Polyporales","f__Steccherinaceae","g__Steccherinum",NA
 "OTU.81","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Pleosporales","f__Pleosporaceae","g__Alternaria","s__dactylidicola"
 "OTU.263","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Hymenochaetales","f__Rickenellaceae","g__Sidera","s__lenis"
 "OTU.142","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Atheliales","f__Atheliaceae","g__Athelia","s__epiphylla"
@@ -95,15 +95,15 @@
 "OTU.287","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Venturiales","f__Venturiaceae","g__Tothia","s__fuscella"
 "OTU.138","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales","f__Cordycipitaceae","g__Lecanicillium","s__psalliotae"
 "OTU.353","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Trichocomaceae","g__Talaromyces","s__bacillisporus"
-"OTU.298","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Agaricales",NA,NA,NA
+"OTU.298","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Agaricales","f__Tricholomataceae","g__Mycena",NA
 "OTU.194","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Russulales","f__Stereaceae","g__Stereum","s__sanguinolentum"
 "OTU.230","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Agaricales","f__Strophariaceae","g__Deconica","s__citrispora"
 "OTU.179","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Agaricales","f__Tricholomataceae","g__Panellus","s__mitis"
 "Trichoderma","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales","f__Hypocreaceae","g__Trichoderma",NA
 "OTU.107","k__Fungi","p__Basidiomycota","c__Malasseziomycetes","o__Malasseziales",NA,NA,NA
 "OTU.342","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Tremellales","f__Bulleribasidiaceae","g__Dioszegia","s__fristingensis"
 "OTU.262","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Tremellales","f__Bulleribasidiaceae","g__Dioszegia","s__hungarica"
-"OTU.177","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Cladosporiaceae","g__Rachicladosporium","s__antarcticum"
+"OTU.177","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Cladosporiaceae","g__Rachicladosporium",NA
 "OTU.225","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Diaporthales","f__Gnomoniaceae","g__Gnomoniopsis","s__fructicola"
 "OTU.233","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Cantharellales","f__Botryobasidiaceae","g__Botryobasidium",NA
 "OTU.318","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Russulales","f__Stereaceae",NA,NA
@@ -161,7 +161,7 @@
 "OTU.47","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Filobasidiales","f__Filobasidiaceae","g__Filobasidium","s__stepposum"
 "OTU.113","k__Fungi","p__Basidiomycota","c__Microbotryomycetes","o__Sporidiobolales","f__Sporidiobolaceae","g__Sporobolomyces","s__roseus"
 "OTU.96","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Aspergillus",NA
-"OTU.195","k__Fungi","p__Basidiomycota","c__Malasseziomycetes","o__Malasseziales","f__Malasseziaceae","g__Malassezia",NA
+"OTU.195","k__Fungi","p__Basidiomycota","c__Malasseziomycetes","o__Malasseziales","f__Malasseziaceae","g__Malassezia","s__restricta"
 "OTU.58","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Mycosphaerellaceae","g__Ramularia","s__mali"
 "Dioszegia","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Tremellales","f__Bulleribasidiaceae","g__Dioszegia","s__butyracea"
 "OTU.169","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Hymenochaetales","f__Schizoporaceae","g__Hyphodontia","s__pallidula"
@@ -176,7 +176,7 @@
 "OTU.123","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Polyporales","f__Meruliaceae","g__Phlebia","s__rufa"
 "OTU.244","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Russulales","f__Peniophoraceae","g__Peniophora",NA
 "OTU.305","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Russulales","f__Peniophoraceae","g__Peniophora",NA
-"OTU.106","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Aspergillus","s__montevidensis"
+"OTU.106","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Aspergillus",NA
 "OTU.250","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Russulales","f__Peniophoraceae","g__Peniophora",NA
 "OTU.77","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Hymenochaetales","f__Tubulicrinaceae","g__Tubulicrinis","s__glebulosus"
 "OTU.277","k__Fungi","p__Basidiomycota","c__Tremellomycetes","o__Cystofilobasidiales","f__Mrakiaceae","g__Itersonilia","s__pannonica"
@@ -192,7 +192,7 @@
 "OTU.265","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Penicillium",NA
 "OTU.66","k__Fungi","p__Ascomycota","c__Eurotiomycetes","o__Eurotiales","f__Aspergillaceae","g__Penicillium","s__simplicissimum"
 "OTU.228","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales",NA,NA,NA
-"OTU.46","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales",NA,NA,NA
+"OTU.46","k__Fungi","p__Ascomycota","c__Sordariomycetes","o__Hypocreales","f__Hypocreaceae","g__Monocillium",NA
 "OTU.334","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Auriculariales","f__Exidiaceae","g__Exidia",NA
 "OTU.279","k__Fungi","p__Basidiomycota","c__Agaricomycetes","o__Auriculariales","f__Exidiaceae","g__Exidia","s__glandulosa"
 "OTU.126","k__Fungi","p__Ascomycota","c__Dothideomycetes","o__Capnodiales","f__Mycosphaerellaceae",NA,NA