# qgammaAppr

##### Compute (Approximate) Quantiles of the Gamma Distribution

Compute approximations to the quantile (i.e., inverse cumulative) function of the Gamma distribution.

- Keywords
- distribution

##### Usage

```
qgammaAppr(p, shape, lower.tail = TRUE, log.p = FALSE,
tol = 5e-07)
qgamma.R (p, alpha, scale = 1, lower.tail = TRUE, log.p = FALSE,
EPS1 = 0.01, EPS2 = 5e-07, epsN = 1e-15, maxit = 1000,
pMin = 1e-100, pMax = (1 - 1e-14),
verbose = getOption("verbose"))
```qgammaApprKG(p, shape, lower.tail = TRUE, log.p = FALSE)

qgammaApprSmallP(p, shape, lower.tail = TRUE, log.p = FALSE)

##### Arguments

- p
numeric vector (possibly log tranformed) probabilities.

- shape, alpha
shape parameter, non-negative.

- scale
scale parameter, non-negative, see

`qgamma`

.- lower.tail, log.p
logical, see, e.g.,

`qgamma()`

; must have length 1.- tol
tolerance of maximal approximation error.

- EPS1
small positive number. ...

- EPS2
small positive number. ...

- epsN
small positive number. ...

- maxit
maximal number of iterations. ...

- pMin, pMax
boundaries for

`p`

. ...- verbose
logical indicating if the algorithm should produce “monitoring” information.

##### Details

`qgammaApprSmallP(p, a)`

should be a good approximation in the
following situation when both `p`

and `shape`

\(= \alpha =:
a\) are small :

If we look at Abramowitz&Stegun \(gamma*(a,x) = x^-a * P(a,x)\) and its series \(g*(a,x) = 1/gamma(a) * (1/a - 1/(a+1) * x + ...)\),

then the first order approximation \(P(a,x) = x^a * g*(a,x) ~= x^a/gamma(a+1)\) and hence its inverse \(x = qgamma(p, a) ~= (p * gamma(a+1)) ^ (1/a)\) should be good as soon as \(1/a >> 1/(a+1) * x\)

<==> x << (a+1)/a = (1 + 1/a)

<==> x < eps *(a+1)/a

<==> log(x) < log(eps) + log( (a+1)/a ) = log(eps) + log((a+1)/a) ~ -36 - log(a) where log(x) ~= log(p * gamma(a+1)) / a = (log(p) + lgamma1p(a))/a

such that the above

<==> (log(p) + lgamma1p(a))/a < log(eps) + log((a+1)/a)

<==> log(p) + lgamma1p(a) < a*(-log(a)+ log(eps) + log1p(a))

<==> log(p) < a*(-log(a)+ log(eps) + log1p(a)) - lgamma1p(a) =: bnd(a)

Note that `qgammaApprSmallP()`

indeed also builds on `lgamma1p()`

.

`.qgammaApprBnd(a)`

provides this bound \(bnd(a)\);
it is simply `a*(logEps + log1p(a) - log(a)) - lgamma1p(a)`

, where
`logEps`

is \(\log(\epsilon)\) = `log(eps)`

where ```
eps <-
.Machine$double.eps
```

, i.e. typically (always?) `logEps`

\(= \log
\epsilon = -52 * \log(2) = -36.04365\).

##### Value

numeric

##### References

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.

##### See Also

`qgamma`

for R's Gamma distribution functions.

##### Examples

```
# NOT RUN {
## TODO : Move some of the curve()s from ../tests/qgamma-ex.R !!
# }
```

*Documentation reproduced from package DPQ, version 0.3-3, License: GPL (>= 2)*