three.point.compute: Computations with a (generalized) three-point structured tree

In phylolm: Phylogenetic Linear Regression

Computations with a (generalized) three-point structured tree

Description

Computes P'V^{-1}Q and the \log(\det V) of a (generalized) three-point structured matrix.

Usage

three.point.compute(phy, P, Q = NULL, diagWeight = NULL, 
            check.pruningwise = TRUE, check.names = TRUE, check.precision = TRUE)

Arguments

phy	a rooted phylogenetic tree of type phylo with branch lengths, to represent the 3-point structured matrix V_0. Note that the matrix of interest is V = D V_0 D.
P, Q	two matrices.
diagWeight	a vector containing the diagonal elements of the diagonal matrix D if V has a generalized 3-point structure: V = D V_0 D
check.pruningwise	If FALSE, the tree is assumed to be in pruningwise order.
check.names	if FALSE, the row names of P, Q, and the names of diagWeight are assumed to be the same as the labels of the tips in the tree.
check.precision	if FALSE, diagWeight will be allowed to be below Machine epsilon. Recommended to remain TRUE.

Value

vec11	1'V^{-1}1.
P1	P'V^{-1}1.
PP	P'V^{-1}P.
Q1	Q'V^{-1}1.
QQ	Q'V^{-1}Q.
PQ	P'V^{-1}Q.
logd	\log(\det V).

Note

The matrix V is assumed to be V = D V_0 D where D is the diagonal matrix with non-zero diagonal elements in diagWeight, and where V_0 is the 3-point structured covariance matrix determined by phy and its branch lengths. Note that D do not correspond to measurement error terms.

The number of rows in P and Q and the length of diagWeight need to be the same as the number of tips in the tree. When Q = NULL, the function only returns 1'V^{-1}1, P'V^{-1}1 and P'V^{-1}P.

Author(s)

Lam Si Tung Ho, Robert Lachlan

References

Ho, L. S. T. and An?, C. (2014). "A linear-time algorithm for Gaussian and non-Gaussian trait evolution models". Systematic Biology 63(3):397-408.

See Also

transf.branch.lengths.

Examples

tre1 = rtree(500)
tre2 = transf.branch.lengths(phy=tre1, model="OUrandomRoot",
                             parameters = list(alpha = 0.5))
Q = rTrait(n=2,tre1)
y = rTrait(n=1,tre1)
P = cbind(1,y)
three.point.compute(tre2$tree,P,Q,tre2$diagWeight)

phylolm documentation built on Oct. 1, 2024, 1:09 a.m.