ratematrix (version 1.2)

samplePrior: Take samples from the prior distribution

Description

Take samples from the prior distribution.

Usage

samplePrior(n, prior, sample.sd = TRUE, rebuild.R = FALSE)

Arguments

n

number of samples to be generated.

prior

the object with the prior function. See 'makePrior' for more information.

sample.sd

whether the function should sample the vector of standard deviations independently from the correlation matrices. See 'Details'.

rebuild.R

whether the prior sample should return an evolutionary rate matrix rather than a correlation matrix and a vector of standard deviations (default is FALSE). See 'Details'.

Value

A list with samples from the prior distribution. The structure of this list is the same as required by the parameter 'start' of the 'ratematrixMCMC'.

Details

The prior samples from this function can be used to start the MCMC sampler. See the examples below.

If 'sample.sd' is set to FALSE the samples from the standard deviations will be derived from the covariance matrices. If 'sample.sd' is set to TRUE (default) then standard deviations are independently sampled from their own prior distribution and are not derived from the samples of the correlation matrix. Option 'sample.sd = TRUE' is the one used during the MCMC.

The option 'rebuild.R' controls if the samples from the posterior distribution should return the standard deviation separated from the correlation matrix or if these elements should be used to rebuild the covariance matrix. Set 'rebuild.R' to TRUE if you want to obtain the covariance matrices. Otherwise, the 'plotPrior' function works better when 'rebuild.R' is set to FALSE.

Examples

# NOT RUN {
data( centrarchidae )
dt.range <- t( apply( centrarchidae$data, 2, range ) )
## The step size for the root value can be set given the range we need to sample from:
w_mu <- ( dt.range[,2] - dt.range[,1] ) / 10
par.sd <- cbind(c(0,0), sqrt( c(10,10) ))
prior <- makePrior(r=2, p=2, den.mu="unif", par.mu=dt.range, den.sd="unif", par.sd=par.sd)
prior.samples <- samplePrior(n = 1000, prior = prior)
start.point <- samplePrior(n=1, prior=prior)
## Plot the prior. Red line shows the sample from the prior that will set the starting 
##        point for the MCMC.
plotRatematrix(prior.samples, point.matrix = start.point$matrix, point.color = "red"
               , point.wd = 2)
handle <- ratematrixMCMC(data=centrarchidae$data, phy=centrarchidae$phy.map, prior=prior
                         , gen=10000, w_mu=w_mu, dir=tempdir())
posterior <- readMCMC(handle, burn = 0.2, thin = 10)
## Again, here the red line shows the starting point of the MCMC.
plotRatematrix( posterior, point.matrix = start.point$matrix, point.color = "red"
               , point.wd = 2)
# }