Ancestral character estimation with a trend
This function estimates the evolutionary parameters and ancestral states for Brownian evolution with directional trend.
anc.trend(tree, x, maxit=2000)
tree |
an object of class |
x |
a vector of tip values for species; |
maxit |
an optional integer value indicating the maximum number of iterations for optimization. |
Note that this will generally only work and produce sensible results for a phylogeny with some non-contemporary tips (i.e., a tree with some fossil species). The function uses optim with method= "L-BFGS-B"; however optimization is only constrained for the sig2 which must be >0.
An object of class "anc.trend" with the following components:
ace |
a vector with the ancestral states. |
mu |
a trend parameter per unit time. |
sig2 |
the variance of the BM process, σ^2. |
logL |
the log-likelihood. |
convergence |
the value of |
Liam Revell liam.revell@umb.edu
Revell, L. J. (2012) phytools: An R package for phylogenetic comparative biology (and other things). Methods Ecol. Evol., 3, 217-223.
## simulate tree & data using fastBM with a trend (m!=0)
tree<-rtree(n=26,tip.label=LETTERS)
x<-fastBM(tree,mu=4,internal=TRUE)
a<-x[as.character(1:tree$Nnode+Ntip(tree))]
x<-x[tree$tip.label]
## fit no trend model
fit.bm<-anc.ML(tree,x,model="BM")
print(fit.bm)
## fit trend model
fit.trend<-anc.trend(tree,x)
print(fit.trend)
## compare trend vs. no-trend models & estimates
AIC(fit.bm,fit.trend)
layout(matrix(c(3,4,1,2,5,6),3,2,byrow=TRUE),
heights=c(1.5,3,1.5),widths=c(3,3))
xlim<-ylim<-range(c(a,fit.bm$ace,
fit.trend$ace))
plot(a,fit.bm$ace,pch=19,
col=make.transparent("blue",0.5),
xlab="true ancestral states",
ylab="ML estimates",
main=paste("Comparison of true and estimated",
"\nstates under a no-trend model"),font.main=3,
cex.main=1.2,bty="l",cex=1.5,
xlim=xlim,ylim=ylim)
lines(xlim,ylim,lty="dotted")
plot(a,fit.trend$ace,pch=19,
col=make.transparent("blue",0.5),
xlab="true ancestral states",
ylab="ML estimates",
main=paste("Comparison of true and estimated",
"\nstates under a trend model"),font.main=3,
cex.main=1.2,bty="l",cex=1.5,
xlim=xlim,ylim=ylim)
lines(xlim,ylim,lty="dotted")
par(mfrow=c(1,1))Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.