DPQmpfr (version 0.3-3)

pbeta_ser: Beta Distribution Function -- ‘BPSER’ Series Expansion from TOMS 708

Description

Compute a version of the Beta cumulative distribution function (pbeta() in R), namely using the series expansion, named BPSER(), from “TOMS 708”, i.e., Didonato and Morris (1992).

This “pure R” function exists for didactical or documentational reasons on one hand, as R's own pbeta() uses this expansion when appropriate and other algorithms otherwise. On the other hand, using high precision q and MPFR arithmetic (via package Rmpfr) may allow to get highly accurate pbeta() values.

Usage

pbeta_ser(q, shape1, shape2, log.p = FALSE, eps = 1e-15, errPb = 0, verbose = FALSE)

Value

An approximation to the Beta probability \(P[X \le q]\)

for \( X \sim B(a,b),\) (where \(a=\)shape1, and \(b=\)shape2).

Arguments

q, shape1, shape2

quantiles and shape parameters of the Beta distribution, q typically in \([0,1]\), see pbeta. Here, q must be scalar, i.e., of length one, and may inherit from class "mpfr", in order to be more accurate (than with the double precision computations).

log.p

if TRUE, probabilities p are given as log(p).

eps

non-negative number; tol <- eps/shape1 will be used for convergence checks in the series computations.

errPb

an integer code, typically in -2, -1, 0 to determine how warnings on convergence failures are handled.

verbose

logical indicating if console output about intermediate results should be printed.

Author

Didonato and Morris and R Core team; separate packaging by Martin Maechler.

Details

pbeta_ser() crucially needs three auxiliary functions which we “mpfr-ized” as well: gam1M(), lgamma1pM(), and algdivM.

References

Didonato, A. and Morris, A., Jr, (1992) Algorithm 708: Significant digit computation of the incomplete beta function ratios, ACM Transactions on Mathematical Software 18, 360--373; tools:::Rd_expr_doi("10.1145/131766.131776").

See Also

pbeta, DPQmpfr's own pbetaD94; even more pbeta() approximations in package DPQ, e.g., pnbetaAS310, or pbetaRv1.

In addition, for integer shape parameters, the potentially “fully accurate” finite sum base pbetaI() in package Rmpfr.

Examples

(p. <- pbeta_ser(1/2, shape1 = 2, shape2 = 3, verbose=TRUE))
(lp <- pbeta_ser(1/2, shape1 = 2, shape2 = 3, log.p = TRUE))
          all.equal(lp, log(p.), tolerance=0) # 1.48e-16
stopifnot(all.equal(lp, log(p.), tolerance = 1e-13))

## Using  Vectorize() in order to allow vector 'q' e.g. for curve():
str(pbetaSer <- Vectorize(pbeta_ser, "q"))
curve(pbetaSer(x, 1.5, 4.5)); abline(h=0:1, v=0:1, lty=2, col="gray")
curve(pbeta   (x, 1.5, 4.5), add=TRUE, col = adjustcolor(2, 1/4), lwd=3)

## now using mpfr-numbers:
half <- 1/Rmpfr::mpfr(2, 256)
(p2 <- pbeta_ser(half, shape1 = 1, shape2 = 123))