View source: R/simulate_rate.R
simulate_rate | R Documentation |
simulate_rate
Simulates three different evolutionary rates models. In
the two first, 'predictor_BM' and 'predictor_gBM', the evolution of y follows
a Brownian motion with variance linear in x, while the evolution of x either
follows a Brownian motion or a geometric Brownian motion, respectively. In
the third model, 'recent_evol', the residuals of the macroevolutionary
predictions of y have variance linear in x.
simulate_rate( tree, startv_x = NULL, sigma_x = NULL, a, b, sigma_y = NULL, x = NULL, model = "predictor_BM" )
tree |
A |
startv_x |
The starting value for x (usually 0 for 'predictor_BM' and 'recent_evol', and 1 for 'predictor_gBM'). |
sigma_x |
The evolutionary rate of x. |
a |
A parameter of the evolutionary rate of y. |
b |
A parameter of the evolutionary rate of y. |
sigma_y |
The evolutionary rate for the macroevolution of y (Brownian motion process) used in the 'recent_evolution' model. |
x |
Optional fixed values of x (only for the 'recent_evol' model), must equal number of tips in the phylogeny, must correspond to the order of the tip labels. |
model |
Either a Brownian motion model of x 'predictor_BM', geometric Brownian motion model of x 'predictor_gBM', or 'recent_evol'. |
See the vignette 'Analyzing rates of evolution' for an explanation
of the evolutionary models. The data of the tips can be analyzed with the
function rate_gls
. Note that a large part of parameter space
will cause negative roots in the rates of y (i.e. negative a+bx). In these
cases the rates are set to 0. A warning message is given if the number of
such instances is larger than 0.1%. For model 1 and 2, it is possible to
set ‘a=’scaleme'', if this chosen then 'a' will be given the lowest
possible value constrained by a+bx>0.
An object of class
'simulate_rate'
, which is a list
with the following components:
tips
| A data frame of x and y values for the tips. | |||
percent_negative_roots | The percent of iterations with negative roots in the rates of y (not given for model = 'recent_evol'). | |||
compl_dynamics | A list with the output of the complete dynamics (not given for model = 'recent_evol'). |
Geir H. Bolstad
# Also see the vignette 'Analyzing rates of evolution'. ## Not run: # Generating a tree with 50 species set.seed(102) tree <- ape::rtree(n = 50) tree <- ape::chronopl(tree, lambda = 1, age.min = 1) ### model = 'predictor_BM' ### sim_data <- simulate_rate(tree, startv_x = 0, sigma_x = 0.25, a = 1, b = 1, model = "predictor_BM") head(sim_data$tips) par(mfrow = c(1, 3)) plot(sim_data) plot(sim_data, response = "y") plot(sim_data, response = "x") par(mfrow = c(1, 1)) ### model = 'predictor_gBM' ### sim_data <- simulate_rate(tree, startv_x = 1, sigma_x = 1, a = 1, b = 0.1, model = "predictor_gBM") head(sim_data$tips) par(mfrow = c(1, 3)) plot(sim_data) plot(sim_data, response = "y") plot(sim_data, response = "x") par(mfrow = c(1, 1)) ### model = 'recent_evol' ### sim_data <- simulate_rate(tree, startv_x = 0, sigma_x = 1, a = 1, b = 1, sigma_y = 1, model = "recent_evol" ) head(sim_data$tips) ## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.