PoissonInverseGaussian {actuar}R Documentation

The Poisson-Inverse Gaussian Distribution

Description

Density function, distribution function, quantile function and random generation for the Poisson-inverse Gaussian discrete distribution with parameters mean and shape.

Usage

dpoisinvgauss(x, mean, shape = 1, dispersion = 1/shape,
              log = FALSE)
ppoisinvgauss(q, mean, shape = 1, dispersion = 1/shape,
              lower.tail = TRUE, log.p = FALSE)
qpoisinvgauss(p, mean, shape = 1, dispersion = 1/shape,
              lower.tail = TRUE, log.p = FALSE)
rpoisinvgauss(n, mean, shape = 1, dispersion = 1/shape)

Arguments

x

vector of (positive integer) quantiles.

q

vector of quantiles.

p

vector of probabilities.

n

number of observations. If length(n) > 1, the length is taken to be the number required.

mean, shape

parameters. Must be strictly positive. Infinite values are supported.

dispersion

an alternative way to specify the shape.

log, log.p

logical; if TRUE, probabilities p are returned as log(p).

lower.tail

logical; if TRUE (default), probabilities are P[X ≤ x], otherwise, P[X > x].

Details

The Poisson-inverse Gaussian distribution is the result of the continuous mixture between a Poisson distribution and an inverse Gaussian, that is, the distribution with probability mass function

p(x) = int_0^Inf (y^x exp(-y))/x! g(y; μ, φ) dy,

where g(y; μ, φ) is the density function of the inverse Gaussian distribution with parameters mean = μ and dispersion = φ (see dinvgauss).

The resulting probability mass function is

p(x) = sqrt(2/(π φ)) exp(1/(φ μ))/x! * [√(2 φ (1 + 1/(2 φ μ^2)))]^(-(x-1/2)) * K(√((2/φ) (1 + 1/(2 φ μ^2))); x-1/2),

for x = 0, 1, …, μ > 0, φ > 0 and where K(x; ν) is the modified Bessel function of the third kind implemented by R's besselK() and defined in its help.

The limiting case μ = Inf has well defined probability mass and distribution functions, but has no finite strictly positive, integer moments. The pmf in this case reduces to

p(x) = sqrt(2/(π φ)) 1/x! [√(2 φ)]^(-(x-1/2)) * K(√(2/φ); x-1/2).

The limiting case φ = 0 is a degenerate distribution in x = 0.

If an element of x is not integer, the result of dpoisinvgauss is zero, with a warning.

The quantile is defined as the smallest value x such that F(x) ≥ p, where F is the distribution function.

Value

dpoisinvgauss gives the probability mass function, ppoisinvgauss gives the distribution function, qpoisinvgauss gives the quantile function, and rpoisinvgauss generates random deviates.

Invalid arguments will result in return value NaN, with a warning.

The length of the result is determined by n for rpoisinvgauss, and is the maximum of the lengths of the numerical arguments for the other functions.

Note

[dpqr]pig are aliases for [dpqr]poisinvgauss.

qpoisinvgauss is based on qbinom et al.; it uses the Cornish–Fisher Expansion to include a skewness correction to a normal approximation, followed by a search.

Author(s)

Vincent Goulet vincent.goulet@act.ulaval.ca

References

Holla, M. S. (1966), “On a Poisson-Inverse Gaussian Distribution”, Metrika, vol. 15, p. 377-384.

Johnson, N. L., Kemp, A. W. and Kotz, S. (2005), Univariate Discrete Distributions, Third Edition, Wiley.

Klugman, S. A., Panjer, H. H. and Willmot, G. E. (2012), Loss Models, From Data to Decisions, Fourth Edition, Wiley.

Shaban, S. A., (1981) “Computation of the poisson-inverse gaussian distribution”, Communications in Statistics - Theory and Methods, vol. 10, no. 14, p. 1389-1399.

See Also

dpois for the Poisson distribution, dinvgauss for the inverse Gaussian distribution.

Examples

## Tables I and II of Shaban (1981)
x <- 0:2
sapply(c(0.4, 0.8, 1), dpoisinvgauss, x = x, mean = 0.1)
sapply(c(40, 80, 100, 130), dpoisinvgauss, x = x, mean = 1)

qpoisinvgauss(ppoisinvgauss(0:10, 1, dis = 2.5), 1, dis = 2.5)

x <- rpoisinvgauss(1000, 1, dis = 2.5)
y <- sort(unique(x))
plot(y, table(x)/length(x), type = "h", lwd = 2,
     pch = 19, col = "black", xlab = "x", ylab = "p(x)",
     main = "Empirical vs theoretical probabilities")
points(y, dpoisinvgauss(y, 1, dis = 2.5),
       pch = 19, col = "red")
legend("topright", c("empirical", "theoretical"),
       lty = c(1, NA), pch = c(NA, 19), col = c("black", "red"))

[Package actuar version 3.1-4 Index]