NAG Library Function Document
nag_gamma_pdf_vector (g01kkc)
1 Purpose
nag_gamma_pdf_vector (g01kkc) returns a number of values of the probability density function (PDF), or its logarithm, for the gamma distribution.
2 Specification
#include <nag.h> |
#include <nagg01.h> |
void |
nag_gamma_pdf_vector (Nag_Boolean ilog,
Integer lx,
const double x[],
Integer la,
const double a[],
Integer lb,
const double b[],
double pdf[],
Integer ivalid[],
NagError *fail) |
|
3 Description
The gamma distribution with shape parameter
and scale parameter
has PDF
If
then an algorithm based directly on the gamma distribution's PDF is used. For values outside this range, the function is calculated via the Poisson distribution's PDF as described in
Loader (2000) (see
Section 9).
The input arrays to this function are designed to allow maximum flexibility in the supply of vector arguments by re-using elements of any arrays that are shorter than the total number of evaluations required. See
Section 2.6 in the g01 Chapter Introduction for further information.
4 References
Loader C (2000) Fast and accurate computation of binomial probabilities (not yet published)
5 Arguments
- 1:
ilog – Nag_BooleanInput
On entry: the value of
ilog determines whether the logarithmic value is returned in
pdf.
- , the probability density function is returned.
- , the logarithm of the probability density function is returned.
- 2:
lx – IntegerInput
On entry: the length of the array
x.
Constraint:
.
- 3:
x[lx] – const doubleInput
On entry: , the values at which the PDF is to be evaluated with , , for .
- 4:
la – IntegerInput
On entry: the length of the array
a.
Constraint:
.
- 5:
a[la] – const doubleInput
On entry: , the shape parameter with , .
Constraint:
, for .
- 6:
lb – IntegerInput
On entry: the length of the array
b.
Constraint:
.
- 7:
b[lb] – const doubleInput
On entry: , the scale parameter with , .
Constraint:
, for .
- 8:
pdf[] – doubleOutput
-
Note: the dimension,
dim, of the array
pdf
must be at least
.
On exit: or .
- 9:
ivalid[] – IntegerOutput
-
Note: the dimension,
dim, of the array
ivalid
must be at least
.
On exit:
indicates any errors with the input arguments, with
- No error.
- .
- .
- overflows, the value returned should be a reasonable approximation.
- 10:
fail – NagError *Input/Output
-
The NAG error argument (see
Section 3.6 in the Essential Introduction).
6 Error Indicators and Warnings
- NE_ARRAY_SIZE
-
On entry, .
Constraint: .
On entry, .
Constraint: .
On entry, .
Constraint: .
- NE_BAD_PARAM
-
On entry, argument had an illegal value.
- NE_INTERNAL_ERROR
-
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact
NAG for assistance.
- NW_IVALID
-
On entry, at least one value of
x,
a or
b was invalid.
Check
ivalid for more information.
7 Accuracy
Not applicable.
8 Parallelism and Performance
Not applicable.
Due to the lack of a stable link to
Loader (2000) paper, we give a brief overview of the method, as applied to the Poisson distribution. The Poisson distribution has a continuous mass function given by,
The usual way of computing this quantity would be to take the logarithm and calculate,
For large and , and are very large, of the same order of magnitude and when calculated have rounding errors. The subtraction of these two terms can therefore result in a number, many orders of magnitude smaller and hence we lose accuracy due to subtraction errors. For example for and , and . But calculated with the method shown later we have . The difference between these two results suggests a loss of about 7 significant figures of precision.
Loader introduces an alternative way of expressing
(1) based on the saddle point expansion,
where
, the deviance for the Poisson distribution is given by,
and
For
close to
,
can be evaluated through the series expansion
otherwise
can be evaluated directly. In addition, Loader suggests evaluating
using the Stirling–De Moivre series,
where the error
is given by
Finally
can be evaluated by combining equations
(1)–
(4) to get,
10 Example
This example prints the value of the gamma distribution PDF at six different points with differing and .
10.1 Program Text
Program Text (g01kkce.c)
10.2 Program Data
Program Data (g01kkce.d)
10.3 Program Results
Program Results (g01kkce.r)