#Nelson et al.—American Journal of Botany 2020—Appendix S1 # Statistical analyses and figures for Nelson et al. "Testing an invasion mechanism for Eucalyptus globulus: is there evidence of allelopathy?" library(car) #for Anova function library(emmeans) #for posthoc comparisons library(multcomp) #for independent contrasts library(sciplot) #for barcharts #### set working directory setwd("~/Dropbox/euciedata/Dryad") ################################# # stats: water soluble experiment ################################# df.water <- read.csv("Euc_watersoluble.csv", header = T) df.water$Treatment <- factor(df.water$Treatment, levels=c('Control', 'Oak', 'Euc', 'WhSage')) mod1 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ Treatment*Species, data = df.water, family = binomial(logit)) Anova(mod1) #post hoc test of main effects (drop interaction term since it wasnt significant) mod2 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ Treatment+Species, data = df.water, family = binomial(logit)) summary(glht(mod2, mcp(Treatment="Tukey"))) ############################# # stats: volatiles experiment ############################# df.vol <- read.csv("Euc_volatiles.csv", header = T) df.vol$Treatment <- factor(df.vol$Treatment, levels=c('Control', 'Oak', 'Euc', 'WhSage')) mod3 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ Treatment*Species, data = df.vol, family = binomial(logit)) Anova(mod3) summary(glht(mod3, mcp(Treatment="Tukey"))) #preplanned contrasts df.vol$trt.sp <- interaction(df.vol$Treatment, df.vol$Species) mod4 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ trt.sp, family=binomial(logit), data = df.vol) AG_Con_Oak <- c(-1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) AG_Con_Euc <- c(-1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) AG_Con_WSa <- c(-1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) AG_Oak_Euc <- c(0,-1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) AG_Oak_WSa <- c(0,-1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) AG_Euc_WSa <- c(0,0,-1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) EC_Con_Oak <- c(0,0,0,0,-1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0) EC_Con_Euc <- c(0,0,0,0,-1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0) EC_Con_WSa <- c(0,0,0,0,-1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0) EC_Oak_Euc <- c(0,0,0,0,0,-1,1,0,0,0,0,0,0,0,0,0,0,0,0,0) EC_Oak_WSa <- c(0,0,0,0,0,-1,0,1,0,0,0,0,0,0,0,0,0,0,0,0) EC_Euc_WSa <- c(0,0,0,0,0,0,-1,1,0,0,0,0,0,0,0,0,0,0,0,0) FM_Con_Oak <- c(0,0,0,0,0,0,0,0,-1,1,0,0,0,0,0,0,0,0,0,0) FM_Con_Euc <- c(0,0,0,0,0,0,0,0,-1,0,1,0,0,0,0,0,0,0,0,0) FM_Con_WSa <- c(0,0,0,0,0,0,0,0,-1,0,0,1,0,0,0,0,0,0,0,0) FM_Oak_Euc <- c(0,0,0,0,0,0,0,0,0,-1,1,0,0,0,0,0,0,0,0,0) FM_Oak_WSa <- c(0,0,0,0,0,0,0,0,0,-1,0,1,0,0,0,0,0,0,0,0) FM_Euc_WSa <- c(0,0,0,0,0,0,0,0,0,0,-1,1,0,0,0,0,0,0,0,0) LS_Con_Oak <- c(0,0,0,0,0,0,0,0,0,0,0,0,-1,1,0,0,0,0,0,0) LS_Con_Euc <- c(0,0,0,0,0,0,0,0,0,0,0,0,-1,0,1,0,0,0,0,0) LS_Con_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,-1,0,0,1,0,0,0,0) LS_Oak_Euc <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,-1,1,0,0,0,0,0) LS_Oak_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,-1,0,1,0,0,0,0) LS_Euc_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,1,0,0,0,0) SP_Con_Oak <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,1,0,0) SP_Con_Euc <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,0,1,0) SP_Con_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,0,0,1) SP_Oak_Euc <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,1,0) SP_Oak_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,0,1) SP_Euc_WSa <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-1,1) cntrMat <- rbind(AG_Con_Oak,AG_Con_Euc,AG_Con_WSa,AG_Oak_Euc,AG_Oak_WSa,AG_Euc_WSa, EC_Con_Oak,EC_Con_Euc,EC_Con_WSa,EC_Oak_Euc,EC_Oak_WSa,EC_Euc_WSa, FM_Con_Oak,FM_Con_Euc,FM_Con_WSa,FM_Oak_Euc,FM_Oak_WSa,FM_Euc_WSa, LS_Con_Oak,LS_Con_Euc,LS_Con_WSa,LS_Oak_Euc,LS_Oak_WSa,LS_Euc_WSa, SP_Con_Oak,SP_Con_Euc,SP_Con_WSa,SP_Oak_Euc,SP_Oak_WSa,SP_Euc_WSa) summary(glht(mod4, linfct=mcp(trt.sp=cntrMat), alternative="two.sided"), test=adjusted("BH")) ################################### # stats: soil experiment germination ################################### df.soil.germ <- read.csv("Euc_soil_germination.csv", header = T) df.soil.germ$Soil <- factor(df.soil.germ$Soil, levels=c('Scrub', 'Euc')) mod5 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ Soil*Species*Season + Site, data = df.soil.germ, family = binomial(logit)) Anova(mod5) #preplanned contrasts df.soil.germ$soil.sp.season <- interaction(df.soil.germ$Soil, df.soil.germ$Species,df.soil.germ$Season) levels(df.soil.germ$soil.sp.season) mod6 <- glm(cbind(No_Germ, Seeds-No_Germ) ~ soil.sp.season, family=binomial(logit), data = df.soil.germ) AG_Scrub_Euc <- c(0,0,0,0,0,0,-1,1,0,0,0,0,0,0) EC_Scrub_Euc <- c(0,0,0,0,0,0,0,0,-1,1,0,0,0,0) FM_Scrub_Euc <- c(-1,1,0,0,0,0,0,0,0,0,0,0,0,0) LS_Scrub_Euc <- c(0,0,-1,1,0,0,0,0,0,0,-1,1,0,0) SP_Scrub_Euc <- c(0,0,0,0,-1,1,0,0,0,0,0,0,-1,1) cntrMat <- rbind(AG_Scrub_Euc,EC_Scrub_Euc,FM_Scrub_Euc,LS_Scrub_Euc,SP_Scrub_Euc) summary(glht(mod6, linfct=mcp(soil.sp.season=cntrMat), alternative="two.sided"), test=adjusted("BH")) ################################ # stats: soil experiment biomass ################################ df.soil.biomass <- read.csv("Euc_soil_biomass.csv", header = T) df.soil.biomass$Soil <- factor(df.soil.biomass$Soil, levels=c('Scrub', 'Euc')) mod7 <- glm(df.soil.biomass$Biomass_g ~ Soil*Species + Site, family = Gamma(link = "log"),data=df.soil.biomass) Anova(mod7) #preplanned contrasts df.soil.biomass$soil.sp.site <- interaction(df.soil.biomass$Soil, df.soil.biomass$Species, df.soil.biomass$Site) levels(df.soil.biomass$soil.sp.site) mod8 <- aov(log(df.soil.biomass$Biomass_g) ~ soil.sp.site, data=df.soil.biomass) FM_Scrub_Euc <- c(-1,1,0,0,0,0,-1,1,0,0,0,0) LS_Scrub_Euc <- c(0,0,-1,1,0,0,0,0,-1,1,0,0) SP_Scrub_Euc <- c(0,0,0,0,-1,1,0,0,0,0,-1,1) cntrMat <- rbind(FM_Scrub_Euc,LS_Scrub_Euc,SP_Scrub_Euc) summary(glht(mod8, linfct=mcp(soil.sp.site=cntrMat), alternative="two.sided"), test=adjusted("BH")) ############################## ######### Figure 1 ########### ############################## pdf(file="Fig1.pdf",width=9.5, height= 8.5) par(mfrow=c(2,1), mar=c(3,4,2,4)) #Panel A. Volatile compounds # make "pooled' species for plotting main effect of treatment df.vol.main <- df.vol df.vol.main$Species <- "pooled" df.vol.fig <- rbind(df.vol, df.vol.main) #order treatment levels df.vol.fig$Species <- factor(df.vol.fig$Species, levels=c('pooled','AG', 'EC', "FM", "LS", "SP")) #boxplots of treatment by species x <- boxplot(df.vol.fig$No_Germ/df.vol.fig$Seeds ~ df.vol.fig$Treatment:df.vol.fig$Species, xlab="", ylab="prop germination", border="gray20", col=c("gray31","gray47","gray77","gray88"), at=c(1,2,3,4,6,7,8,9,11,12,13,14,16,17,18,19,21,22,23,24,26,27,28,29), frame=F, xaxt="n", ylim=c(0,1.13), outline=F) df.vol.fig$sptrt <- df.vol.fig$Species:df.vol.fig$Treatment df.vol.fig$value <- df.vol.fig$No_Germ/df.vol.fig$Seeds # Add data points mylevels <- levels(df.vol.fig$sptrt) mylevels <- c(mylevels[1:4],"",mylevels[5:8],"",mylevels[9:12],"",mylevels[13:16],"",mylevels[17:20],"", mylevels[21:24]) levelProportions <- summary(df.vol.fig$sptrt)/nrow(df.vol.fig) for(i in 1:length(mylevels)){ thislevel <- mylevels[i] thisvalues <- df.vol.fig[df.vol.fig$sptrt==thislevel, "value"] # take the x-axis indices and add a jitter, proportional to the N in each level myjitter <- jitter(rep(i, length(thisvalues)), amount=0.15) points(myjitter, thisvalues, pch=20, col="black", cex=0.3) } axis(1,at=c(0.5,4.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(5.5,9.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(10.5,14.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(15.5,19.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(20.5,24.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(25.5,29.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1, at=c(2.5), labels=c("pooled"), tick=F, font=1, cex.axis=0.8) axis(1, at=c(7.5,12.5,17.5,22.5,27.5), labels=c("Acmispon glaber","Eschscholzia californica","Festuca microstachys","Lupinus succulentus","Stipa pulchra"), tick=F, font=3, cex.axis=0.65) #add legend legend(x=24, y=1.2, legend=c("air", "Quercus agrifolia", "Eucalyptus globulus", "Salvia apiana"), pt.bg=c("gray31","gray47","gray77","gray88"), pch=22, bty="n", cex=0.7, pt.cex=1.5, y.intersp=1, x.intersp=0.9, text.font=c(1,3,3,3)) #add preplanned contrast results for each species myy <- x$stats[5,1:24]+0.05 myy[14] <- 1.07 myy[21] <- 0.77 text( x=c(1,2,3,4,6,7,8,9,11,12,13,14,16,17,18,19,21,22,23,24,26,27,28,29),y=myy, cex=0.8,labels=c( "","","","", "A","A","A","A", "A","A,B","A","B", "A","A","A","B", "A","A","A","A", "A","A","A","B") ) mtext("A. Volatile compounds",side=3, line=0,adj=0) #Panel B. Water soluble compounds # make "pooled' species for plotting main effect of treatment df.water.main <- df.water df.water.main$Species <- "pooled" df.water.fig <- rbind(df.water, df.water.main) #order treatment levels df.water.fig$Species <- factor(df.water.fig$Species, levels=c('pooled','AG', 'EC', "FM", "LS", "SP")) # bargraphs of treatment by species #boxplots of treatment by species x <- boxplot(df.water.fig$No_Germ/df.water.fig$Seeds ~ df.water.fig$Treatment:df.water.fig$Species, xlab="", ylab="prop germination", border="gray20", col=c("gray31","gray47","gray77","gray88"), at=c(1,2,3,4,6,7,8,9,11,12,13,14,16,17,18,19,21,22,23,24,26,27,28,29), frame=F, xaxt="n", ylim=c(0,1.12), outline=F) df.water.fig$sptrt <- df.water.fig$Species:df.water.fig$Treatment df.water.fig$value <- df.water.fig$No_Germ/df.water.fig$Seeds # Add data points mylevels <- levels(df.water.fig$sptrt) mylevels <- c(mylevels[1:4],"",mylevels[5:8],"",mylevels[9:12],"",mylevels[13:16],"",mylevels[17:20],"", mylevels[21:24]) for(i in 1:length(mylevels)){ thislevel <- mylevels[i] thisvalues <- df.water.fig[df.water.fig$sptrt==thislevel, "value"] # take the x-axis indices and add a jitter, proportional to the N in each level myjitter <- jitter(rep(i, length(thisvalues)), amount=0.25) points(myjitter, thisvalues, pch=20, col="black", cex=0.3) } axis(1,at=c(0.5,4.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(5.5,9.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(10.5,14.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(15.5,19.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(20.5,24.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(25.5,29.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1, at=c(2.5), labels=c("pooled"), tick=F, font=1, cex.axis=0.8) axis(1, at=c(7.5,12.5,17.5,22.5,27.5), labels=c("Acmispon glaber","Eschscholzia californica","Festuca microstachys","Lupinus succulentus","Stipa pulchra"), tick=F, font=3, cex.axis=0.65) #add legend legend(x=24, y=1.2, legend=c("water", "Quercus agrifolia", "Eucalyptus globulus", "Salvia apiana"), pt.bg=c("gray31","gray47","gray77","gray88"), pch=22, bty="n", cex=0.7, pt.cex=1.5, y.intersp=1, x.intersp=0.9, text.font=c(1,3,3,3)) #add treatment comparisons text( x=c(1,2,3,4),y=1.05, labels=c("A","A,B","B","C") , cex=0.8) mtext("B. Water-soluble compounds",side=3, line=0,adj=0) dev.off() ############################## ######### Figure 2 ########### ############################## pdf(file="Fig2.pdf",width=7, height= 8.5) par(mfrow=c(2,1), mar=c(3,4,2,4)) #SOIL GERMINATION PANEL # make "pooled' species for plotting main effect of treatment df.soil.main <- df.soil.germ df.soil.main$Species <- "pooled" df.soil.fig <- rbind(df.soil.germ, df.soil.main) #order treatment levels df.soil.fig$Species <- factor(df.soil.fig$Species, levels=c('pooled','AG', 'EC', "FM", "LS", "SP")) #boxplots of treatment by species x <- boxplot(df.soil.fig$No_Germ/df.soil.fig$Seeds ~ df.soil.fig$Soil:df.soil.fig$Species, xlab="", ylab="prop germination", border="gray20", col=c("gray47","gray88"), at=c(1,2,4,5,7,8,10,11,13,14,16,17), frame=F, xaxt="n", ylim=c(0,1.15), outline=F) #add preplanned contrast results myy <- x$stats[5,1:12]+0.04 text(x=c(1,2,4,5,7,8,10,11,13,14,16,17), y=myy,cex=0.8,labels=c( "","", "A","A", "A","A", "A","B", "A","A", "A","A")) mtext("A.",side=3, line=-1,adj=0) axis(1,at=c(0.5,2.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(3.5,5.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(6.5,8.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(9.5,11.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(12.5,14.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(15.5,17.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) df.soil.fig$sptrt <- df.soil.fig$Species:df.soil.fig$Soil df.soil.fig$value <- df.soil.fig$No_Germ/df.soil.fig$Seeds # Add data points mylevels <- levels(df.soil.fig$sptrt) mylevels <- c(mylevels[1:2],"",mylevels[3:4],"",mylevels[5:6],"",mylevels[7:8],"",mylevels[9:10],"", mylevels[11:12]) for(i in 1:length(mylevels)){ thislevel <- mylevels[i] thisvalues <- df.soil.fig[df.soil.fig$sptrt==thislevel, "value"] # take the x-axis indices and add a jitter, proportional to the N in each level myjitter <- jitter(rep(i, length(thisvalues)), amount=0.25) points(myjitter, thisvalues, pch=20, col="black", cex=0.3) } axis(1, at=c(1.5), labels=c("pooled"), tick=F, font=1, cex.axis=0.8) axis(1, at=c(4.5,7.5,10.5,13.5,16.5), labels=c("A. glaber","E. californica","F. microstachys","L. succulentus","S. pulchra"), tick=F, font=3, cex.axis=0.65) #add legend legend(x=15, y=1.20, legend=c("coastal scrub", "Eucalyptus"), pt.bg=c("gray47","gray88"), pch=22, bty="n", cex=0.6, pt.cex=1.3, y.intersp=1, x.intersp=0.9, text.font=c(1,3)) #SOIL BIOMASS PANEL # make "pooled' species for plotting main effect of treatment df.soil.biomass.main <- df.soil.biomass df.soil.biomass.main$Species <- "pooled" df.soil.biomass.fig <- rbind(df.soil.biomass, df.soil.biomass.main) #order treatment levels df.soil.biomass.fig$Species <- factor(df.soil.biomass.fig$Species, levels=c('pooled','AG', 'EC', "FM", "LS", "SP")) #boxplots of treatment by species x <- boxplot(df.soil.biomass.fig$Biomass_g ~ df.soil.biomass.fig$Soil:df.soil.biomass.fig$Species, xlab="", ylab="biomass (g)", border="gray20", col=c("gray47","gray88"), at=c(1,2,4,5,7,8,10,11,13,14,16,17), frame=F, xaxt="n", outline=F, ylim=c(0,6.2)) #add preplanned contrast results myy <- x$stats[5,1:12]+0.3 text(x=c(1,2,4,5,7,8,10,11,13,14,16,17), y=myy,cex=0.8,labels=c( "","", "","", "","", "A","A", "A","B", "A","B")) mtext("B.",side=3, line=0,adj=0) text(x=c(4.5,7.5),y=0.1,labels="no data", cex=0.8) axis(1,at=c(0.5,2.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(3.5,5.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(6.5,8.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(9.5,11.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(12.5,14.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) axis(1,at=c(15.5,17.5),line=0.5,tick=T,labels=c("",""),lwd=1,lwd.ticks=0) df.soil.biomass.fig$sptrt <- df.soil.biomass.fig$Species:df.soil.biomass.fig$Soil df.soil.biomass.fig$value <- df.soil.biomass.fig$Biomass_g # Add data points mylevels <- levels(df.soil.biomass.fig$sptrt) mylevels <- c(mylevels[1:2],"","","","","","","",mylevels[7:8],"",mylevels[9:10],"", mylevels[11:12]) for(i in 1:length(mylevels)){ thislevel <- mylevels[i] thisvalues <- df.soil.biomass.fig[df.soil.biomass.fig$sptrt==thislevel, "value"] # take the x-axis indices and add a jitter, proportional to the N in each level myjitter <- jitter(rep(i, length(thisvalues)), amount=0.25) points(myjitter, thisvalues, pch=20, col="black", cex=0.3) } axis(1, at=c(1.5), labels=c("pooled"), tick=F, font=1, cex.axis=0.8) axis(1, at=c(4.5,7.5,10.5,13.5,16.5), labels=c("A. glaber","E. californica","F. microstachys","L. succulentus","S. pulchra"), tick=F, font=3, cex.axis=0.65) #add legend legend(x=15, y=6.1, legend=c("coastal scrub", "Eucalyptus"), pt.bg=c("gray47","gray88"), pch=22, bty="n", cex=0.6, pt.cex=1.3, y.intersp=1, x.intersp=0.9, text.font=c(1,3)) dev.off()