Description Usage Arguments Value Author(s) References Examples
Run the Bayesian inference of the clade age-richness index alpha
1 2 3 4 |
tree |
A phylo object |
epsilon |
Minimum size of unsampled splits (see appendix 1) |
beta |
Imbalance index |
niter |
Number of iterations in the mcmc |
ini |
Initial alpha value (default to 0) |
V |
Initial scaling value for the mcmc proposal (default to 0.1) |
chain |
Former mcmc chain (if NULL (the default), a new one is started) |
verbose |
Number of iterations after which the state of the mcmc is printed (if silent == FALSE) |
silent |
If TRUE (the default) the state of the mcmc is not printed |
Nadapt |
Number of iterations between each proposal scalling (default to Inf) |
NadaptMin |
Minimum nmber of iterations before the first proposal scalling (default to 10) |
NadaptMax |
Number of iterations after which the proposal stops being scalled (default to Inf) |
ma |
Minimal alpha value (default to -4) |
Ma |
Maximal alpha value (default to 4) |
proposal |
Shape of the proposal. Can be "bactrian" (the default, ref), "uniform", or "normal" |
accOpt |
Optimal acceptance value (default to 0.3) |
Vmin |
Minimal scaling value for the mcmc proposal (default to 0.001) |
A list with a mcmc field contening the resulting chain. The other fields are only used to resume runing the inference if the chain has to be completed.
Odile Maliet, Fanny Gascuel & Amaury Lambert
Maliet O., Gascuel F., Lambert A. (2018) Ranked tree shapes, non-random extinctions and the loss of phylogenetic diversity, bioRxiv 224295, doi: https://doi.org/10.1101/224295
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 | ntip=30
set.seed(123)
tree=simulate_tree(epsilon = 0.01,alpha = -1,beta = 0,N = ntip,equal.ab = TRUE)
beta=maxlik.betasplit(tree,up=10)$max_lik
plot(tree)
niter=1000
## Not run:
chain=mcmc_alpha(tree,epsilon=0.01,beta=beta,niter=600,V = c(0.1),ini=c(0),
verbose = 100,silent = FALSE,Nadapt = 100,NadaptMin = 100)
# Continue the same chain
chain=mcmc_alpha(tree,epsilon=0.01,beta=beta,niter=400,verbose = 100,silent = FALSE,
chain = chain,Nadapt = 100,NadaptMin = 500,NadaptMax = 700)
thinned=mcmc(chain$mcmc[seq(200,1000,10),])
plot(thinned)
da=density(thinned[,"alpha"])
MPa=da$x[which.max(da$y)]
print(MPa)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.