--- title: "1. Bivariate Allometry" author: "Nathalie Feiner" date: "16/01/2021" output: html_document: toc: yes toc_depth: 3 toc_float: yes pdf_document: toc: yes toc_depth: '3' --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) knitr::opts_chunk$set(message = FALSE) ``` The input dataset contains 374 individual embryos belonging to 15 Anolis species and 693 individual adults belonging to 267 species. Note that relative pathes have to be adjusted for reproducing the results. ```{r load libraries, include=F} library(ggplot2) library(openxlsx) library(pcaMethods) library(tidyr) library(dplyr) library(naniar) library(VIM) library(missMDA) library(car) library(tinytex) library(phangorn) library(geiger) library(phylobase) library(ggtree) library(phytools) library(lmerTest) library(emmeans) library(jmv) #load ontogenetic data Data <- read.csv("C:/Users/Nathalie/Dropbox/stuff/MSc_Projects/ElianeVanDerCruyssen_2020/Data_Embryos_Final.csv", row.names = 1) Data$logrelFL <- log(Data$FL)/log(Data$spine) Data$logrelHL <- log(Data$HL)/log(Data$spine) #load adult data Adult_data <- read.csv("C:/Users/Nathalie/Dropbox/stuff/MSc_Projects/ElianeVanDerCruyssen_2020/Data_Adult_Individuals.csv", header=T, row.names = 1, na.strings=c("","NA")) Adult_data$FL <- Adult_data$humerus + Adult_data$ulna + Adult_data$digit_1 + Adult_data$digit_2 + Adult_data$digit_3 + Adult_data$digit_4 Adult_data$HL <- Adult_data$femur + Adult_data$tibia + Adult_data$toe_1 + Adult_data$toe_2 + Adult_data$toe_3 + Adult_data$toe_4 + Adult_data$toe_5 Adult_data$logrelFL <- log(Adult_data$FL)/log(Adult_data$Csize) Adult_data$logrelHL <- log(Adult_data$HL)/log(Adult_data$Csize) Adult_data_selec <- subset(Adult_data, species %in% Data$species) Adult_data_selec$species <- droplevels(Adult_data_selec$species) ``` # The datasets ```{r , echo=F} paste("Individuals per species for embryos") summary(Data$species) paste("Average individuals per species for embryos") mean(table(Data$species)) paste("Standard deviations of individuals per species for embryos") sd(table(Data$species)) paste("Individuals per species for adults") summary(Adult_data$species) paste("Average individuals per species for adults") mean(table(Adult_data$species)) paste("Standard deviations of individuals per species for adults") sd(table(Adult_data$species)) ``` # Part I: Differences in relative limb length between species (Table S2) ## Are relative limb length different between species when adult? ***Adult - Forelimb*** ```{r , echo=T} model <- lm(logrelFL ~ species, data=Adult_data_selec) Anova(model, type="II") ``` ***Adult - Hindlimb*** ```{r , echo=T} model <- lm(logrelHL ~ species, data=Adult_data_selec) Anova(model, type="II") ``` Result: Yes, species have different relative limb length when adult. ## Are relative limb length different between species already at hatching? **Forelimb:** ```{r , echo=T} Data_hatch <- subset(Data, stage == "hatch") model <- lm(logrelFL ~ species, data=Data_hatch) Anova(model, type="II") ``` **Hindlimb:** ```{r , echo=T} model <- lm(logrelHL ~ species, data=Data_hatch) Anova(model, type="II") ``` Result: Yes, species have different relative limb length at hatching. ## Are relative limb length different between species already at stage 11? **Forelimb:** ```{r , echo=T} Data_st11 <- subset(Data, stage == "11") model <- lm(logrelFL ~ species, data=Data_st11) Anova(model, type="II") ``` **Hindlimb:** ```{r , echo=T} model <- lm(logrelHL ~ species, data=Data_st11) Anova(model, type="II") ``` Result: Yes for forelimb, borderline no for hindlimb at the earliest stage. ## Are relative limb length correlated between adults and hatchlings? ```{r, echo=F, warning = FALSE, fig.show="hold", out.width="50%"} #forelimb: RelFL_adult <- aggregate(Adult_data_selec[,c("logrelFL")], list(Adult_data_selec$species), mean) colnames(RelFL_adult) <- c("Species","relFL_adult") RelFL_hatch <- aggregate(Data_hatch[,c("logrelFL")], list(Data_hatch$species), mean) colnames(RelFL_hatch) <- c("Species","relFL_hatch") Comp_FL <- merge(RelFL_adult, RelFL_hatch, by="Species") ggplot(Comp_FL,aes(x = relFL_adult, y = relFL_hatch)) + geom_point(aes(colour=Species)) + ggtitle("Relative forelimb length per species") + theme_bw() #hindlimb: RelHL_adult <- aggregate(Adult_data_selec[,c("logrelHL")], list(Adult_data_selec$species), mean) colnames(RelHL_adult) <- c("Species","relHL_adult") RelHL_hatch <- aggregate(Data_hatch[,c("logrelHL")], list(Data_hatch$species), mean) colnames(RelHL_hatch) <- c("Species","relHL_hatch") Comp_HL <- merge(RelHL_adult, RelHL_hatch, by="Species") ggplot(Comp_HL,aes(x = relHL_adult, y = relHL_hatch)) + geom_point(aes(colour=Species)) + ggtitle("Relative hindlimb length per species") + theme_bw() paste("Pearson's correlation between relative forelimb length") cor.test(Comp_FL$relFL_adult,Comp_FL$relFL_hatch) paste("Pearson's correlation between relative hindlimb length") cor.test(Comp_HL$relHL_adult,Comp_HL$relHL_hatch) ``` Answer: Yes, adult and hatchling relative hindlimb length do correlate (not sign. for forelimb), showing that species-specific patterns are largely established before hatching. # Part II: Graphically explore allometric relationships (Figure 2) ```{r graphical allometry, echo=F, results="hide"} ggplot(Data,aes(y = log(FL), x = log(spine), colour = species)) + geom_point(size=0.3) + geom_smooth(method = "lm", fill = NA) + ggtitle("Allometry of forelimbs per species") + theme_bw() ggplot(Data,aes(y = log(HL), x = log(spine), colour = species)) + geom_point(size=0.3) + geom_smooth(method = "lm", fill = NA) + ggtitle("Allometry of hindlimbs per species") + theme_bw() ``` # Part III: Test if total growth is different from isometric (Table S4) ```{r , echo=F} Data_hend <- subset(Data, species == "hendersoni") Data_marm <- subset(Data, species == "marmoratus") Data_sagr <- subset(Data, species == "sagrei") Data_cybo <- subset(Data, species == "cybotes") Data_cris <- subset(Data, species == "cristatellus") Data_saba <- subset(Data, species == "sabanus") Data_bahu <- subset(Data, species == "bahorucoensis") Data_bart <- subset(Data, species == "bartschi") Data_chri <- subset(Data, species == "christophei") Data_coel <- subset(Data, species == "coelestinus") Data_eque <- subset(Data, species == "equestris") Data_roqu <- subset(Data, species == "roquet") Data_porc <- subset(Data, species == "porcus") Data_ferr <- subset(Data, species == "ferreus") Data_reje <- subset(Data, species == "rejectus") species.list <- list(hend = Data_hend, marm = Data_marm, sagr = Data_sagr, cybo = Data_cybo, cris = Data_cris, saba = Data_saba, bahu = Data_bahu, bart = Data_bart, chri = Data_chri, eque = Data_eque, coel = Data_coel, roqu = Data_roqu, ferr = Data_ferr, porc = Data_porc, reje = Data_reje) Slopes_FL <- Slopes_HL <- array(NA, dim=c(15,4)) colnames(Slopes_FL) <- colnames(Slopes_HL) <- c("Species","mean","2.5 %","97.5 %") for (i in 1:length(species.list)) { X <- species.list[[i]] model <- lm(log(FL) ~ log(spine), data=X) Slopes_FL[i,1] <- paste((X$species)[1]) Slopes_FL[i,2] <- round(unname(coef((model))[2]),3) Slopes_FL[i,3] <- round(unname(confint(model, level=0.95)[2,1]),3) Slopes_FL[i,4] <- round(unname(confint(model, level=0.95)[2,2]),3) } for (i in 1:length(species.list)) { X <- species.list[[i]] model <- lm(log(HL) ~ log(spine), data=X) Slopes_HL[i,1] <- paste((X$species)[1]) Slopes_HL[i,2] <- round(unname(coef((model))[2]),3) Slopes_HL[i,3] <- round(unname(confint(model, level=0.95)[2,1]),3) Slopes_HL[i,4] <- round(unname(confint(model, level=0.95)[2,2]),3) } Slopes_FL <- as.data.frame(Slopes_FL) Slopes_HL <- as.data.frame(Slopes_HL) paste("Table shows the slopes of regressions against spine with confidence intervals for Forelimb") Slopes_FL paste("Table shows the slopes of regressions against spine with confidence intervals for Hindlimb") Slopes_HL paste("Mean slope for forelimbs") mean(as.numeric(as.character(Slopes_FL$mean))) paste("Mean slope for hindlimbs") mean(as.numeric(as.character(Slopes_HL$mean))) ``` Results: Growth is mostly positively allometric. Only coelestinus (FL & HL), porcus (FL) and rejectus (FL) cannot be distinguished from isometric growth. The mean slope across all species is 1.44 for Forelimbs and 1.76 for Hindlimbs. # Part IV: Statistically explore allometric relationships ## ANCOVA for total limb length by species (Table S3 and S5) First, test for interaction of species and log(spine), i.e., is there a difference in slope. **Forelimb:** ```{r , echo=T} model.1 <- lm(log(FL) ~ species * log(spine), data=Data) Anova(model.1, type="III") ``` ```{r , echo=F} contrasts <- as.data.frame(emmeans(model.1, pairwise ~ species)$contrasts) contrasts_sig <- subset(contrasts, p.value <= 0.05) paste("These are the significant contrasts (species pairs)") contrasts_sig ``` **Hindlimb:** ```{r Anolis_HL, echo=T} model.1 <- lm(log(HL) ~ species * log(spine), data=Data) Anova(model.1, type="III") ``` ```{r , echo=F} contrasts <- as.data.frame(emmeans(model.1, pairwise ~ species)$contrasts) contrasts_sig <- subset(contrasts, p.value <= 0.05) paste("These are the significant contrasts (species pairs)") contrasts_sig ``` Result: Yes, both Fore- and Hindlimb have significant interaction term, i.e. different slopes per species. ## MANCOVA for individual limb elements by species (Table S3 and S6) **Forelimb:** ```{r, echo=T} Data_FL <- as.matrix(Data[,c(9:14)]) model.2 <- manova(log(Data_FL) ~ species * log(spine), data = Data) summary((model.2)) summary.aov((model.2)) ``` **Hindlimb:** ```{r, echo=T} Data_HL <- as.matrix(Data[,c(15:21)]) model.2 <- manova(log(Data_HL) ~ species * log(spine), data = Data) summary((model.2)) summary.aov((model.2)) ```