PEM-class: Class and Methods for Phylogenetic Eigenvector Maps (PEM)

In guenardg/MPSEM: Modelling Phylogenetic Signals using Eigenvector Maps

Class and Methods for Phylogenetic Eigenvector Maps (PEM)

Description

Class and methods to calculate and manipulate Phylogenetic Eigenvector Maps (PEM).

Usage

PEM(
  x,
  ...,
  a = -Inf,
  mm_a = matrix(1, nedge(x), 1L),
  psi = NULL,
  mm_psi = matrix(NA, nedge(x), 0L),
  d = "distance",
  sp = "species",
  tol = .Machine$double.eps^0.5
)

## S3 method for class 'PEM'
print(x, ...)

## S3 method for class 'PEM'
as.matrix(x, ...)

## S3 method for class 'PEM'
as.data.frame(x, row.names = NULL, optional = FALSE, ...)

## S3 method for class 'PEM'
update(object, a, psi, ...)

## S3 method for class 'PEM'
evolution.model(object, y, ..., x = NULL)

## S3 method for class 'PEM'
predict(object, newdata, ...)

Arguments

x	A graph-class object, a PEM-class object, or a matrix containing auxiliary traits.
...	Additional parameters to be passed to the method. Currently ignored.
a	A numeric vector of linear coefficients to estimate the steepness in the trait evolution model; it must have as many values as the number of columns in argument mm_a (default: -Inf).
mm_a	A numeric matrix with as many rows as the number of edges and as many columns as the number of values of argument a (default: an all-ones single-column matrix).
psi	A numeric vector of linear coefficients to estimate the evolution rate in the trait evolution model; it must have as many values as the number of columns in argument mm_psi (default: NULL, meaning an homogeneous evolution rate throughout the graph).
mm_psi	A numeric matrix with as many rows as the number of edges and as many columns as the number of values of argument psi (default: a zero-column matrix).
d	A character string; the name of the (numeric) edge property containing the phylogenetic distances (default: "distance").
sp	A character string; the name of the (logical) vertex property containing vertex property (default: "species").
tol	A numeric; the lowest singular value above which to retain a singular vector as part of a PEM (default: .Machine$double.eps^0.5).
row.names	Included for method consistency reason; ignored.
optional	Included for method consistency reason; ignored.
object	A PEM-class object.
y	A vector or matrix of (numeric) trait values to estimate the parameters of the trait evolution model.
newdata	A three-column data frame giving the locations of a set of targets in the graph such as the one provided by function locate.graph (element $location; see Details below).

Format

A PEM-class object contains:

d

the name of the edge property containing the phylogenetic distances (see argument d above),

sp

the name of the vertex property identifying which which of the vertices are bearing trait values (see argument sp above),

tol

the tolerance value (see argument tol above),

mm

a list with the model matrices (see arguments mm_a and mm_psi above) that are used to estimate the steepness and relative evolution rate throughout the graph manifold,

nsp

the number of vertices that bear trait values,

B

the influence matrix (see inflMat),

print

function; prints the PEM object,

graph

function; returns the graph object used to calculate the PEM,

pem

function; returns a list of components (ie., vectors and matrices) involved in PEM calculation,

par

function; returns the actual model parameter values as a list,

nodal

function; returns the steepness and evolution rates at any given vertex of the graph,

update

function; updates the PEM with a new set of parameters,

var

function; returns the PVCV matrix: phylogenetic variance and covariances matrix among the vertices with trait values,

inv_var

function; returns the inverse PVCV matrix,

logdet

function; returns the natural logarithm of the PVCV matrix determinant,

dev

function; calculates the deviance of a target trait, or many such traits, with respect to the PVCV matrix, optionally conditional to one or more auxiliary trait(s),

S2

function; calculates the variance(s) of one or more trait(s) with respect to the PVCV matrix or the residual variance(s) of said trait(s) conditioned to one or more auxiliary trait(s),

predictVertex

function; calculates the PEM prediction scores at the graph's vertices, and

predictEdge

function; calculates the PEM prediction scores at the graph's edge.

Details

The class is implemented using function PEM. It provides methods as.matrix and as.data.frame to extract the phylogenetic vectors in the forms of a matrix or a data frame, respectively, an update method to recalculate an existing PEM with new parameters (ie., arguments a and psi), an evolution.model method to empirically estimate parameter values from a data set of traits (ie., argument y and x), and a predict method to estimate PEM scores for arbitrary graph locations (using argument newdata).

Graph locations are given as a three-column data frame:

ref

the index of the vertex or edge where each target can be found,

dist

the phylogenetic distance along the edge where each target can be found, or NA for a target located at a vertex,

lca

the phylogenetic distance between the latest common ancestor of the target in the graph and the target itself.

One such data frame is one of the second element of the list produced by method locate.graph (called "location").

Functions

• PEM(): Create PEM-class

  A function creating a PEM-class object.

• print(PEM): Print PEM-class

  A print method for PEM-class objects.

• as.matrix(PEM): Method as.matrix for PEM-class Objects

  A method to extract the phylogenetic eigenvectors from a PEM-class object.

• as.data.frame(PEM): Method as.data.frame for PEM-class Objects

  A method to extract the phylogenetic eigenvectors from a PEM-class object.

• update(PEM): Update Method for PEM-class Objects

  An updating method to recalculate a PEM with new smoothness and evolution rate parameters.

• evolution.model(PEM): Evolution Model Method for PEM-class Objects

  A method to estimate the smoothness and evolution rate parameters of a PEM-class object empirically on the basis of an observed quantitative trait and, optionally, optional trait(s).

• predict(PEM): Predict Method for PEM-class Objects

  A predict method to obtain PEM prediction scores for arbitrary graph locations.

Author(s)

Guillaume Guénard [aut, cre] (<https://orcid.org/0000-0003-0761-3072>), Pierre Legendre [ctb] (<https://orcid.org/0000-0002-3838-3305>) Maintainer: Guillaume Guénard <[e-mail]>

References

Guénard, G., Legendre, P., and Peres-Neto, P. 2013. Phylogenetic eigenvector maps: a framework to model and predict species traits. Methods in Ecology and Evolution 4: 1120-1131

See Also

PEM-functions

Examples

## Synthetic example

## This example describes the phyogeny of 7 species (A to G) in a tree with 6
## nodes, presented in Newick format, read by function
## read.tree of package ape.

t1 <- read.tree(text=paste(
  "(((A:0.15,B:0.2)N4:0.15,C:0.35)N2:0.25,((D:0.25,E:0.1)N5:0.3,",
  "(F:0.15,G:0.2)N6:0.3)N3:0.1)N1:0.1;",sep=""))
t1
summary(t1)

## Example of a made up set of trait values for the 7 species
y <- c(A=-1.1436265,B=-0.3186166,C=1.9364105,D=1.7164079,E=1.0013993,
       F=-1.8586351,G=-2.0236371)

## The phylogenetic tree is turned into an evolutionary graph as follows:
x <- as.graph(t1)

## The graph:
x

## This is the edge information of that simple graph:
edge(x)

## Calculate the (binary) influence matrix; E1 to E12 are the tree edges:
IM1 <- InflMat(x)
IM1

## This is the edges associated with the influence matrix:
edge(IM1)

## This is the influence matrix for the tree leaves (ie., the graph's 
## terminal vertices):
IM1[x$species,]

## Here is a suite of weighting function profiles:
seq(0,1.5,0.001) %>%
  plot(y=PEMweights(., a=0), ylim=c(0,1.7), type="l", xlab="distance",
       ylab="weight")

seq(0,1.5,0.001) %>%
  lines(y=PEMweights(., a=0.5), col="red")

seq(0,1.5,0.001) %>%
  lines(y=PEMweights(., a=0.5, psi=1.5), col="green")

seq(0,1.5,0.001) %>%
  lines(y=PEMweights(., a=0.9), col="blue")

## A Phylogenetic eigenvector maps (PEM) is build with the default parameters
## as follows:
PEM1 <- PEM(x)

## The PEM object:
PEM1

## The as.matrix methods returns the phylogenetic eigenvectors as a matrix:
as.matrix(PEM1)

## The as.data.frame methods returns the phylogenetic eigenvectors as a data
## frame:
as.data.frame(PEM1)

## The update method recalculates the PEM with new parameters:
update(PEM1, a=log(0.25/(1 - 0.25)), psi=NULL)

## Note: The way 'a' is estimated is based on an inverse-logit-transformed
## linear sub-model. This approach is used to facilitate the bounding of the
## steepness parameter in the [0,1] interval. Therefore, to set a given
## steepness value, one has to provide its logit-transformed value. As a
## reminder, the logit transformation is f(x) = log(x/(1 - x)).

## To show the edges of the graph within the PEM object:
edge(PEM1$graph())

## A second PEM is initialized with a = 0.2 as follows:
PEM2 <- PEM(x, a=log(0.2/(1 - 0.2)))
PEM2

## The edges of PEM2 graph:
edge(PEM2$graph())

## A third PEM is initialized with a = 1 as follows:
PEM3 <- PEM(x, a=Inf)
PEM3

## The edges of PEM3 graph:
edge(PEM3$graph())

## Estimation on the steepness parameter empirically using data set y
## Initial value must be finite (not -Inf nor Inf, meaning a > 0 and a < 1):
PEM4 <- PEM(x, a=0)

## Method evolution.model carry out the estimation:
opt <- evolution.model(PEM4, y=y)
opt

## The edges of PEM4 graph after the estimation:
edge(PEM4$graph())

## Graph locations for target species X, Y, and Z not found in the original
## data set:
read.tree(
  text = paste(
    "((X:0.45,((A:0.15,B:0.2)N4:0.15,(C:0.25,Z:0.2)NZ:0.1)N2:0.05)NX:0.2,",
    "(((D:0.25,E:0.1)N5:0.05,Y:0.25)NY:0.25,(F:0.15,G:0.2)N6:0.3)N3:0.1)N1;",
    sep=""
  )
) -> tree
tree

## The tree is converted into a graph as follows:
x2 <- as.graph(tree)
x2

## The target X, Y, and Z are extracted from graph x2 as follows:
loc <- locate(x2, target = c("X","Y","Z"))
loc

## The list contains the residual graph and a data frame with the locations
## of the target in the residual graph.

## Building a PEM on the residual graph (assuming a = 0):
PEM5 <- PEM(loc$x)

## For making predictions, we need to do any of the following:
## 1. run evolution.model(); doing this operation will change the parameter
##    values,
## 2  estimate the deviance using $dev(y, ...)
## 3. estimate the variance using $S2(y, ...)

## Before any of these operations is done, no variance estimate is available:
PEM5$S2()

## Presenting a variable to the PEM makes its variance estimate available
## later on. For instance, presenting 'y' to the PEM5 while estimating its
## deviance as follows:
PEM5$dev(y)

## This operation makes the variance estimate for 'y' available for later
## access:
PEM5$S2()

## The predict method is used to generate the predictors as follows:
scr <- predict(PEM5, loc$location)
scr

## The 'vf' attribute, which can be accessed as follows:
attr(scr,"vf")
## is the amount of variance associated with the evolutionary distance
## between each target species and its latest common ancestor in the
## evolutionary graph.

## Given a model built using the PEM as follows:
lm1 <- lm(y ~ U_2 + U_3 + U_5, data=as.data.frame(PEM5))

## Predictions are obtained as follows:
ypred <- predict(lm1, as.data.frame(scr))
ypred

## Estimate the ancestral trait values

## The ancestors of the species of the tree are the vertices without extant
## species:
data.frame(
  ref = which(!x$species),
  dist = rep(NA_real_, sum(!x$species)),
  lca = rep(0, sum(!x$species)),
  row.names = rownames(x)[!x$species]
) -> anc

## The scores of the ancestors are obtained as follows:
src_anc <- predict(PEM1, anc)
## Note the absence of the 'vf' attributes since no variable has been
## presented to PEM1.

## Predictions are obtained from the previous linear model as follows:
predict(lm1, as.data.frame(src_anc))

guenardg/MPSEM documentation built on April 14, 2025, 3:53 p.m.