## Simulated dataset
set.seed(14)
# Generating a random tree
tree<-pbtree(n=50)
# Setting the regime states of tip species
sta<-as.vector(c(rep("Forest",20),rep("Savannah",30))); names(sta)<-tree$tip.label
# Making the simmap tree with mapped states
tree<-make.simmap(tree,sta , model="ER", nsim=1)
col<-c("blue","orange"); names(col)<-c("Forest","Savannah")
# Plot of the phylogeny for illustration
plotSimmap(tree,col,fsize=0.6,node.numbers=FALSE,lwd=3, pts=FALSE)
# Simulate two correlated traits evolving along the phylogeny
traits<-mvSIM(tree,nsim=1, model="BMM", param=list(sigma=list(matrix(c(2,1,1,1.5),2,2),
matrix(c(4,1,1,4),2,2)), names_traits=c("head.size","mouth.size")))
# Introduce some missing cases (NA values)
data<-traits
data[8,2]<-NA
data[25,1]<-NA
# Fit of model 1
fit<-mvBM(tree,data,model="BMM")
# Estimate the missing cases
imp<-estim(tree, data, fit)
# Check the imputed data
imp$estim[1:10,]
## We want the ancestral states values at each nodes:
nodelabels() # To see where the nodes are situated
imp2<-estim(tree, data, fit, asr=TRUE)
# Check the 10 firsts ancestral states
imp2$estim[1:10,]
Run the code above in your browser using DataLab