#graphing the environmental data #nutrients at 44, 42, 25 and HIMB from 2014-2015 library(scales); install.packages (zoo) library(zoo) reefcols <- c("#8dd3c7", "#8dbdd3", "#bebada", "#d9d9d9") # Import nutrient data nutrients <- read.csv("~/Desktop/Kbay_nutrient_data_2014_2016.csv") head (nutrients) # set date to date class nutrients$date <- as.Date(nutrients$date, format="%Y-%m-%d") head (nutrients) #reduce file to the 4 reefs rf44.nutr <- subset(nutrients, site=="reef44") rf25.nutr <- subset(nutrients, site=="reef25") rf42.nutr <- subset(nutrients, site=="reef42") HIMB.nutr <- subset(nutrients, site=="HIMB") head (rf42.nutr) head (HIMB.nutr) # plot phosphate values for each reef pdf (file="~/Desktop/Fig_phosphate.pdf", width= 7, height= 3) par(mar=c(2,3,2,1), mgp=c(2,0.5,0)) plot(phosphate ~ date, rf44.nutr, type="n", ylab="phosphate", xaxt="n", xlab="", ylim=c(0.02,.25)) lines (rf44.nutr$date, rf44.nutr$phosphate, type="o", col=reefcols[1]) lines (rf25.nutr$date, rf25.nutr$phosphate, type="o", col=reefcols[3]) lines (HIMB.nutr$date, HIMB.nutr$phosphate, type="o", col=reefcols[4]) lines (rf42.nutr$date, rf42.nutr$phosphate, type="o", col=reefcols[2]) title ("Phosphate", adj=0) #mtext(expression(bold("B")), 2, adj=4.5, las=1, padj=-8) axis.Date(1, at=rf44.nutr$date, format="%m/%Y") legend("topleft", lty=1, col=reefcols, legend=c("44","42","25","HIMB"), lwd=2, bty="n", inset=c(0.07,0)) dev.off() #plot for Silicate at all reefs pdf (file="~/Desktop/Fig_silicate.pdf", width= 7, height= 3) par(mar=c(2,3,2,1), mgp=c(2,0.5,0)) plot(silicate ~ date, rf44.nutr, type="n", ylab="Silicate", xaxt="n", xlab="", ylim=c(0,38)) lines (rf44.nutr$date, rf44.nutr$silicate, type="o", col=reefcols[1]) lines (rf25.nutr$date, rf25.nutr$silicate, type="o", col=reefcols[3]) lines (HIMB.nutr$date, HIMB.nutr$silicate, type="o", col=reefcols[4]) lines (rf42.nutr$date, rf42.nutr$silicate, type="o", col=reefcols[2]) title ("Silicate", adj=0) #mtext(expression(bold("B")), 2, adj=4.5, las=1, padj=-8) axis.Date(1, at=rf44.nutr$date, format="%m/%Y") legend("topleft", lty=1, col=reefcols, legend=c("44","42","25","HIMB"), lwd=2, bty="n", inset=c(0.07,0)) dev.off() #plot for Nitrate+Nitrite at all reefs pdf (file="~/Desktop/Fig_N_N.pdf", width= 7, height= 3) par(mar=c(2,3,2,1), mgp=c(2,0.5,0)) plot(nitrate_nitrite ~ date, rf44.nutr, type="n", ylab="Nitrate+Nitrite", xaxt="n", xlab="", ylim=c(0,3.8)) lines (rf44.nutr$date, rf44.nutr$nitrate_nitrite, type="o", col=reefcols[1]) lines (rf25.nutr$date, rf25.nutr$nitrate_nitrite, type="o", col=reefcols[3]) lines (HIMB.nutr$date, HIMB.nutr$nitrate_nitrite, type="o", col=reefcols[4]) lines (rf42.nutr$date, rf42.nutr$nitrate_nitrite, type="o", col=reefcols[2]) title ("Nitrate+Nitrite", adj=0) #mtext(expression(bold("B")), 2, adj=4.5, las=1, padj=-8) axis.Date(1, at=rf44.nutr$date, format="%m/%Y") legend("topleft", lty=1, col=reefcols, legend=c("44","42","25","HIMB"), lwd=2, bty="n", inset=c(0.07,0)) dev.off() #plot for Ammonia at all reefs pdf (file="~/Desktop/Fig_ammonia.pdf", width= 7, height= 3) par(mar=c(2,3,2,1), mgp=c(2,0.5,0)) plot(ammonia ~ date, rf44.nutr, type="n", ylab="Ammonia", xaxt="n", xlab="", ylim=c(0,5.6)) lines (rf44.nutr$date, rf44.nutr$ammonia, type="o", col=reefcols[1]) lines (rf25.nutr$date, rf25.nutr$ammonia, type="o", col=reefcols[3]) lines (HIMB.nutr$date, HIMB.nutr$ammonia, type="o", col=reefcols[4]) lines (rf42.nutr$date, rf42.nutr$ammonia, type="o", col=reefcols[2]) title ("Ammonia", adj=0) #mtext(expression(bold("B")), 2, adj=4.5, las=1, padj=-8) axis.Date(1, at=rf44.nutr$date, format="%m/%Y") legend("topleft", lty=1, col=reefcols, legend=c("44","42","25","HIMB"), lwd=2, bty="n", inset=c(0.07,0)) dev.off() #plot for total nitrogen at all reefs pdf (file="~/Desktop/Fig_tot_N.pdf", width= 7, height= 3) par(mar=c(2,3,2,1), mgp=c(2,0.5,0)) plot(total_n ~ date, rf44.nutr, type="n", ylab="Total_N", xaxt="n", xlab="", ylim=c(0,31)) lines (rf44.nutr$date, rf44.nutr$total_n, type="o", col=reefcols[1]) lines (rf25.nutr$date, rf25.nutr$total_n, type="o", col=reefcols[3]) lines (HIMB.nutr$date, HIMB.nutr$total_n, type="o", col=reefcols[4]) lines (rf42.nutr$date, rf42.nutr$total_n, type="o", col=reefcols[2]) title ("Total_N", adj=0) #mtext(expression(bold("B")), 2, adj=4.5, las=1, padj=-8) axis.Date(1, at=rf44.nutr$date, format="%m/%Y") legend("topleft", lty=1, col=reefcols, legend=c("44","42","25","HIMB"), lwd=2, bty="n", inset=c(0.07,0)) dev.off() #analysis for 1-way ANOVA in Phosphate among all the reefs fit <- aov(phosphate ~ site, data= nutrients) summary(fit) # difference, p=0.004 TukeyHSD (fit) #44 is different from 25 #analysis for 1-way ANOVA in Silica among all the reefs fit2 <- aov(silicate ~ site, data= nutrients) summary(fit2) # difference, p=0.338 #TukeyHSD (fit2) #analysis for 1-way ANOVA in N_N among all the reefs fit3 <- aov(nitrate_nitrite ~ site, data= nutrients) summary(fit3) # difference, p=0.000 TukeyHSD (fit3) #44 is different from 25 and HIMB but not 42 #analysis for 1-way ANOVA in Ammonia among all the reefs fit4 <- aov(ammonia ~ site, data= nutrients) summary(fit4) # difference, p=0.005 TukeyHSD (fit4) #44 and #42 are different from HIMB #analysis for 1-way ANOVA in Total_N among all the reefs fit5 <- aov(total_n ~ site, data= nutrients) summary(fit5) # difference, p=0.0001 TukeyHSD (fit5) #44 is different from 42, 25 and HIMB