Hi Emmanuel,
Great point that it will be faster to simulate from the root up rather than
compute the matrix for large trees. In that case, how about simulate from
the equation
x(t'') = x(t') e^{-alpha (t''-t') }+theta (1-e^{-alpha (t''-t') }) + sqrt(
sigma^2( 1 -e^{-2\alpha (t''-t') )/2 alpha
Dear list,
The simulation of OU processes using rTraitCont seems to return
erroneous results for non-small values of alpha:
Consider this example:
require(geiger)
require(TreeSim)
tree - sim.bd.taxa(100,1,1,0)[[1]]
x - rTraitCont(tree, model=OU, sigma=1, alpha=10, theta=0, root.value=0)