Description

Usage

1

Arguments

alpha

scalar, first shape parameter of the Beta density. Must be
greater than 1, see details

beta

scalar, second shape parameter of the Beta density. Must be
greater than 1, see details

p

scalar, content of HPD, must lie between 0 and 1

plot

logical flag, if TRUE then plot the density and
show the HDR

xlim

numeric vector of length 2, the limits of the density's
support to show when plotting; the default is NULL, in which
case the function will confine plotting to where the density is
non-neglible

debug

logical flag, if TRUE produce messages to the
console

Details

The Beta density arises frequently in Bayesian models of
binary events, rates, and proportions, which take on values in the
open unit interval. For instance, the Beta density is a conjugate prior
for the unknown success probability in binomial trials. With shape
parameters α > 1 and β > 1, the Beta density is
unimodal.

In general, suppose θ \in Θ \subseteq R^k
is a random variable with density f(θ). A highest
density region (HDR) of f(θ) with content p \in
(0,1] is a set \mathcal{Q} \subseteq Θ with the
following properties:

\int_\mathcal{Q} f(θ) dθ = p

and

f(θ) > f(θ^*) \, \forall\
θ \in \mathcal{Q},
θ^* \not\in \mathcal{Q}.

For a unimodal
Beta density (the class of Beta densities handled by this function),
a HDR of content 0 < p < 1 is simply an interval \mathcal{Q} \in (0,1).

This function uses numerical methods to solve for the
end points of a HDR for a Beta density with user-specified shape
parameters, via repeated calls to the functions dbeta,
pbeta and qbeta. The function
optimize is used to find points v and w
such that

f(v) = f(w)

subject to the constraint

\int_v^w f(θ; α, β) dθ = p,

where f(θ; α, β) is a Beta density with shape
parameters α and β.

In the special case of α = β > 1, the end points
of a HDR with content p are given by the (1 \pm p)/2
quantiles of the Beta density, and are computed with the
qbeta function.

Again note that the function will only compute a HDR for a unimodal
Beta density, and exit with an error if alpha<=1 | beta <=1.
Note that the uniform density results with α = β = 1,
which does not have a unique HDR with content 0 < p <
1. With shape parameters α<1 and β>1 (or
vice-versa, respectively), the Beta density is infinite at 0 (or 1,
respectively), but still integrates to one, and so a HDR is still
well-defined (but not implemented here, at least not yet).
Similarly, with 0 < α, β < 1 the Beta density is
infinite at both 0 and 1, but integrates to one, and again a HDR of
content p<1 is well-defined in this case, but will be a set of
two disjoint intervals (again, at present, this function does not
cover this case).

Value

If the numerical optimization is successful an vector of length 2,
containing v and w, defined above. If the optimization
fails for whatever reason, a vector of NAs is returned.

The function will also produce a plot of the density with area under
the density supported by the HDR shaded, if the user calls the
function with plot=TRUE; the plot will appear on the current
graphics device.

Debugging messages are printed to the console if the debug
logical flag is set to TRUE.

Author(s)

Simon Jackman [email protected]. Thanks to John
Bullock who discovered a bug in an earlier version.