ewenss.sampler: Ewens's sampling formula with optional gamma prior on the...

In wesleycrouse/TIMBR: Tree-based Inference of Multiallelism via Bayesian Regression

Description

Sample allelic series (partitions) from Ewen's sampling formula, optionally informed by user-specified tree(s). Trees must be in coalescent units for appropriate inference.

Usage

ewenss.sampler(samples, trees, prior.alpha, verbose = T)

Arguments

samples	number of samples
trees	either a user-specified tree(s) of class "phylo" ("multiPhylo"), detailed in the 'ape' package, or an integer with the number of leaves to be partitioned
prior.alpha	prior type c("fixed","gamma") for the concentration parameter, see examples for format. type "beta.prime" supported but sampler appears unstable, use ewenss.calc function instead
verbose	optionally report function progress

Value

list of allelic series IDs and probabilities, formatted as prior.M object for TIMBR function

Examples

#specifying hyperparameters for gamma prior using calc.concentration.prior
hyperparam <- calc.concentration.prior(8, 0.05, 0.01)
prior.alpha <- list(type="gamma", shape=hyperparam[1], rate=hyperparam[2])

#running the sampler without user-specified trees; compare with target prior probabilities
prior.M <- ewenss.sampler(100000, 8, prior.alpha)
exp(prior.M$ln.probs[prior.M$M.IDs=="0,0,0,0,0,0,0,0"])
exp(prior.M$ln.probs[prior.M$M.IDs=="0,1,2,3,4,5,6,7"])

#running the sampler with a user-specified tree and fixed concentration parameter; compare with tree structure
tree <- ape::rcoal(8, LETTERS[1:8])
ape::plot.phylo(tree)
prior.alpha <- list(type="fixed", alpha=1)
prior.M <- ewenss.sampler(100000, tree, prior.alpha)
head(prior.M$M.IDs)
head(exp(prior.M$ln.probs))

wesleycrouse/TIMBR documentation built on Feb. 19, 2021, 7:31 a.m.