############### #Variance in plant functional traits in Eucalyptus baxteri and E. obliqua – a pilot study to assess utility for inferring climate vulnerability #Bradshaw, H. & Guerin, G. (2022) #SSD stem specific density and SLA specific leaf area to indicate general vulnerability to water stress and cavitation, also leaf area and L:W ratio #generally x3 or x5 individual trees per population #x1 stem sample / x5 leaf samples per individual ### #load raw data #stem ssd_data <- read.csv("SSD_data_Figshare.csv") #leaf sla_data <- read.csv("SLA_data_Figshare.csv") ############# #variance partitioning car::Anova(lm(sla ~ Species/Location/Tree_ID, sla_data)) #anova output - sum of squares (% of variance) adding to 100% #Species: 21.5% #Locality:18.9% #Individual: 31.6% #Within-individual (residual): 28.1% car::Anova(lm(LWR ~ Species/Location/Tree_ID, sla_data)) #Species: 8.5% #Locality:18.4% #Individual: 30.1% #Within-individual (residual): 42.9% car::Anova(lm(Leaf.area_mm2 ~ Species/Location/Tree_ID, sla_data)) #Species: 1.2% #Locality: 23.3% #Individual: 35.2% #Within-individual (residual): 40.2% #### car::Anova(lm(ssd ~ Species/Location, ssd_data)) #Species: 4.1% (note very little replication for baxteri) #Location: 26.9% #Within location (between individual): 69% #################### #compare tree-level data to rainfall #mean of SLA reps per individual ind_means <- aggregate(sla_data$sla, by=list(sla_data$Tree_ID), FUN=mean) names(ind_means) <- c("Tree_ID", "SLA_mean") #add back to main data sla_data <- merge(sla_data, ind_means, by="Tree_ID") #mean of LWR reps per individual ind_means2 <- aggregate(sla_data$LWR, by=list(sla_data$Tree_ID), FUN=mean) names(ind_means2) <- c("Tree_ID", "LWR_mean") #add back to main data sla_data <- merge(sla_data, ind_means2, by="Tree_ID") #mean of leaf area reps per individual ind_means3 <- aggregate(sla_data$Leaf.area_mm2, by=list(sla_data$Tree_ID), FUN=mean) names(ind_means3) <- c("Tree_ID", "Leaf.area_mean") #add back to main data sla_data <- merge(sla_data, ind_means3, by="Tree_ID") ### pdf("Trat_rainfall_lm.pdf", height=6, width=8) par(mfrow=c(2,2)) par(mar=c(5,3.6,1,0.6)) #set margins par(cex.axis=0.8) #set axis label text size par(mgp=c(2,0.7,0)) #position axis labels plot(SLA_mean ~ Rainfall.annual..mean.ml., subset(sla_data, Species == "obliqua"), ylim=c(3,8), xlim=c(600,1200), col="blue", pch=20, cex=1.2, xlab="Mean annual rainfall (mm)", ylab="SLA (mm2/mg)") points(subset(sla_data, Species == "baxteri")$Rainfall.annual..mean.ml., subset(sla_data, Species == "baxteri")$SLA_mean, col="red", pch=1, cex=1.2) abline(lm(SLA_mean ~ Rainfall.annual..mean.ml., sla_data[-which(duplicated(sla_data$SLA_mean)),])) legend("topleft", pch=c(1,20), col=c("red", "blue"), legend=c("baxteri", "obliqua"), bty="n", cex=1.2) # plot(LWR_mean ~ Rainfall.annual..mean.ml., subset(sla_data, Species == "obliqua"), ylim=c(2,6), xlim=c(600,1200), col="blue", pch=20, cex=1.2, xlab="Mean annual rainfall (mm)", ylab="L:W ratio") points(subset(sla_data, Species == "baxteri")$Rainfall.annual..mean.ml., subset(sla_data, Species == "baxteri")$LWR_mean, col="red", pch=1, cex=1.2) abline(lm(LWR_mean ~ Rainfall.annual..mean.ml., sla_data[-which(duplicated(sla_data$SLA_mean)),]), lty=2) # plot(Leaf.area_mean ~ Rainfall.annual..mean.ml., subset(sla_data, Species == "obliqua"), ylim=c(1000,6000), xlim=c(600,1200), col="blue", pch=20, cex=1.2, xlab="Mean annual rainfall (mm)", ylab="Leaf area (mm2)") points(subset(sla_data, Species == "baxteri")$Rainfall.annual..mean.ml., subset(sla_data, Species == "baxteri")$Leaf.area_mean, col="red", pch=1, cex=1.2) abline(lm(Leaf.area_mean ~ Rainfall.annual..mean.ml., sla_data[-which(duplicated(sla_data$SLA_mean)),])) # plot(ssd ~ Rainfall.annual..mean.ml., subset(ssd_data, Species == "obliqua"), ylim=c(0.3,1), xlim=c(600,1200), col="blue", pch=20, cex=1.2, xlab="Mean annual rainfall (mm)", ylab="Stem specific density (mg/mm3)") points(subset(ssd_data, Species == "baxteri")$Rainfall.annual..mean.ml., subset(ssd_data, Species == "baxteri")$ssd, col="red", pch=1, cex=1.2) abline(lm(ssd ~ Rainfall.annual..mean.ml., ssd_data), lty=2) #slight decreasing trend, non-significant dev.off() ### #################### library(RColorBrewer) Set3 <- brewer.pal(10, "Set3") #colour vector my.order <- aggregate(sla_data, by=list(sla_data$Species), FUN=function(x) {length(unique(x))})$Tree_ID #vector of number of trees sampled per species (for arranging colours) ############################# #arrange sites by sla or rainfall for plotting, e.g. new_order <- with(sla_data, reorder(Location, sla, mean, na.rm=T)) # new_order2 <- with(sla_data, reorder(Location, MAP_mm, mean , na.rm=T)) # new_order4 <- with(subset(sla_data, Species=="obliqua"), reorder(Location, MAP_mm, mean, na.rm=T)) # new_order3 <- sla_data$Tree_ID[order(sla_data$Species)] levels(new_order3) <- unique(new_order3) # new_order5 <- with(sla_data, reorder(site_spp, MAP_mm, mean, na.rm=T)) # new_order6 <- with(ssd_data, reorder(site_spp, MAP_mm, mean , na.rm=T)) ######################## #plots: ###### #individuals boxplot (sep by species) pdf("Individuals_traits_boxplots.pdf", width=5, height=12) #write to pdf par(mfrow=c(3,1)) par(mar=c(8,3.6,1,0.6)) par(cex.axis=0.8) par(mgp=c(2.5,0.7,0)) #sla boxplot(sla ~ new_order3, sla_data, col=rep(Set3[1:2], my.order), las=2, cex.lab=0.9, pch=20, xlab="", ylab="SLA (mm2/mg)") #lwr boxplot(LWR ~ new_order3, sla_data, col=rep(Set3[1:2], my.order), las=2, cex.lab=0.9, pch=20, xlab="", ylab="L:W ratio") #area boxplot(Leaf.area_mm2 ~ new_order3, sla_data, col=rep(Set3[1:2], my.order), las=2, cex.lab=0.9, pch=20, xlab="", ylab="Leaf area (mm2)") # dev.off() ## #sites boxoplots pdf("Sites_traits_boxplots,pdf", width=6, height=7) par(mfrow=c(2,2)) par(mar=c(7,3.6,1,0.6)) par(cex.axis=0.7) par(mgp=c(2.5,0.7,0)) #order by rainfall boxplot(sla ~ new_order5, sla_data, las=2, cex.lab=0.9, pch=20, xlab="", col=brewer.pal(8, "Set3")[c(1,2,3,3,4,5,5,6,7,8)], ylab="SLA (mm2/mg)") boxplot(LWR ~ new_order5, sla_data, las=2, cex.lab=0.9, pch=20, xlab="", col=brewer.pal(8, "Set3")[c(1,2,3,3,4,5,5,6,7,8)], ylab="L:W ratio") boxplot(Leaf.area_mm2 ~ new_order5, sla_data, las=2, cex.lab=0.9, pch=20, xlab="", col=brewer.pal(8, "Set3")[c(1,2,3,3,4,5,5,6,7,8)], ylab="Leaf area (mm2)") boxplot(ssd ~ new_order6, ssd_data, las=2, cex.lab=0.9, pch=20, xlab="", col=brewer.pal(8, "Set3")[c(1,2,3,3,4,5,5,6,7,8)], ylab="Stem specific density (mg/mm3)") dev.off() #species boxplot pdf("Species_trait_boxplots.pdf", width=6, height=6) par(mar=c(4,4,1,1)) par(cex.axis=0.9) par(mgp=c(3,1,0)) par(mfrow=c(2,2)) boxplot(sla ~ Species, sla_data, las=2, cex.lab=0.9, col=Set3, pch = 20, xlab="", ylab="SLA (mm2/mg)") boxplot(LWR ~ Species, sla_data, las=2, cex.lab=0.9, col=Set3, pch = 20, xlab="", ylab="L:W ratio") boxplot(Leaf.area_mm2 ~ Species, sla_data, las=2, cex.lab=0.9, col=Set3, pch = 20, xlab="", ylab="Leaf area (mm2)") boxplot(ssd ~ Species, ssd_data, las=2, cex.lab=0.9, col=Set3, pch = 20, xlab="", ylab="Stem specific density (mg/mm3)") # dev.off()