simulateMRF: Simulates a single Markov random field trajectory.

In RevGadgets: Visualization and Post-Processing of 'RevBayes' Analyses

Simulates a single Markov random field trajectory.

Description

This function simulates a draw from a HSMRF or GMRF distribution given a user-specified global scale parameter. The MRF can be taken to be on the log-scale (such as for a birth rate) or the real-scale. The first value must be specified

Usage

simulateMRF(
  n_episodes,
  model,
  global_scale_hyperprior,
  initial_value = NULL,
  exponentiate = TRUE
)

Arguments

n_episodes	(numeric; no default) The number of episodes in the random field (the parameter vector will be this long).
model	(character; no default) What model should the global scale parameter be set for? Options are "GMRF" and "HSMRF".
global_scale_hyperprior	(numeric; no default) The hyperprior on the global scale parameter.
initial_value	(numeric; NULL) The first value in the MRF. If no value is specified, the field is assumed to start at 0 (if exponentiate=FALSE) or 1 (if exponentiate=TRUE).
exponentiate	(logical; TRUE) If TRUE, the MRF model is taken to be on the log-scale and the values are returned on the real-scale (note this means that the specified initial value will be the log of the true initial value). If FALSE, the model is taken to be on the real scale.

Value

A vector drawn from the specified MRF model on the specified (log- or real-) scale.

References

Magee et al. (2020) Locally adaptive Bayesian birth-death model successfully detects slow and rapid rate shifts. PLoS Computational Biology, 16 (10): e1007999.

Faulkner, James R., and Vladimir N. Minin. Locally adaptive smoothing with Markov random fields and shrinkage priors. Bayesian analysis, 13 (1), 225.

Examples

# Simulate a 100-episode HSMRF model for a speciation-rate through time
trajectory <- simulateMRF(n_episodes = 100,
                          model = "HSMRF",
                          global_scale_hyperprior = 0.0021)
plot(1:100,
     rev(trajectory),
     type = "l",
     xlab = "time",
     ylab = "speciation rate")

RevGadgets documentation built on May 29, 2024, 10:03 a.m.