########################### ### GAMs on r, R0 and T ### ########################### library(bbmle) library(multcomp) library(mgcv) Life_Table<-read.csv("Life Table Summary.csv") bootstraps<-read.csv("r boostraps.csv") Births<-read.csv("Births_DS.csv", na.string="NA") BP<-Births[Births$NonDeath=="No",] Birth_Pared<-BP[BP$LastCheck>0,] P_Birth<-Births[Births$Total>0,] gam_b<-gam(Total~s(Temp,k=3,bs='cr')+TRT+s(Temp,by=TRT,k=3,bs='cr'),data=P_Birth) par(mfrow=c(1,4)) gam.check(gam_b) summary.gam(gam_b) anova.gam(gam_b) gam_r<- gam(r ~ s(Temp,k=5,bs='cr')+TRT+s(Temp,by=TRT,k=3,bs='cr'),data=Life_Table) par(mfrow=c(1,4)) gam.check(gam_r) summary.gam(gam_r) anova.gam(gam_r) #par(mfrow=c(1,5)) #plot.gam(gam_r) #res.r<-residuals(gam_r) #par(mfrow=c(1,3)) #qq.gam(gam_r,rep=100) #plot(fitted(gam_r),res.r) #plot(Life_Table$Temp,res.r) gam_R0<- gam(R0 ~ s(Temp,k=3,bs='cr')+TRT+s(Temp,by=TRT,k=3,bs='cr'),data=Life_Table) par(mfrow=c(1,4)) gam.check(gam_R0) summary.gam(gam_R0) anova.gam(gam_R0) #par(mfrow=c(1,5)) #plot.gam(gam_R0) gam_T<- gam(T ~ s(Temp,k=3,bs='cr')+TRT+s(Temp,by=TRT,k=3,bs='cr'),data=Life_Table) par(mfrow=c(1,4)) gam.check(gam_T) summary.gam(gam_T) anova.gam(gam_T) #par(mfrow=c(1,5)) #plot.gam(gam_T) ########################################################### ########################################################### ############# Graphs of Temp-dependence by TRT ############ ########################################################### ########################################################### ########################################################### library(plyr) mA1<-ddply(Life_Table[Life_Table$TRT=="A",],.(Temp)) mB1<-ddply(Life_Table[Life_Table$TRT=="B",],.(Temp)) mC1<-ddply(Life_Table[Life_Table$TRT=="C",],.(Temp)) mD1<-ddply(Life_Table[Life_Table$TRT=="D",],.(Temp)) bA1<-ddply(P_Birth[P_Birth$TRT=="A",],.(Temp)) bB1<-ddply(P_Birth[P_Birth$TRT=="B",],.(Temp)) bC1<-ddply(P_Birth[P_Birth$TRT=="C",],.(Temp)) bD1<-ddply(P_Birth[P_Birth$TRT=="D",],.(Temp)) ## Create an empty dataframe for the predict function ## Do for each response variable and Treatment combination ## Smnaller increment is to allow better placement of maxima on lines predictorA_r<-data.frame(TRT="A",Temp=seq(from = 10, to = 34, by = 0.1)) predictorB_r<-data.frame(TRT="B",Temp=seq(from = 10, to = 34, by = 0.1)) predictorC_r<-data.frame(TRT="C",Temp=seq(from = 10, to = 34, by = 0.1)) predictorD_r<-data.frame(TRT="D",Temp=seq(from = 10, to = 34, by = 0.1)) predictorA_T<-data.frame(TRT="A",Temp=seq(from = 10, to = 34, by = 0.1)) predictorB_T<-data.frame(TRT="B",Temp=seq(from = 10, to = 34, by = 0.1)) predictorC_T<-data.frame(TRT="C",Temp=seq(from = 10, to = 34, by = 0.1)) predictorD_T<-data.frame(TRT="D",Temp=seq(from = 10, to = 34, by = 0.1)) predictorA_R0<-data.frame(TRT="A",Temp=seq(from = 10, to = 34, by = 0.1)) predictorB_R0<-data.frame(TRT="B",Temp=seq(from = 10, to = 34, by = 0.1)) predictorC_R0<-data.frame(TRT="C",Temp=seq(from = 10, to = 34, by = 0.1)) predictorD_R0<-data.frame(TRT="D",Temp=seq(from = 10, to = 34, by = 0.1)) predictorA_b<-data.frame(TRT="A",Temp=seq(from = 10, to = 34, by = 0.1)) predictorB_b<-data.frame(TRT="B",Temp=seq(from = 10, to = 34, by = 0.1)) predictorC_b<-data.frame(TRT="C",Temp=seq(from = 10, to = 34, by = 0.1)) predictorD_b<-data.frame(TRT="D",Temp=seq(from = 10, to = 34, by = 0.1)) ## Use best model for predict gam. # Use type="response" to get predictions on the response scale! pred_A_b<-predict.gam(gam_b,type="response",newdata = predictorA_b,se.fit=TRUE) pred_B_b<-predict.gam(gam_b,type="response",newdata = predictorB_b,se.fit=TRUE) pred_C_b<-predict.gam(gam_b,type="response",newdata = predictorC_b,se.fit=TRUE) pred_D_b<-predict.gam(gam_b,type="response",newdata = predictorD_b,se.fit=TRUE) pred_A_r<-predict.gam(gam_r,type="response",newdata = predictorA_r,se.fit=TRUE) pred_B_r<-predict.gam(gam_r,type="response",newdata = predictorB_r,se.fit=TRUE) pred_C_r<-predict.gam(gam_r,type="response",newdata = predictorC_r,se.fit=TRUE) pred_D_r<-predict.gam(gam_r,type="response",newdata = predictorD_r,se.fit=TRUE) pred_A_T<-predict.gam(gam_T,type="response",newdata = predictorA_T,se.fit=TRUE) pred_B_T<-predict.gam(gam_T,type="response",newdata = predictorB_T,se.fit=TRUE) pred_C_T<-predict.gam(gam_T,type="response",newdata = predictorC_T,se.fit=TRUE) pred_D_T<-predict.gam(gam_T,type="response",newdata = predictorD_T,se.fit=TRUE) pred_A_R0<-predict.gam(gam_R0,type="response",newdata = predictorA_R0,se.fit=TRUE) pred_B_R0<-predict.gam(gam_R0,type="response",newdata = predictorB_R0,se.fit=TRUE) pred_C_R0<-predict.gam(gam_R0,type="response",newdata = predictorC_R0,se.fit=TRUE) pred_D_R0<-predict.gam(gam_R0,type="response",newdata = predictorD_R0,se.fit=TRUE) #### Polygon info for 95% CIs fit_A_b <-predict(gam_b,se=TRUE,newdata=predictorA_b)$fit se_A_b <-predict(gam_b,se=TRUE,newdata=predictorA_b)$se.fit lcl_A_b <-fit_A_b - 1.96*se_A_b ucl_A_b <-fit_A_b + 1.96*se_A_b fit_B_b <-predict(gam_b,se=TRUE,newdata=predictorB_b)$fit se_B_b <-predict(gam_b,se=TRUE,newdata=predictorB_b)$se.fit lcl_B_b <-fit_B_b - 1.96*se_B_b ucl_B_b <-fit_B_b + 1.96*se_B_b fit_C_b <-predict(gam_b,se=TRUE,newdata=predictorC_b)$fit se_C_b <-predict(gam_b,se=TRUE,newdata=predictorC_b)$se.fit lcl_C_b <-fit_C_b - 1.96*se_C_b ucl_C_b <-fit_C_b + 1.96*se_C_b fit_D_b <-predict(gam_b,se=TRUE,newdata=predictorD_b)$fit se_D_b <-predict(gam_b,se=TRUE,newdata=predictorD_b)$se.fit lcl_D_b <-fit_D_b - 1.96*se_D_b ucl_D_b <-fit_D_b + 1.96*se_D_b fit_A_r <-predict(gam_r,se=TRUE,newdata=predictorA_r)$fit se_A_r <-predict(gam_r,se=TRUE,newdata=predictorA_r)$se.fit lcl_A_r <-fit_A_r - 1.96*se_A_r ucl_A_r <-fit_A_r + 1.96*se_A_r fit_B_r <-predict(gam_r,se=TRUE,newdata=predictorB_r)$fit se_B_r <-predict(gam_r,se=TRUE,newdata=predictorB_r)$se.fit lcl_B_r <-fit_B_r - 1.96*se_B_r ucl_B_r <-fit_B_r + 1.96*se_B_r fit_C_r <-predict(gam_r,se=TRUE,newdata=predictorC_r)$fit se_C_r <-predict(gam_r,se=TRUE,newdata=predictorC_r)$se.fit lcl_C_r <-fit_C_r - 1.96*se_C_r ucl_C_r <-fit_C_r + 1.96*se_C_r fit_D_r <-predict(gam_r,se=TRUE,newdata=predictorD_r)$fit se_D_r <-predict(gam_r,se=TRUE,newdata=predictorD_r)$se.fit lcl_D_r <-fit_D_r - 1.96*se_D_r ucl_D_r <-fit_D_r + 1.96*se_D_r fit_A_T <-predict(gam_T,se=TRUE,newdata=predictorA_T)$fit se_A_T <-predict(gam_T,se=TRUE,newdata=predictorA_T)$se.fit lcl_A_T <-fit_A_T - 1.96*se_A_T ucl_A_T <-fit_A_T + 1.96*se_A_T fit_B_T <-predict(gam_T,se=TRUE,newdata=predictorB_T)$fit se_B_T <-predict(gam_T,se=TRUE,newdata=predictorB_T)$se.fit lcl_B_T <-fit_B_T - 1.96*se_B_T ucl_B_T <-fit_B_T + 1.96*se_B_T fit_C_T <-predict(gam_T,se=TRUE,newdata=predictorC_T)$fit se_C_T <-predict(gam_T,se=TRUE,newdata=predictorC_T)$se.fit lcl_C_T <-fit_C_T - 1.96*se_C_T ucl_C_T <-fit_C_T + 1.96*se_C_T fit_D_T <-predict(gam_T,se=TRUE,newdata=predictorD_T)$fit se_D_T <-predict(gam_T,se=TRUE,newdata=predictorD_T)$se.fit lcl_D_T <-fit_D_T - 1.96*se_D_T ucl_D_T <-fit_D_T + 1.96*se_D_T fit_A_R0 <-predict(gam_R0,se=TRUE,newdata=predictorA_R0)$fit se_A_R0 <-predict(gam_R0,se=TRUE,newdata=predictorA_R0)$se.fit lcl_A_R0 <-fit_A_R0 - 1.96*se_A_R0 ucl_A_R0 <-fit_A_R0 + 1.96*se_A_R0 fit_B_R0 <-predict(gam_R0,se=TRUE,newdata=predictorB_R0)$fit se_B_R0 <-predict(gam_R0,se=TRUE,newdata=predictorB_R0)$se.fit lcl_B_R0 <-fit_B_R0 - 1.96*se_B_R0 ucl_B_R0 <-fit_B_R0 + 1.96*se_B_R0 fit_C_R0 <-predict(gam_R0,se=TRUE,newdata=predictorC_R0)$fit se_C_R0 <-predict(gam_R0,se=TRUE,newdata=predictorC_R0)$se.fit lcl_C_R0 <-fit_C_R0 - 1.96*se_C_R0 ucl_C_R0 <-fit_C_R0 + 1.96*se_C_R0 fit_D_R0 <-predict(gam_R0,se=TRUE,newdata=predictorD_R0)$fit se_D_R0 <-predict(gam_R0,se=TRUE,newdata=predictorD_R0)$se.fit lcl_D_R0 <-fit_D_R0 - 1.96*se_D_R0 ucl_D_R0 <-fit_D_R0 + 1.96*se_D_R0 # Creating polygons for shaded area (CI's) i.for_A_b <-order(predictorA_b$Temp) i.back_A_b <-order(predictorA_b$Temp,decreasing = TRUE) x.poly_A_b <-c(predictorA_b$Temp[i.for_A_b], predictorA_b$Temp[i.back_A_b]) y.poly_A_b <-c(ucl_A_b[i.for_A_b],lcl_A_b[i.back_A_b]) i.for_B_b <-order(predictorB_b$Temp) i.back_B_b <-order(predictorB_b$Temp,decreasing = TRUE) x.poly_B_b <-c(predictorB_b$Temp[i.for_B_b], predictorB_b$Temp[i.back_B_b]) y.poly_B_b <-c(ucl_B_b[i.for_B_b],lcl_B_b[i.back_B_b]) i.for_C_b <-order(predictorC_b$Temp) i.back_C_b <-order(predictorC_b$Temp,decreasing = TRUE) x.poly_C_b <-c(predictorC_b$Temp[i.for_C_b], predictorC_b$Temp[i.back_C_b]) y.poly_C_b <-c(ucl_C_b[i.for_C_b],lcl_C_b[i.back_C_b]) i.for_D_b <-order(predictorD_b$Temp) i.back_D_b <-order(predictorD_b$Temp,decreasing = TRUE) x.poly_D_b <-c(predictorD_b$Temp[i.for_D_b], predictorD_b$Temp[i.back_D_b]) y.poly_D_b <-c(ucl_D_b[i.for_D_b],lcl_D_b[i.back_D_b]) i.for_A_r <-order(predictorA_r$Temp) i.back_A_r <-order(predictorA_r$Temp,decreasing = TRUE) x.poly_A_r <-c(predictorA_r$Temp[i.for_A_r], predictorA_r$Temp[i.back_A_r]) y.poly_A_r <-c(ucl_A_r[i.for_A_r],lcl_A_r[i.back_A_r]) i.for_B_r <-order(predictorB_r$Temp) i.back_B_r <-order(predictorB_r$Temp,decreasing = TRUE) x.poly_B_r <-c(predictorB_r$Temp[i.for_B_r], predictorB_r$Temp[i.back_B_r]) y.poly_B_r <-c(ucl_B_r[i.for_B_r],lcl_B_r[i.back_B_r]) i.for_C_r <-order(predictorC_r$Temp) i.back_C_r <-order(predictorC_r$Temp,decreasing = TRUE) x.poly_C_r <-c(predictorC_r$Temp[i.for_C_r], predictorC_r$Temp[i.back_C_r]) y.poly_C_r <-c(ucl_C_r[i.for_C_r],lcl_C_r[i.back_C_r]) i.for_D_r <-order(predictorD_r$Temp) i.back_D_r <-order(predictorD_r$Temp,decreasing = TRUE) x.poly_D_r <-c(predictorD_r$Temp[i.for_D_r], predictorD_r$Temp[i.back_D_r]) y.poly_D_r <-c(ucl_D_r[i.for_D_r],lcl_D_r[i.back_D_r]) i.for_A_R0 <-order(predictorA_R0$Temp) i.back_A_R0 <-order(predictorA_R0$Temp,decreasing = TRUE) x.poly_A_R0 <-c(predictorA_R0$Temp[i.for_A_R0], predictorA_R0$Temp[i.back_A_R0]) y.poly_A_R0 <-c(ucl_A_R0[i.for_A_R0],lcl_A_R0[i.back_A_R0]) i.for_B_R0 <-order(predictorB_R0$Temp) i.back_B_R0 <-order(predictorB_R0$Temp,decreasing = TRUE) x.poly_B_R0 <-c(predictorB_R0$Temp[i.for_B_R0], predictorB_R0$Temp[i.back_B_R0]) y.poly_B_R0 <-c(ucl_B_R0[i.for_B_R0],lcl_B_R0[i.back_B_R0]) i.for_C_R0 <-order(predictorC_R0$Temp) i.back_C_R0 <-order(predictorC_R0$Temp,decreasing = TRUE) x.poly_C_R0 <-c(predictorC_R0$Temp[i.for_C_R0], predictorC_R0$Temp[i.back_C_R0]) y.poly_C_R0 <-c(ucl_C_R0[i.for_C_R0],lcl_C_R0[i.back_C_R0]) i.for_D_R0 <-order(predictorD_R0$Temp) i.back_D_R0 <-order(predictorD_R0$Temp,decreasing = TRUE) x.poly_D_R0 <-c(predictorD_R0$Temp[i.for_D_R0], predictorD_R0$Temp[i.back_D_R0]) y.poly_D_R0 <-c(ucl_D_R0[i.for_D_R0],lcl_D_R0[i.back_D_R0]) i.for_A_T <-order(predictorA_T$Temp) i.back_A_T <-order(predictorA_T$Temp,decreasing = TRUE) x.poly_A_T <-c(predictorA_T$Temp[i.for_A_T], predictorA_T$Temp[i.back_A_T]) y.poly_A_T <-c(ucl_A_T[i.for_A_T],lcl_A_T[i.back_A_T]) i.for_B_T <-order(predictorB_T$Temp) i.back_B_T <-order(predictorB_T$Temp,decreasing = TRUE) x.poly_B_T <-c(predictorB_T$Temp[i.for_B_T], predictorB_T$Temp[i.back_B_T]) y.poly_B_T <-c(ucl_B_T[i.for_B_T],lcl_B_T[i.back_B_T]) i.for_C_T <-order(predictorC_T$Temp) i.back_C_T <-order(predictorC_T$Temp,decreasing = TRUE) x.poly_C_T <-c(predictorC_T$Temp[i.for_C_T], predictorC_T$Temp[i.back_C_T]) y.poly_C_T <-c(ucl_C_T[i.for_C_T],lcl_C_T[i.back_C_T]) i.for_D_T <-order(predictorD_T$Temp) i.back_D_T <-order(predictorD_T$Temp,decreasing = TRUE) x.poly_D_T <-c(predictorD_T$Temp[i.for_D_T], predictorD_T$Temp[i.back_D_T]) y.poly_D_T <-c(ucl_D_T[i.for_D_T],lcl_D_T[i.back_D_T]) ####################################################################### ####################################################################### ####################################################################### ####################### ################################ ####################### 4-Panel Graphs ################################ ####################### ################################ ####################################################################### ####################################################################### ##################################### ############ R0 #################### ##################################### # Columns are Treatments tiff(filename = "Fig 3 Net Rep Rate.tiff",width = 8,height = 4,units = "in",res = 600) ### R0 ### ymin=0 ymax=32 xmin=10 xmax=34 par(mfrow=c(1,4)) par(cex = 0.6) par(mar=c(0,0,0,0),oma=c(4,4,0.5,0.5)) par(tcl=-0.25) par(mgp=c(2,0.6,0)) ## (Overlay Control (black) within column for reference ## Column 1-Clutch Size plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_A_R0,y.poly_A_R0,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mA1$Temp,mA1$R0, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("No Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorA_R0$Temp,fit_A_R0,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_B_R0,y.poly_B_R0,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mB1$Temp,mB1$R0, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Early Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorB_R0$Temp,fit_B_R0,lwd=2) lines(predictorA_R0$Temp,fit_A_R0,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_C_R0,y.poly_C_R0,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mC1$Temp,mC1$R0, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Late Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorC_R0$Temp,fit_C_R0,lwd=2) lines(predictorA_R0$Temp,fit_A_R0,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_D_R0,y.poly_D_R0,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mD1$Temp,mD1$R0, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Constant Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorD_R0$Temp,fit_D_R0,lwd=2) lines(predictorA_R0$Temp,fit_A_R0,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) mtext("Temperature (C)",side = 1, outer=TRUE,cex=1,line = 2.2,col="black") mtext("Net Reproductive Rate (R0)",side = 2, outer=TRUE,cex=1,line = 2.2,col="black") dev.off() ##################################### ############ T #################### ##################################### # Columns are Treatments tiff(filename = "Fig 2 Generation Time.tiff",width = 8,height = 4,units = "in",res = 600) ### T ### ymin=5 ymax=25 xmin=10 xmax=34 par(mfrow=c(1,4)) par(cex = 0.6) par(mar=c(0,0,0,0),oma=c(4,4,0.5,0.5)) par(tcl=-0.25) par(mgp=c(2,0.6,0)) ## (Overlay Control (black) within column for reference ## Column 1-Clutch Size plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_A_T,y.poly_A_T,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mA1$Temp,mA1$T, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',col="black") mtext("No Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorA_T$Temp,fit_A_T,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_B_T,y.poly_B_T,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mB1$Temp,mB1$T, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Early Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorB_T$Temp,fit_B_T,lwd=2) lines(predictorA_T$Temp,fit_A_T,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_C_T,y.poly_C_T,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mC1$Temp,mC1$T, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Late Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorC_T$Temp,fit_C_T,lwd=2) lines(predictorA_T$Temp,fit_A_T,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_D_T,y.poly_D_T,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mD1$Temp,mD1$T, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Constant Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorD_T$Temp,fit_D_T,lwd=2) lines(predictorA_T$Temp,fit_A_T,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) mtext("Temperature (C)",side = 1, outer=TRUE,cex=1,line = 2.2,col="black") mtext("Generation Time (T)",side = 2, outer=TRUE,cex=1,line = 2.2,col="black") dev.off() ##################################### ############ r #################### ##################################### # Columns are Treatments tiff(filename = "Fig 1 r.tiff",width = 8,height = 4,units = "in",res = 600) ### r ### ymin=0 ymax=0.25 xmin=10 xmax=34 par(mfrow=c(1,4)) par(cex = 0.6) par(mar=c(0,0,0,0),oma=c(4,4,0.5,0.5)) par(tcl=-0.25) par(mgp=c(2,0.6,0)) par(mfrow=c(1,4)) ## (Overlay Control (black) within column for reference ## Column 1-Clutch Size plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_A_r,y.poly_A_r,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mA1$Temp,mA1$r, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',col="black") mtext("No Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorA_r$Temp,fit_A_r,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_B_r,y.poly_B_r,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mB1$Temp,mB1$r, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Early Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorB_r$Temp,fit_B_r,lwd=2) lines(predictorA_r$Temp,fit_A_r,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_C_r,y.poly_C_r,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mC1$Temp,mC1$r, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Late Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorC_r$Temp,fit_C_r,lwd=2) lines(predictorA_r$Temp,fit_A_r,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) plot(0,typ="n",bty="n",xlab="", ylab="", xaxt='n',yaxt='n', xlim=c(xmin,xmax),ylim=c(ymin,ymax)) polygon(x.poly_D_r,y.poly_D_r,col=rgb(0.5,0.5,0.5,0.5),border="black") par(new=TRUE) plot(mD1$Temp,mD1$r, xlim=c(xmin,xmax),ylim=c(ymin,ymax), xlab='',ylab='',main='', xaxt='n',yaxt='n',col="black") mtext("Constant Cue",side=3,line=-2,adj=0.5,cex=1,col="black") par(new=TRUE) lines(predictorD_r$Temp,fit_D_r,lwd=2) lines(predictorA_r$Temp,fit_A_r,lty=2,lwd=2) xtick<-c(11,17,23,27,29,31,33) axis(side=1,at=xtick) mtext("Temperature (C)",side = 1, outer=TRUE,cex=1,line = 2.2,col="black") mtext("Intrinsic Rate of Increase (r)",side = 2, outer=TRUE,cex=1,line = 2.2,col="black") dev.off()