Compute log(beta(a,b))
in a simple (fast) or asymptotic way.
lbetaM (a, b, k.max = 5, give.all = FALSE)
lbeta_asy(a, b, k.max = 5, give.all = FALSE)
lbetaMM (a, b, cutAsy = 1e-2, verbose = FALSE) betaI(a, n)
lbetaI(a, n)
logQab_asy(a, b, k.max = 5, give.all = FALSE)
Qab_terms(a, k)
a fast or simple (approximate) computation of lbeta(a,b)
.
the Beta parameters, see beta
; n
must
be a positive integer and “small”.
..
logical
..
cutoff value from where to switch to asymptotic formula.
logical (or integer) indicating if and how much monitoring information should be printed to the console.
the number of terms in the series expansion of
Qab_terms()
, currently must be in \({0, 1, .., 5}\).
Martin Maechler
All lbeta*()
functions compute log(beta(a,b))
.
We use \(Qab = Qab(a,b)\) for
$$Q_{a,b} := \frac{\Gamma(a + b)}{\Gamma(b)},$$
which is numerically challenging when \(b\) becomes large compared to
a
, or \(a \ll b\).
With the beta function
$$B(a,b) = \frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)} = \frac{\Gamma(a)}{Qab},$$
and hence
$$\log B(a,b) = \log\Gamma(a) + \log\Gamma(b) - \log\Gamma(a+b) = \log\Gamma(a) - \log Qab,$$
or in R, lBeta(a,b) := lgamma(a) - logQab(a,b)
.
Indeed, typically everything has to be computed in log scale, as both \(\Gamma(b)\)
and \(\Gamma(a+b)\) would overflow numerically for large \(b\).
Consequently, we use logQab*()
, and for the large \(b\) case
logQab_asy()
specifically,
$$\code{logQab(a,b)} := \log( Qab(a,b) ).$$
Note this is related to trying to get asymptotic formula for \(\Gamma\) ratios, notably formula (6.1.47) in Abramowitz and Stegun.
Note how this is related to computing qbeta()
in boundary
cases, and see algdiv()
‘Details’ about this.
We also have a vignette
about this, but really the problem has been adressed pragmatically
by the authors of TOMS 708, see the ‘References’ in
pbeta
,
by their routine algdiv()
which also is available in our
package DPQ.
Abramowitz, M. and Stegun, I. A. (1972)
Handbook of Mathematical Functions. New York: Dover.
https://en.wikipedia.org/wiki/Abramowitz_and_Stegun provides
links to the full text which is in public domain.
Formula (6.1.47), p.257
R's beta
function; algdiv()
.