| 42333 |
ripley |
1 |
% File src/library/stats/man/GammaDist.Rd
|
| 68948 |
ripley |
2 |
% Part of the R package, https://www.R-project.org
|
| 74265 |
hornik |
3 |
% Copyright 1995-2018 R Core Team
|
| 42333 |
ripley |
4 |
% Distributed under GPL 2 or later
|
|
|
5 |
|
| 27442 |
ripley |
6 |
\name{GammaDist}
|
| 56186 |
murdoch |
7 |
\alias{GammaDist}
|
| 27442 |
ripley |
8 |
\alias{dgamma}
|
|
|
9 |
\alias{pgamma}
|
|
|
10 |
\alias{qgamma}
|
|
|
11 |
\alias{rgamma}
|
| 28897 |
maechler |
12 |
\concept{incomplete gamma function}
|
| 27442 |
ripley |
13 |
\title{The Gamma Distribution}
|
|
|
14 |
\description{
|
|
|
15 |
Density, distribution function, quantile function and random
|
|
|
16 |
generation for the Gamma distribution with parameters \code{shape} and
|
|
|
17 |
\code{scale}.
|
|
|
18 |
}
|
|
|
19 |
\usage{
|
|
|
20 |
dgamma(x, shape, rate = 1, scale = 1/rate, log = FALSE)
|
| 30912 |
ripley |
21 |
pgamma(q, shape, rate = 1, scale = 1/rate, lower.tail = TRUE,
|
|
|
22 |
log.p = FALSE)
|
|
|
23 |
qgamma(p, shape, rate = 1, scale = 1/rate, lower.tail = TRUE,
|
|
|
24 |
log.p = FALSE)
|
| 27442 |
ripley |
25 |
rgamma(n, shape, rate = 1, scale = 1/rate)
|
|
|
26 |
}
|
|
|
27 |
\arguments{
|
|
|
28 |
\item{x, q}{vector of quantiles.}
|
|
|
29 |
\item{p}{vector of probabilities.}
|
|
|
30 |
\item{n}{number of observations. If \code{length(n) > 1}, the length
|
|
|
31 |
is taken to be the number required.}
|
|
|
32 |
\item{rate}{an alternative way to specify the scale.}
|
| 43895 |
maechler |
33 |
\item{shape, scale}{shape and scale parameters. Must be positive,
|
|
|
34 |
\code{scale} strictly.}
|
| 36220 |
ripley |
35 |
\item{log, log.p}{logical; if \code{TRUE}, probabilities/densities \eqn{p}
|
|
|
36 |
are returned as \eqn{log(p)}.}
|
| 27442 |
ripley |
37 |
\item{lower.tail}{logical; if TRUE (default), probabilities are
|
| 50805 |
ripley |
38 |
\eqn{P[X \le x]}, otherwise, \eqn{P[X > x]}.}
|
| 27442 |
ripley |
39 |
}
|
|
|
40 |
\value{
|
|
|
41 |
\code{dgamma} gives the density,
|
| 34848 |
ripley |
42 |
\code{pgamma} gives the distribution function,
|
| 27442 |
ripley |
43 |
\code{qgamma} gives the quantile function, and
|
|
|
44 |
\code{rgamma} generates random deviates.
|
| 37274 |
ripley |
45 |
|
|
|
46 |
Invalid arguments will result in return value \code{NaN}, with a warning.
|
| 61975 |
murdoch |
47 |
|
|
|
48 |
The length of the result is determined by \code{n} for
|
|
|
49 |
\code{rgamma}, and is the maximum of the lengths of the
|
| 65782 |
ripley |
50 |
numerical arguments for the other functions.
|
| 61975 |
murdoch |
51 |
|
| 65782 |
ripley |
52 |
The numerical arguments other than \code{n} are recycled to the
|
| 61975 |
murdoch |
53 |
length of the result. Only the first elements of the logical
|
| 65782 |
ripley |
54 |
arguments are used.
|
| 27442 |
ripley |
55 |
}
|
|
|
56 |
\details{
|
|
|
57 |
If \code{scale} is omitted, it assumes the default value of \code{1}.
|
|
|
58 |
|
|
|
59 |
The Gamma distribution with parameters \code{shape} \eqn{=\alpha}{= a}
|
|
|
60 |
and \code{scale} \eqn{=\sigma}{= s} has density
|
|
|
61 |
\deqn{
|
|
|
62 |
f(x)= \frac{1}{{\sigma}^{\alpha}\Gamma(\alpha)} {x}^{\alpha-1} e^{-x/\sigma}%
|
|
|
63 |
}{f(x)= 1/(s^a Gamma(a)) x^(a-1) e^-(x/s)}
|
| 50805 |
ripley |
64 |
for \eqn{x \ge 0}, \eqn{\alpha > 0}{a > 0} and \eqn{\sigma > 0}{s > 0}.
|
| 34851 |
ripley |
65 |
(Here \eqn{\Gamma(\alpha)}{Gamma(a)} is the function implemented by \R's
|
| 59979 |
ripley |
66 |
\code{\link{gamma}()} and defined in its help. Note that \eqn{a = 0}
|
| 43899 |
maechler |
67 |
corresponds to the trivial distribution with all mass at point 0.)
|
| 34848 |
ripley |
68 |
|
| 27442 |
ripley |
69 |
The mean and variance are
|
|
|
70 |
\eqn{E(X) = \alpha\sigma}{E(X) = a*s} and
|
|
|
71 |
\eqn{Var(X) = \alpha\sigma^2}{Var(X) = a*s^2}.
|
| 28231 |
maechler |
72 |
|
| 34848 |
ripley |
73 |
The cumulative hazard \eqn{H(t) = - \log(1 - F(t))}{H(t) = - log(1 - F(t))}
|
| 61173 |
ripley |
74 |
is
|
|
|
75 |
\preformatted{-pgamma(t, ..., lower = FALSE, log = TRUE)
|
|
|
76 |
}
|
| 39015 |
ripley |
77 |
|
| 46325 |
ripley |
78 |
Note that for smallish values of \code{shape} (and moderate
|
|
|
79 |
\code{scale}) a large parts of the mass of the Gamma distribution is
|
|
|
80 |
on values of \eqn{x} so near zero that they will be represented as
|
| 61177 |
ripley |
81 |
zero in computer arithmetic. So \code{rgamma} may well return values
|
| 46325 |
ripley |
82 |
which will be represented as zero. (This will also happen for very
|
|
|
83 |
large values of \code{scale} since the actual generation is done for
|
| 61160 |
ripley |
84 |
\code{scale = 1}.)
|
| 50414 |
ripley |
85 |
}
|
| 50399 |
maechler |
86 |
% Have caught all currently known problems; hence no longer say:
|
|
|
87 |
% Similarly, \code{qgamma} has a very hard job for
|
|
|
88 |
% small \code{scale}, and warns of potential unreliability for
|
|
|
89 |
% \code{scale < 1e-10}.
|
| 27442 |
ripley |
90 |
\note{
|
| 74363 |
maechler |
91 |
The S (Becker \emph{et al}, 1988) parametrization was via \code{shape}
|
| 71883 |
ripley |
92 |
and \code{rate}: S had no \code{scale} parameter. It is an error
|
|
|
93 |
to supply and \code{scale} and \code{rate}.
|
| 28231 |
maechler |
94 |
|
| 28789 |
ripley |
95 |
\code{pgamma} is closely related to the incomplete gamma function. As
|
| 42706 |
ripley |
96 |
defined by Abramowitz and Stegun 6.5.1 (and by \sQuote{Numerical
|
|
|
97 |
Recipes}) this is
|
| 61177 |
ripley |
98 |
\deqn{P(a,x) = \frac{1}{\Gamma(a)} \int_0^x t^{a-1} e^{-t} dt}{P(a,x) = 1/Gamma(a) integral_0^x t^(a-1) exp(-t) dt}
|
| 28789 |
ripley |
99 |
\eqn{P(a, x)} is \code{pgamma(x, a)}. Other authors (for example
|
|
|
100 |
Karl Pearson in his 1922 tables) omit the normalizing factor,
|
| 61042 |
maechler |
101 |
defining the incomplete gamma function \eqn{\gamma(a,x)} as
|
|
|
102 |
\eqn{\gamma(a,x) = \int_0^x t^{a-1} e^{-t} dt,}{gamma(a,x) =
|
|
|
103 |
integral_0^x t^(a-1) exp(-t) dt,} i.e., \code{pgamma(x, a) * gamma(a)}.
|
|
|
104 |
Yet other use the \sQuote{upper} incomplete gamma function,
|
| 61177 |
ripley |
105 |
\deqn{\Gamma(a,x) = \int_x^\infty t^{a-1} e^{-t} dt,}{Gamma(a,x) = integral_x^Inf t^(a-1) exp(-t) dt,}
|
| 61042 |
maechler |
106 |
which can be computed by
|
| 61160 |
ripley |
107 |
\code{pgamma(x, a, lower = FALSE) * gamma(a)}.
|
| 61042 |
maechler |
108 |
|
|
|
109 |
Note however that \code{pgamma(x, a, ..)} currently requires \eqn{a > 0},
|
|
|
110 |
whereas the incomplete gamma function is also defined for negative
|
|
|
111 |
\eqn{a}. In that case, you can use \code{gamma_inc(a,x)} (for
|
| 62120 |
hornik |
112 |
\eqn{\Gamma(a,x)}) from package \CRANpkg{gsl}.
|
| 61042 |
maechler |
113 |
|
|
|
114 |
See also
|
| 68939 |
ripley |
115 |
\url{https://en.wikipedia.org/wiki/Incomplete_gamma_function}, or
|
| 61042 |
maechler |
116 |
\url{http://dlmf.nist.gov/8.2#i}.
|
| 37222 |
ripley |
117 |
}
|
|
|
118 |
\source{
|
| 37226 |
ripley |
119 |
\code{dgamma} is computed via the Poisson density, using code contributed
|
| 37222 |
ripley |
120 |
by Catherine Loader (see \code{\link{dbinom}}).
|
| 43895 |
maechler |
121 |
|
| 42987 |
ripley |
122 |
\code{pgamma} uses an unpublished (and not otherwise documented)
|
|
|
123 |
algorithm \sQuote{mainly by Morten Welinder}.
|
| 34848 |
ripley |
124 |
|
| 37222 |
ripley |
125 |
\code{qgamma} is based on a C translation of
|
|
|
126 |
|
|
|
127 |
Best, D. J. and D. E. Roberts (1975).
|
|
|
128 |
Algorithm AS91. Percentage points of the chi-squared distribution.
|
|
|
129 |
\emph{Applied Statistics}, \bold{24}, 385--388.
|
|
|
130 |
|
|
|
131 |
plus a final Newton step to improve the approximation.
|
| 43895 |
maechler |
132 |
|
| 37222 |
ripley |
133 |
\code{rgamma} for \code{shape >= 1} uses
|
|
|
134 |
|
|
|
135 |
Ahrens, J. H. and Dieter, U. (1982).
|
|
|
136 |
Generating gamma variates by a modified rejection technique.
|
|
|
137 |
\emph{Communications of the ACM}, \bold{25}, 47--54,
|
|
|
138 |
|
|
|
139 |
and for \code{0 < shape < 1} uses
|
|
|
140 |
|
|
|
141 |
Ahrens, J. H. and Dieter, U. (1974).
|
|
|
142 |
Computer methods for sampling from gamma, beta, Poisson and binomial
|
|
|
143 |
distributions. \emph{Computing}, \bold{12}, 223--246.
|
| 27442 |
ripley |
144 |
}
|
|
|
145 |
\references{
|
| 74265 |
hornik |
146 |
Becker, R. A., Chambers, J. M. and Wilks, A. R. (1988).
|
| 27442 |
ripley |
147 |
\emph{The New S Language}.
|
| 47262 |
ripley |
148 |
Wadsworth & Brooks/Cole.
|
| 28231 |
maechler |
149 |
|
| 74265 |
hornik |
150 |
Shea, B. L. (1988).
|
|
|
151 |
Algorithm AS 239: Chi-squared and incomplete Gamma integral,
|
|
|
152 |
\emph{Applied Statistics (JRSS C)}, \bold{37}, 466--473.
|
|
|
153 |
\doi{10.2307/2347328}.
|
| 28789 |
ripley |
154 |
|
|
|
155 |
Abramowitz, M. and Stegun, I. A. (1972)
|
|
|
156 |
\emph{Handbook of Mathematical Functions.} New York: Dover.
|
|
|
157 |
Chapter 6: Gamma and Related Functions.
|
| 61042 |
maechler |
158 |
|
|
|
159 |
NIST Digital Library of Mathematical Functions.
|
|
|
160 |
\url{http://dlmf.nist.gov/}, section 8.2.
|
| 27442 |
ripley |
161 |
}
|
|
|
162 |
\seealso{
|
| 52770 |
ripley |
163 |
\code{\link{gamma}} for the gamma function.
|
|
|
164 |
|
|
|
165 |
\link{Distributions} for other standard distributions, including
|
|
|
166 |
\code{\link{dbeta}} for the Beta distribution and \code{\link{dchisq}}
|
|
|
167 |
for the chi-squared distribution which is a special case of the Gamma
|
|
|
168 |
distribution.
|
| 27442 |
ripley |
169 |
}
|
|
|
170 |
\examples{
|
| 61160 |
ripley |
171 |
-log(dgamma(1:4, shape = 1))
|
| 27442 |
ripley |
172 |
p <- (1:9)/10
|
| 61168 |
ripley |
173 |
pgamma(qgamma(p, shape = 2), shape = 2)
|
| 61160 |
ripley |
174 |
1 - 1/exp(qgamma(p, shape = 1))
|
| 39015 |
ripley |
175 |
|
| 55251 |
ripley |
176 |
\donttest{# even for shape = 0.001 about half the mass is on numbers
|
| 39015 |
ripley |
177 |
# that cannot be represented accurately (and most of those as zero)
|
|
|
178 |
pgamma(.Machine$double.xmin, 0.001)
|
| 52417 |
maechler |
179 |
pgamma(5e-324, 0.001) # on most machines 5e-324 is the smallest
|
| 39015 |
ripley |
180 |
# representable non-zero number
|
|
|
181 |
table(rgamma(1e4, 0.001) == 0)/1e4
|
| 55251 |
ripley |
182 |
}}
|
| 27442 |
ripley |
183 |
\keyword{distribution}
|