#################################### #Load Packages #################################### library(geiger) library(phytools) library(ape) #################################### #Set working directory & load files #################################### setwd("your_working_directory") MHCdata<-read.csv("diversity_data.txt") MHCtree<-read.nexus(file="tree.nex") #################################### #Data Checks #################################### #Distribution of data hist(MHCdata$Diversity_Estmiate) #Check for polytomies "true" means it has binary burfication is.binary.tree(MHCtree) #MHCtree<- multi2di(MHCtree, random = TRUE) #Check ultrametric i.e. has the tree been scaled so all brached end at same point. "Ture" if ultrametric is.ultrametric(MHCtree) max(MHCtree$edge.length) min(MHCtree$edge.length) # to ensure there are no minus values at nodes #################################### #Test best model of evolution using Geiger: BM, OU, versus D #################################### MHCgeiger<-as.matrix(log(MHCdata[,"Diversity_Estmiate"])) # puts it into correct format for Geiger rownames(MHCgeiger)<-MHCdata$Spp #Gives row names to new doc MHC.BM<-fitContinuous(MHCtree,MHCgeiger,model="BM") MHC.OU<-fitContinuous(MHCtree,MHCgeiger,model="OU") MHC.DT<-fitContinuous(MHCtree,MHCgeiger,model="delta") MHC.BM MHC.OU MHC.DT #################################### #Estimate phylogenetic signal using Phytools #################################### MHCrow<-MHCdata$Diversity_Estmiate names(MHCrow)<-MHCdata$Spp # to Phytools format phylosig(MHCtree,MHCrow,method="lambda",test=TRUE)