library(ape)
library(geiger)
library(picante)
library(ade4)
library(phytools)
library(car)
### Example from Lynch 1991
# This code writes the phylogeny in "tre" format to a file, then loads it with read.tree
cat("((((Homo:0.21,Pongo:0.21):0.28,","Macaca:0.49):0.13,Ateles:0.62):0.38,Galago:1.00);",file = "ex.tre", sep = "\n")
tree.primates <- read.tree("ex.tre")
# Weight and longevity data for the primates in the tree
weight <- c(60, 37, 10.7, 7.6, 0.23)
longevity <- c(115, 28, 29, 18, 10)
# names of the primates in the dataset
names(weight) <- names(longevity) <- c("Homo", "Pongo", "Macaca", "Ateles", "Galago")
# log transform data
weight<-log(weight)
longevity<-log(longevity)
# plot data
plot(weight,longevity,pch=16)
# plot tree
plot(tree.primates)
# plot species weight on phylogeny
contMap(tree.primates,weight)
# plot species longevity on phylogeny
contMap(tree.primates,longevity)
# Check both for strength of phylogenetic signal
phylosignal(longevity, tree.primates, reps=1000, checkdata=TRUE)
phylosignal(weight, tree.primates, reps=1000)
# calculate independent contrasts for weight and longevity
pic.weight <- pic(weight, tree.primates)
pic.longevity <- pic(longevity, tree.primates)
# form models for weight and longevity with and without phylogeny
weight_long<-lm(weight~longevity)
weight_long_phylo<-lm(pic.weight~pic.longevity)
# Anova tables for each model
Anova(weight_long)
summary(weight_long)
Anova(weight_long_phylo)
summary(weight_long_phylo)
# Tests for phylogenetic signal with four models
# Brownian motion
brownFit <-fitContinuous(tree.primates, longevity)
# Pagel's lambda
lambdaFit<-fitContinuous(tree.primates, longevity, model="lambda")
# Ornstein-Uhlenbeck model
ouFit<-fitContinuous(tree.primates, longevity, model="OU")
# White noise
whiteFit<-fitContinuous(tree.primates, longevity, model="white")
# Get the AIC values from models
model.aic<-c(brownFit$opt$aicc,lambdaFit$opt$aicc,ouFit$opt$aicc,whiteFit$opt$aicc)
# Print delta AIC in order: Brownian, lambda, ou, white noise
model.aic-min(model.aic)
# Do the same but with weight
brownFit <-fitContinuous(tree.primates, weight)
lambdaFit<-fitContinuous(tree.primates, weight, model="lambda")
ouFit<-fitContinuous(tree.primates, weight, model="OU")
whiteFit<-fitContinuous(tree.primates, weight, model="white")
# Get the AIC values from models
model.aic<-c(brownFit$opt$aicc,lambdaFit$opt$aicc,ouFit$opt$aicc,whiteFit$opt$aicc)
# Print delta AIC in order: Brownian, lambda, ou, white noise
model.aic-min(model.aic)