PEM: Phylogenetic Eigenvector Map

Description Usage Arguments Details Value Author(s) References Examples

Description

Functions to calculate and manipulate Phylogenetic Eigenvector Maps (PEM).

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
  PEMInfluence(x, mroot = TRUE)
  PEMweights(d, a = 0, psi = 1)
  PEM.build(x, d = "distance", sp = "species", a = 0, psi = 1,
  tol=.Machine$double.eps^0.5)
  PEM.updater(object, a, psi = 1, tol = .Machine$double.eps^0.5)
  PEM.fitSimple(y, x, w, d = "distance", sp = "species", lower = 0,
  upper = 1, tol = .Machine$double.eps^0.5)
  PEM.forcedSimple(y, x, w, d = "distance", sp = "species", a = 0, psi = 1,
  tol = .Machine$double.eps^0.5)
  getGraphLocations(tpall, targets)
  getAncGraphLocations(x, tpall)
  Locations2PEMscores(object, gsc)

Arguments

x

A graph-class object containing a phylogenetic graph.

w

A graph-class object containing a phylogenetic graph.

object

A PEM-class object containing a Phylogenetic Eigenvector Map.

y

One or many response variable(s) in the form of a numeric vector or a matrix, respectively.

mroot

Boolean (TRUE or FALSE) specifying whether multiple rooting is allowed.

d

The name of the member of x$edge where the phylogenetic distances (edge lengths) can be found.

a

The steepness parameter describing whether changes occur, on average, progressively long edges (a close to 0) or abruptly at vertices (a close to 1).

psi

Relative evolution rate along the edges (default: 1). This parameter is only relevant when multiple values are assigned to different portions of the phylogeny.

sp

Name of the member of x$vertex where a logical vertex property specifying which vertices are species can be found. (see graph-class).

tol

Eigenvalue threshold to regard eigenvectors as usable.

lower

Lower limit for the L-BFGS-B optimization algorithm as implemented in optim.

upper

Upper limit for the L-BFGS-B optimization algorithm as implemented in optim.

tpall

Parameter of function getGraphLocations: Phylogenetic tree object of class ‘phylo’ (package ape) containing all species (model and target) used in the study.

targets

Name of the target species to extract using the tpall.

gsc

The output of getGraphLocations.

Details

Functions PEMInfluence and PEMweights are used internally by PEM.build to create a binary matrix referred to as an ‘influence matrix’ and weight its columns. That matrix has a row for each vertex of graph ‘x’ and a column for each of its edges. The elements of the influence matrix are 1 whenever the vertex associated with a row is located in the tree either directly or indirectly downward the edge associated with a column. That function is implemented in C language using recursive function calls. Although PEMInfluence allows one to use multiple roots as its default parameter, it is called within PEM.build with mroot = FALSE. User must therefore ensure that the graph provided to PEMap is single-rooted.

Function PEM.build is used to produce a phylogenetic eigenvector map, while function PEM.updater allows one to re-calculate a PEM-class object with new weighting function parameters. Function PEM.fitSimple performs a maximum likelihood estimation of a and psi assuming single values for the whole tree whereas function PEM.forcedSimple allows one the force parameters a and psi to a PEM-class object while adding the same computational details as those PEM.fitSimple would have produced (and which are necessary to make predictions).

Functions getGraphLocations returns the coordinates of a species in terms of its position with respect to the influence matrix while function Locations2PEMscores transforms these coordinates into sets of scores that can be used to make predictions. Function getAncGraphLocations produce the same output as getGraphLocations, but of the ancestral species (i.e. the nodes of the phylogeny) in order to estimate ancestral trait values.

Value

Function PEMInfluence returns the influence matrix of graph x and function PEMweights returns weights corresponding to the distances. Functions PEM.build, PEM.fitSimple, PEM.forcedSimple returns a PEM-class object. Function getGraphLocations returns a list whose first member is an influence coordinates matrix whose rows refer to the target species and columns refer to the edges and second member is the lengths of the terminal edges connecting each target species to the rest of the phylogeny. Function Locations2PEMscores returns a list whose first member is a PEM score matrix whose rows refer to the target species and columns refer to the eigenvectors and second member is the variance associated with the terminal edges connecting the target species to the phylogeny.

Author(s)

Guillaume Guénard, Département de sciences biologiques Université de Montréal, Montréal, QC, Canada.

References

Guénard, G., Legendre, P., and Peres-Neto, P. 2013. Phylogenetic eigenvector maps (PEM): a framework to model and predict species traits. Meth. Ecol. Evol. In press.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
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;",sep=""))
x <- Phylo2DirectedGraph(t1)
#
## Calculates the (binary) influence matrix
PEMInfluence(x)
PEMInfluence(x)[x$vertex$species,]
#
## Building phylogenetic eigenvector maps
PEM1 <- PEM.build(x)
print(PEM1)
PEM2 <- PEM.build(x, a = 0.2)
PEM3 <- PEM.build(x, a = 1)
PEM4 <- PEM.updater(PEM3,a=0.5)
#
## Extracts the eigenvectors
as.data.frame(PEM4)
#
### Example of an hypothetical set of trait values
y <- c(A=-1.1436265,B=-0.3186166,C=1.9364105,D=1.7164079,E=1.0013993,
       F=-1.8586351,G=-2.0236371)
#
## Estimate single steepness parameter for the whole tree.
PEMfs1 <- PEM.fitSimple(y=y,x=NULL,w=x,d="distance",sp="species",lower=0,upper=1)
PEMfs1$optim       # Results of the optimization.
#
## Force neutral evolution for the whole tree.
PEMfrc1 <- PEM.forcedSimple(y=y,x=NULL,w=x,d="distance",sp="species",a=0)
PEMfrc1$x$edge$a   # Steepness parameter forced for each individual edge.
#
## Get graph locations for target species X, Y, and Z
tpAll <- 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=""))
grloc <- getGraphLocations(tpAll, c("X","Y","Z"))
#
PEMfs2 <- PEM.fitSimple(y=y,x=NULL,w=grloc$x,d="distance",sp="species",lower=0,upper=1)
PEMfs2$optim       # Same as for PEMfs1$optim
#
PEMsc1 <- Locations2PEMscores(PEMfs2, grloc)
lm1 <- lm(y~V_2+V_3+V_5,data=PEMfs2)
#
ypred <- predict(object=PEMfs2,targets=grloc,lmobject=lm1,interval="none")
#
tpModel <- drop.tip(tpAll,c("X","Y","Z"))
#
## Plotting the results:
layout(t(c(1,1,2)))
par(mar=c(6,2,2,0.5)+0.1)
plot(tpModel,show.tip.label=TRUE,show.node.label=TRUE,root.edge = TRUE,
     srt = 0,adj=0.5,label.offset=0.08,font=1,cex=1.5,xpd=TRUE)
edgelabels(paste("E",1:nrow(tpModel$edge),sep=""),
           edge=1:nrow(tpModel$edge),bg="white",font=1,cex=1)
points(x=0.20,y=2.25,pch=21,bg="black")
lines(x=c(0.20,0.20,0.65),y=c(2.25,0.55,0.55),xpd=TRUE,lty=2)
text("X",x=0.69,y=0.55,xpd=TRUE,font=1,cex=1.5)
points(x=0.35,y=4.5,pch=21,bg="black")
lines(x=c(0.35,0.35,0.6),y=c(4.5,5.47,5.47),xpd=TRUE,lty=2)
text("Y",x=0.64,y=5.47,xpd=TRUE,font=1,cex=1.5)
points(x=0.35,y=3,pch=21,bg="black")
lines(x=c(0.35,0.35,0.55),y=c(3,3.5,3.5),xpd=TRUE,lty=2)
text("Z",x=0.59,y=3.5,xpd=TRUE,font=1,cex=1.5)
text(c("NX","NY","NZ"),x=c(0.20,0.35,0.35),y=c(2.25,4.5,3)+0.3*c(1,-1,-1),font=1,cex=1)
add.scale.bar(length=0.1,cex=1.25)
par(mar=c(3.75,0,2,2)+0.1)
plot(x=y,y=1:7,ylim=c(0.45,7),xlim=c(-4,4), axes=FALSE, type="n", xlab="")
axis(1,label=c("-4","-2","0","2","4"),at=c(-4,-2,0,2,4)) ; abline(v=0)
points(x=y,y=1:7,xlim=c(-2,2),pch=21,bg="black")             # Observed values
text("B)",x=-3.5,y=7,cex=1.5,xpd=TRUE) ; text("Trait value",x=0,y=-0.5,cex=1.25,xpd=TRUE)
points(x=ypred,y=c(0.5,5.5,3.5),pch=23,bg="white",cex=1.25)  # Predicted values
#
## Estimating ancestral trait values:
ANCloc <- getAncGraphLocations(x)
PEMfsAnc <- PEM.fitSimple(y=y,x=NULL,w=ANCloc$x,d="distance",sp="species",lower=0,upper=1)
PEMfsAnc$optim
#
PEManc1 <- Locations2PEMscores(PEMfsAnc, ANCloc)
y_anc <- predict(object=PEMfsAnc,targets=ANCloc,lmobject=lm1,interval="confidence")
#

Example output

Loading required package: ape
Loading required package: MASS
   E1 E2 E3 E4 E5 E6 E7 E8 E9 E10 E11 E12
A   1  1  1  0  0  0  0  0  0   0   0   0
B   1  1  0  1  0  0  0  0  0   0   0   0
C   1  0  0  0  1  0  0  0  0   0   0   0
D   0  0  0  0  0  1  1  1  0   0   0   0
E   0  0  0  0  0  1  1  0  1   0   0   0
F   0  0  0  0  0  1  0  0  0   1   1   0
G   0  0  0  0  0  1  0  0  0   1   0   1
N1  0  0  0  0  0  0  0  0  0   0   0   0
N2  1  0  0  0  0  0  0  0  0   0   0   0
N4  1  1  0  0  0  0  0  0  0   0   0   0
N3  0  0  0  0  0  1  0  0  0   0   0   0
N5  0  0  0  0  0  1  1  0  0   0   0   0
N6  0  0  0  0  0  1  0  0  0   1   0   0
  E1 E2 E3 E4 E5 E6 E7 E8 E9 E10 E11 E12
A  1  1  1  0  0  0  0  0  0   0   0   0
B  1  1  0  1  0  0  0  0  0   0   0   0
C  1  0  0  0  1  0  0  0  0   0   0   0
D  0  0  0  0  0  1  1  1  0   0   0   0
E  0  0  0  0  0  1  1  0  1   0   0   0
F  0  0  0  0  0  1  0  0  0   1   1   0
G  0  0  0  0  0  1  0  0  0   1   0   1
A phylogenetic eigenvector map (PEM) for  7  species:
A,B,C,D,E,F,G 
obtained from the following phylogenetic graph:

A graph with 12 edges and 13 vertices. 
Edge labels: E1 E2 E3 E4 E5 E6 E7 E8 E9 E10 E11 E12 
Vertex labels: A B C D E F G N1 N2 N4 N3 N5 N6 
Available edge information:  distance 
Available vertex information:  species 

         V_1           V_2         V_3         V_4           V_5          V_6
A -0.4739477  0.0003520888 -0.32451718  0.12046476 -0.7158877081  0.014068974
B -0.4869085  0.0003718940 -0.38520154 -0.11713112  0.6767099858 -0.005147892
C -0.3394274  0.0001256851  0.86047786 -0.01297167  0.0323610484  0.017153514
D  0.3376007 -0.5407991012 -0.04695301 -0.05689084  0.0006053317  0.667280651
E  0.3108811 -0.4577115750 -0.02783727  0.04702704 -0.0005050412 -0.740276438
F  0.3215049  0.4851845416 -0.03473658  0.70471926  0.1227273054  0.073996965
G  0.3302969  0.5124764666 -0.04123228 -0.68521743 -0.1160109221 -0.027075774
$par
[1] 0.5236929

$value
[1] 21.84892

$counts
function gradient 
       6        6 

$convergence
[1] 0

$message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

 [1] 0 0 0 0 0 0 0 0 0 0 0 0
$par
[1] 0.5236929

$value
[1] 21.84892

$counts
function gradient 
       6        6 

$convergence
[1] 0

$message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

$par
[1] 0.5236929

$value
[1] 21.84892

$counts
function gradient 
       6        6 

$convergence
[1] 0

$message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

MPSEM documentation built on June 4, 2019, 1:02 a.m.