gsi_mcmc_fb | R Documentation |
Given a list of key parameters from a genetic dataset, this function samples values of pi and the posteriors for all the individuals. Each MCMC iteration includes a recalculation of the scaled genotype likelihood matrix, with baseline allele frequencies updated based on the previous iteration's allocations. It returns the output in a list.
gsi_mcmc_fb(
par_list,
Pi_init,
lambda,
reps,
burn_in,
sample_int_Pi,
sample_int_PofZ
)
par_list |
genetic data converted to the param_list format by |
Pi_init |
Starting value for the pi (collection mixture proportion) vector. |
lambda |
the prior to be added to the collection allocations, in order to generate pseudo-count Dirichlet parameters for the simulation of a new pi vector |
reps |
total number of reps (sweeps) to do. |
burn_in |
how many reps to discard in the beginning when doing the mean calculation. They will still be returned in the traces if desired |
sample_int_Pi |
the number of reps between samples being taken for Pi traces. If 0 no trace samples are taken |
sample_int_PofZ |
the number of reps between samples being taken for the traces of posterior of each individual's origin. If 0 no trace samples are taken. |
gsi_mcmc_fb
returns a list of three. $mean
lists the posterior
means for collection proportions pi
, for the individual posterior
probabilities of assignment PofZ
, and for the allele frequencies theta
.
$sd
returns the posterior standard deviations for the same values.
If the corresponding sample_int
variables are not 0, $trace
contains
samples taken from the Markov chain at intervals of sample_int_
(variable) steps.
# this demonstrates it with scaled likelihoods computed from
# assignment of the reference samples
# we have to get the ploidies to pass to tcf2param_list
locnames <- names(alewife)[-(1:16)][c(TRUE, FALSE)]
ploidies <- rep(2, length(locnames))
names(ploidies) <- locnames
params <- tcf2param_list(alewife, 17, ploidies = ploidies)
lambda <- rep(1/params$C, params$C)
# use very short run and burn in so it doesn't take too long
# when checking on CRAN
mcmc <- gsi_mcmc_fb(params, lambda, lambda, 20, 5, 4, 4)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.