JClavel / mvMORPH

mvMORPH: an R package for fitting multivariate evolutionary models to morphometric data
17 stars 5 forks source link

Simulated data has counterintuitive values when alpha is really small #3

Closed meireles closed 6 years ago

meireles commented 6 years ago

The variance in traits simulated under an OU model with a small alpha seems to be really large (see example below). Is mvSIM drawing data from the stationary distribution for the OU model (I think that would explain it)? Anyhow, I wonder what the best way to simulate a dataset that mixes BM and OU traits without running each model separately. Thanks!

library("mvMORPH")
library("phytools")

tree     = phytools::pbtree(n = 50, scale = 1)
trait_mu = c(10, 0.2)
sigma    = c(0.2, 0.0001)
sigma_sq = sigma * sigma

alpha    = matrix(data = diag(c(1e-9, 0.1)),
                  nrow = 2, ncol = 2)

params   = list(ntraits      = 2,
                sigma        = diag(sigma_sq),
                alpha        = alpha,
                theta        = trait_mu)

data_ou   = mvMORPH::mvSIM(tree  = tree,
                           param = params,
                           model = "OU1",
                           nsim  = 10)
JClavel commented 6 years ago

Hi Jose,

Yes the data are drawn by default from the stationary distribution. The stationary distribution of the process you're simulating have a very high variance (about 2e07). A simple way to check it is to use the 'stationary' function:

?stationary class(params)<-c("mvmorph","mvmorph.ou") stationary(params)

You can sample from the process with a fixed root instead:

params = list(ntraits = 2, sigma = diag(sigma_sq), alpha = alpha, theta = trait_mu, vcv = "fixedRoot")

data_ou = mvMORPH::mvSIM(tree = tree, param = params, model = "OU1", nsim = 10)

Regards,

Julien


De : Jose Eduardo Meireles notifications@github.com Envoyé : mardi 13 novembre 2018 05:10 À : JClavel/mvMORPH Cc : Subscribed Objet : [JClavel/mvMORPH] Simulated data is counterintuitive when alpha is really small (#3)

The variance in traits simulated under an OU model with a small alpha seems to be really large (see example below). Is mvSIM drawing data from the stationary distribution for the OU model (I think that would explain it)? Anyhow, I wonder what the best way to simulate a dataset that mixes BM and OU traits without running each model separately. Thanks!

library("mvMORPH") library("phytools")

tree = phytools::pbtree(n = 50, scale = 1) trait_mu = c(10, 0.2) sigma = c(0.2, 0.0001) sigma_sq = sigma * sigma

alpha = matrix(data = diag(c(1e-9, 0.1)), nrow = 2, ncol = 2)

params = list(ntraits = 2, sigma = diag(sigma_sq), alpha = alpha, theta = trait_mu)

data_ou = mvMORPH::mvSIM(tree = tree, param = params, model = "OU1", nsim = 10)

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/JClavel/mvMORPH/issues/3, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AKQkyVOYoS7_y3R9l3UEwPZnEGFn4UObks5uukYegaJpZM4Ya4UA.

meireles commented 6 years ago

@JClavel somehow I had missed the "vcv" argument when reading the doc. Thanks!