Density function, cumulative distribution function, quantile function
and random number generation for the generalized inverse Gaussian
distribution with parameter vector Theta. Utility routines are
included for the derivative of the density function and to find
suitable break points for use in determining the distribution function.
for x>0, where K_lambda() is the
modified Bessel function of the third kind with order
lambda.
The generalized inverse Gaussian distribution is investigated in
detail in J<c3><b6>rgensen (1982).
Use gigChangePars to convert from the
(delta,gamma),
(alpha,beta), or
(omega,eta) parameterisations to the
(chi,psi), parameterisation used above.
pgig breaks the real line into eight regions in order to
determine the integral of dgig. The break points determining
the regions are found by gigBreaks, based on the values of
small, tiny, and deriv. In the extreme tails of
the distribution where the probability is tiny according to
gigCalcRange, the probability is taken to be zero. For the
generalized inverse Gaussian distribution the leftmost breakpoint is
not affected by the value of tiny but is always taken as 0. In
the inner part of the distribution, the range is divided in 6 regions,
3 above the mode, and 3 below. On each side of the mode, there are two
break points giving the required three regions. The outer break point
is where the probability in the tail has the value given by the
variable small. The inner break point is where the derivative
of the density function is deriv times the maximum value of the
derivative on that side of the mode. In each of the 6 inner regions
the numerical integration routine safeIntegrate (which
is a wrapper for integrate) is used to integrate the
density dgig.
qgig use the breakup of the real line into the same 8
regions as pgig. For quantiles which fall in the 2 extreme
regions, the quantile is returned as -Inf or Inf as
appropriate. In the 6 inner regions splinefun is used to fit
values of the distribution function generated by pgig. The
quantiles are then found using the uniroot function.
pgig and qgig may generally be expected to be accurate
to 5 decimal places. Unfortunately, when lambda is less than
-0.5, the accuracy may be as little as 3 decimal places.
Generalized inverse Gaussian observations are obtained via the
algorithm of Dagpunar (1989).
Value
dgig gives the density, pgig gives the distribution
function, qgig gives the quantile function, and rgig
generates random variates. rgig1 generates random variates in
the special case where lambda=1. An estimate of the
accuracy of the approximation to the distribution function may be
found by setting accuracy = TRUE in the call to phyperb
which then returns a list with components value and
error.
ddgig gives the derivative of dgig.
gigBreaks returns a list with components:
xTiny
Takes value 0 always.
xSmall
Value such that probability to the left is less than
small.
lowBreak
Point to the left of the mode such that the
derivative of the density is deriv times its maximum value
on that side of the mode
highBreak
Point to the right of the mode such that the
derivative of the density is deriv times its maximum value
on that side of the mode
xLarge
Value such that probability to the right is less than
small.
xHuge
Value such that probability to the right is less than
tiny.
modeDist
The mode of the given generalized inverse Gaussian
distribution.
J<c3><b6>rgensen, B. (1982). Statistical Properties of
the Generalized Inverse Gaussian Distribution. Lecture Notes in
Statistics, Vol. 9, Springer-Verlag, New York.
See Also
safeIntegrate, integrate for its
shortfalls, splinefun, uniroot and
gigChangePars for changing parameters to the
(chi,psi) parameterisation, dghyp for
the generalized hyperbolic distribution.
Examples
Theta <- c(1,2,3)
gigRange <- gigCalcRange(Theta, tol = 10^(-3))
par(mfrow = c(1,2))
curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
n = 1000)
title("Density of the\n Generalized Inverse Gaussian")
curve(pgig(x, Theta), from = gigRange[1], to = gigRange[2],
n = 1000)
title("Distribution Function of the\n Generalized Inverse Gaussian")
dataVector <- rgig(500, Theta)
curve(dgig(x, Theta), range(dataVector)[1], range(dataVector)[2],
n = 500)
hist(dataVector, freq = FALSE, add =TRUE)
title("Density and Histogram\n of the Generalized Inverse Gaussian")
logHist(dataVector, main =
"Log-Density and Log-Histogram\n of the Generalized Inverse Gaussian")
curve(log(dgig(x, Theta)), add = TRUE,
range(dataVector)[1], range(dataVector)[2], n = 500)
par(mfrow = c(2,1))
curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
n = 1000)
title("Density of the\n Generalized Inverse Gaussian")
curve(ddgig(x, Theta), from = gigRange[1], to = gigRange[2],
n = 1000)
title("Derivative of the Density\n of the Generalized Inverse Gaussian")
par(mfrow = c(1,1))
gigRange <- gigCalcRange(Theta, tol = 10^(-6))
curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
n = 1000)
bks <- gigBreaks(Theta)
abline(v = bks)
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(HyperbolicDist)
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/HyperbolicDist/dgig.Rd_%03d_medium.png", width=480, height=480)
> ### Name: Generalized Inverse Gaussian
> ### Title: Generalized Inverse Gaussian Distribution
> ### Aliases: dgig pgig qgig rgig rgig1 ddgig gigBreaks
> ### Keywords: distribution
>
> ### ** Examples
>
> Theta <- c(1,2,3)
> gigRange <- gigCalcRange(Theta, tol = 10^(-3))
> par(mfrow = c(1,2))
> curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
+ n = 1000)
> title("Density of the\n Generalized Inverse Gaussian")
> curve(pgig(x, Theta), from = gigRange[1], to = gigRange[2],
+ n = 1000)
> title("Distribution Function of the\n Generalized Inverse Gaussian")
> dataVector <- rgig(500, Theta)
> curve(dgig(x, Theta), range(dataVector)[1], range(dataVector)[2],
+ n = 500)
> hist(dataVector, freq = FALSE, add =TRUE)
> title("Density and Histogram\n of the Generalized Inverse Gaussian")
> logHist(dataVector, main =
+ "Log-Density and Log-Histogram\n of the Generalized Inverse Gaussian")
> curve(log(dgig(x, Theta)), add = TRUE,
+ range(dataVector)[1], range(dataVector)[2], n = 500)
> par(mfrow = c(2,1))
> curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
+ n = 1000)
> title("Density of the\n Generalized Inverse Gaussian")
> curve(ddgig(x, Theta), from = gigRange[1], to = gigRange[2],
+ n = 1000)
> title("Derivative of the Density\n of the Generalized Inverse Gaussian")
> par(mfrow = c(1,1))
> gigRange <- gigCalcRange(Theta, tol = 10^(-6))
> curve(dgig(x, Theta), from = gigRange[1], to = gigRange[2],
+ n = 1000)
> bks <- gigBreaks(Theta)
> abline(v = bks)
>
>
>
>
>
> dev.off()
null device
1
>