It is possible to fix p to a given value
describing the location of the distribution, that is, how close to the
minimum age is the mode in the truncated Cauchy, and then estimate
c so that the area under the curve attains a certain
maximum age, allowing for some probability to be allocated to the right
of it. The function c_truncauchy can be used provided that
we have both age minimum and maximum, and define some probability
pr to be allocated to the right of the maximum age.
In this example, our minimum age is e.g. \(1\) in units of \(100 Ma\), and our maximum age is \(4.93\) also in units of \(100 Ma\). The minimum is a soft bound
allowing \(0.025\) of the density to be
allocated to the left of it, whereas \(0.975\) is the percentile at which we
observe the maximum age. The quantity p=0.001 has been
chosen so that the mode is closer to the minimum age.
# load the package
library(tbea)
# estimate the c parameter for the L distribution
cparam <- c_truncauchy(tl=0.4094, tr=0.4160, p=0.001, pr=0.975, al=0.025)
cparam## [1] 0.0009289821
The function gives us an estimated of approximately \(2.0\) for c. Thus we can
specify this distribution in MCMCTree as
L(0.4094, 0.1, 0.0009289821, 0.025). We can use the
packages mcmc3r 1 to plot the L density. The
following code should do the trick.
# load the package
library(mcmc3r)
# using the function dL to plot the L density
curve(dL(x, tL=0.4094, p=0.001, c=cparam), from=0.4094, to=0.4160)dos Reis, M. et al. (2018). Using phylogenomic data to explore the effects of relaxed clocks and calibration strategies on divergence time estimation: Primates as a test case. Systematic Biology, 67(4):594-615.↩︎