| trait-simulator {MPSEM} | R Documentation |
Simulates the Evolution of a Quantitative Trait.
Description
Functions to simulate the evolution of a quantitative trait
along a phylogenetic tree inputted as an object of class ‘phylo’
(package ape) or graph-class object.
Usage
EvolveOptimMarkovTree(tp, tw, anc, p = 1, root = tp$edge[1, 1])
TraitOUsimTree(tp, a, sigma, opt, p = 1, root = tp$edge[1, 1])
OUvar(d, a = 0, theta = 1, sigma = 1)
PEMvar(d, a = 0, psi = 1)
TraitVarGraphSim(x, variance, distance = "distance", p = 1, ...)
Arguments
tp |
A rooted phylogenetic tree of class ‘phylo’ (see package ape). |
tw |
Transition matrix giving the probability that the optimum trait value changes from one state to another at vertices. All rows must sum to 1. |
anc |
Ancestral state of a trait (at the root). |
p |
Number of variates to generate. |
root |
Root node of the tree. |
a |
|
sigma |
Neutral evolution rate, i.e. mean trait shift by drift. |
opt |
An index vector of optima at the nodes. |
d |
Phylogenetic distances (edge lengths). |
theta |
Adaptive evolution rate, i.e. mean trait shift by natural selection. |
psi |
Mean evolution rate. |
x |
A |
variance |
Variance function ( |
distance |
The name of the member of ‘x$edge’ where edge lengths can be found. |
... |
Additional parameters for the specified variance function. |
Details
Function EvolveOptimMarkovTree allows one to simulate the
changes of optimum trait values as a Markov process. The index whereby the
process starts, at the tree root, is set by parameter anc; this is the
ancestral character state. From the root onwards to the tips, the optimum is
given the opportunity to change following a multinomial random draw with
transition probabilities given by the rows of matrix tw. The integers
thus obtained can be used as indices of a vector featuring the actual optimum
trait values corresponding to the simulated selection regimes. The resulting
optimum trait values at the nodes are used by TraitOUsimTree as
its parameters opt to simulate trait values at nodes and tips.
Function TraitVarGraphSim uses a graph variance function
(either OUvar or PEMvar) to reconstruct a covariance matrix
that is used to generate covariates drawn from a multi-normal distribution.
Value
Functions EvolveOptimMarkovTree and
TraitOUsimTree return a matrix whose rows represent the
vertices (nodes and tips) of the phylogenetic tree and whose columns stand
for the n different trials the function was asked to perform. For
EvolveQTraitTree, the elements of the matrix are integers,
representing the selection regimes prevailing at the nodes and tips, whereas
for TraitOUsimTree, the elements are simulated quantitative
trait values at the nodes and tips. These functions are implemented in C
language and therefore run swiftly even for large (10000+ species) trees.
Function TraitVarGraphSim returns p phylogenetic signals
and is implemented using a rotation of a matrix of standard normal random
(mean=0, variance=1) deviates. The rotation matrix is itself obtained by
Choleski factorization of the trait covariance matrix expected for a given
set of trees, variance function, and variance function parameters.
Functions
-
EvolveOptimMarkovTree: Simulates the evolution of trait optima along a phylogeny. -
TraitOUsimTree: Simulates the evolution of trait values along a phylogeny. -
OUvar: Describe here... -
PEMvar: Describe here... -
TraitVarGraphSim: Describe here...
Author(s)
Guillaume Guenard, with contribution from Pierre Legendre Maintainer: Guillaume Guenard <guillaume.guenard@gmail.com>
References
Butler, M. A. & King, A. A. 2004. Phylogenetic comparative analysis: a modeling approach for adaptive evolution. Am. Nat. 164: 683-695.
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. 4: 1120–1131
Examples
opt <- c(-2,0,2) # Three trait optima: -2, 0, and 2
## Transition probabilities:
transit <- matrix(c(0.7,0.2,0.2,0.2,0.7,0.1,0.1,0.1,0.7),
length(opt),length(opt),dimnames=list(from=opt,to=opt))
## In this example, the trait has a probability of 0.7 to stay at a given
## optimum, a probability of 0.2 for the optimum to change from -2 to 0,
## from 0 to -2, and from 2 to -2, and a probability of 0.1 for the
## optimum to change from -2 to 2, from 0 to 2, and from 2 to 0.
nsp <- 25 # A random tree for 25 species.
tree2 <- rtree(nsp,tip.label=paste("Species",1:nsp,sep=""))
tree2$node.label=paste("N",1:tree2$Nnode,sep="") # Node labels.
## Simulate 10 trials of optimum change.
reg <- EvolveOptimMarkovTree(tp=tree2,tw=transit,p=10,anc=2)
y1 <- TraitOUsimTree(tp=tree2,a=0,sigma=1,
opt=opt[reg[,1]],p=10) ## Neutral
y2 <- TraitOUsimTree(tp=tree2,a=1,sigma=1,
opt=opt[reg[,1]],p=10) ## Few selection.
y3 <- TraitOUsimTree(tp=tree2,a=10,sigma=1,
opt=opt[reg[,1]],p=10) ## Strong selection.
## Display optimum change with colours.
displayOUprocess <- function(tp,trait,regime,mvalue) {
layout(matrix(1:2,1,2))
n <- length(tp$tip.label)
ape::plot.phylo(tp,show.tip.label=TRUE,show.node.label=TRUE,root.edge=FALSE,
direction="rightwards",adj=0,
edge.color=rainbow(length(trait))[regime[tp$edge[,2]]])
plot(y=1:n,x=mvalue[1:n],type="b",xlim=c(-5,5),ylab="",xlab="Trait value",yaxt="n",
bg=rainbow(length(trait))[regime[1:n]],pch=21)
text(trait[regime[1:n]],y=1:n,x=5,col=rainbow(length(trait))[regime[1:n]])
abline(v=0)
}
displayOUprocess(tree2,opt,reg[,1],y1[,1]) # Trait evolve neutrally,
displayOUprocess(tree2,opt,reg[,1],y2[,1]) # under weak selection,
displayOUprocess(tree2,opt,reg[,1],y3[,1]) # under strong selection.
x <- Phylo2DirectedGraph(tree2)
y4 <- TraitVarGraphSim(x, variance = OUvar, p=10, a=5)
DisplayTreeEvol <- function(tp,mvalue) {
layout(matrix(1:2,1,2))
n <- length(tp$tip.label)
ape::plot.phylo(tp,show.tip.label = TRUE, show.node.label = TRUE,
root.edge = FALSE, direction = "rightwards", adj = 0)
plot(y=1:n, x=mvalue[1:n], type="b", xlim=c(-5,5), ylab="",
xlab="Trait value", yaxt="n", pch=21)
abline(v=0)
}
## Recursively displays the simulated traits.
for(i in 1:10) {
DisplayTreeEvol(tree2,y4[i,])
if(is.null(locator(1)))
break ## Stops recursive display on a mouse right-click.
}