pruning | R Documentation |
This function uses the pruning algorithm (Felsenstein 1973) to efficiently compute the determinant of the phylogenetic covariance matrix as well as the square root of this matrix (or its inverse; Stone 2011, Khabbazian et al. 2016). This algorithm is faster than using "eigen" or "cholesky" function to compute the determinant or the square root (see e.g., Clavel et al. 2015) and can be used to compute independent contrasts scores and the log-likelihood of a model in linear time.
pruning(tree, inv=TRUE, scaled=TRUE, trans=TRUE, check=TRUE)
tree |
Phylogenetic tree (an object of class "phylo" or "simmap"). |
inv |
Return the matrix square root of either the covariance matrix (inv=FALSE) or its inverse (inv=TRUE, the default). This matrix is a "contrasts" matrix. |
scaled |
Indicates whether the contrasts should be scaled with their expected variances (default to TRUE). |
trans |
Return the transpose (trans=TRUE) of the matrix square root/contrasts matrix. (by default - i.e., trans=TRUE - it returns a matrix equivalent to the upper triangular Cholesky factor) |
check |
Check if the input tree is dichotomous and in "postorder" (see ?is.binary.tree and ?reorder.phylo). |
The tree is assumed to be fully dichotomic and in "postorder", otherwise the functions multi2di and reorder.phylo are used internally when check=TRUE.
sqrtMat |
The matrix square root (contrasts matrix) |
varNode |
Variance associated to each node values (similar to "contrasts" variance) |
varRoot |
Variance associated to the root value (similar to the ancestral state variance) |
det |
Log-determinant of the phylogenetic covariance of the tree |
Julien Clavel
Clavel J., Escarguel G., Merceron G. 2015. mvMORPH: an r package for fitting multivariate evolutionary models to morphometric data. Methods Ecol. Evol. 6:1311-1319.
Felsenstein J. 1973. Maximum-likelihood estimation of evolutionary trees from continuous characters. Am. J. Hum. Genet. 25:471-492.
Khabbazian M., Kriebel R., Rohe K., Ane C. 2016. Fast and accurate detection of evolutionary shifts in Ornstein-Uhlenbeck models. Methods Ecol. Evol. 7:811-824.
Stone E.A. 2011. Why the phylogenetic regression appears robust to tree misspecification. Syst. Biol. 60:245-260
mvLL
mvgls
## Simulated dataset
set.seed(14)
# Generating a random tree
tree<-pbtree(n=50)
Y <- mvSIM(tree, model="BM1", param=list(sigma=1, theta=0)) # trait
X <- matrix(1, nrow=Ntip(tree), ncol=1) # design matrix
## Use the GLS trick
# Compute the matrix square root
C <- vcv.phylo(tree)
D <- chol(C)
Cinv <- solve(C)
Di <- chol(Cinv)
# transform the traits
Xi <- Di%*%X
Yi <- Di%*%Y
# Compute the GLS estimate and determinant (see Clavel et al. 2015)
# GLS estimate for the root
print(pseudoinverse(Xi)%*%Yi)
# Determinant of the phylogenetic covariance matrix
print(sum(log(diag(D)^2)))
## Use the pruning algorithm (much faster)
M <- pruning(tree, inv=TRUE)
Xi <- M$sqrtMat%*%X
Yi <- M$sqrtMat%*%Y
# GLS estimate
print(pseudoinverse(Xi)%*%Yi)
# determinant
print(M$det)
## REML determinant (without variance of the root state; see Felsenstein 1973)
# full REML
log(det(C)) + log(det(t(X)%*%Cinv%*%X))
# pruning REML
sum(log(M$varNode))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.