3
$\begingroup$

I have the mean and standard deviation of my data, which I determined follows a gamma distribution. I don't understand the function I found online for the CDF of a gamma distribution because of the gamma functions nested within. How do I write the CDF?

Also, do I have to take the log of mean and standard deviation for a gamma distribution like I would when representing the data with a lognormal distribution?

I'm doing all of this in R, so not sure if there's a simpler way to do this in there.

$\endgroup$
1
  • $\begingroup$ If you know the mean and standard deviation of the gamma distribution, then you can use method-of-moments estimators for the distribution's parameters. Any statistical package will allow you to compute the CDF of a gamma distribution given its parameters. $\endgroup$
    – Sycorax
    Commented Jul 11, 2019 at 15:41

1 Answer 1

4
$\begingroup$

Example: Let's look at a specific example: With R, I take a random sample of size $n = 100$ from $$\mathsf{Gamma}(\text{shape}=\alpha=3, \text{scale}=\theta=10)\\ \equiv \mathsf{Gamma}(\text{shape}=\alpha=3, \text{rate}=\lambda=0.1).$$ rounding to three places. (R uses the latter parameterization, in terms of shape and rate.)

set.seed(711)
x = round(rgamma(100, 3, .1), 3)
summary(x); var(x)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  4.078  16.839  28.058  30.534  38.702  99.846 
[1] 336.9411

Method of Moments. Because the mean of this distribution is $\mu = \alpha\theta$ and its variance is $\sigma^2 = \alpha\theta^2,$ you can get a method-of-moments estimator (MME) of $\theta$ as $\check\theta = 336.9411/30.534 = 11.03495.$ So $\check\lambda = 1/\check\theta - 0.09062118$ and $\check\alpha = \bar X/\check\theta = 30.534/11.03495 = 2.767027.$ (As suggested in the comment by @Sycorax.)

These estimates are not far from the parameter values $\lambda = 0.1$ and $\alpha = 3$ of the population.

CDFs and Gamma functions. No important complications arise from use of the gamma function $\Gamma(\cdot)$ in the constant of integration of gamma distributions.

For positive integer arguments the Gamma function has $\Gamma(k) = (n-1)!.$ In applications, many uses of the gamma distribution have integer values of the shape parameter $\alpha.$ For example $\mathsf{Gamma}(3, .1)$ is the distribution of the sum of three independent random variables $X_i \sim \mathsf{Exp}(\lambda = 0.1).$ However, the Gamma function is also defined for positive real $k$ (and, irrelevant here, for noninteger negative values of $k).$

gamma(5); factorial(4)
[1] 24
[1] 24
gamma(1/2);  sqrt(pi)
[1] 1.772454
[1] 1.772454

The use of the incomplete gamma function $\gamma$ in the CDF, indicates that the CDF is not available in closed form for all choices of parameters. Similarly, the CDF of the normal distribution is not available in closed form for any choice of parameters. R and other statistical software provide values of the gamma CDF (in R, as pgamma) just as they do for the normal CDF.

Applications. I don't know what applications you might have in mind for the CDF after you estimate it. If you want the 80th percentile $(42.8)$ of the population $\mathsf{Gamma}(3, 0.1),$ you can get it with qgamma (inverse CDF) using the population parameters (generally unknown in a real application). If you want to estimate the 80th percentile of the population, using the sample, you can use quantile to find the sample 80th percentile $(43.5)$ or you can use qgamma with the estimated parameters to get $44.0.$

qgamma(.8, 3, .1)
[1] 42.7903
quantile(x, .8)
    80% 
43.4998 
qgamma(.8, 2.77, .0906)
[1] 44.02735

Suppose you want the probability $P(X \le 60),$ for $X \sim \mathsf{Gamma}(3, 0.1).$ The exact value is $P(X \le 60) = 0.9380.$ If you want to estimate this probability directly from the sample, you can note that 93 of the 100 observed values are at or below 60. If you want to estimate this probability from the CDF with estimated values, you find $P(X \le 60) \approx 0.927.$

pgamma(60, 3, .1)
[1] 0.9380312
mean(x <= 60)
[1] 0.93
pgamma(60, 2.77, .0906)
[1] 0.9269133

Moreover, you can plot the CDF of $\mathsf{Gamma}(3, 0.1),$ as shown in both plots below. Superimposed (in red) on the plot at left is the empirical CDF (ECDF) of our sample, which 'jumps up' by $1/100$ at each of the 100 sampled values. (Multiple-sized jumps in case rounding had caused ties, but there are no ties in my x.) In the plot at right we superimpose the CDF with MMEs from the sample instead of the actual population parameters. (R code for the plot is provided in the notes.)

enter image description here

Maximum likelihood estimates. Generally speaking, maximum likelihood estimates (MLEs) are better than method of moments estimates, but they require numerical methods beyond simple arithmetic. The estimate of the shape parameter $\alpha$ is the difficult part. See Wikipedia, where the notation is a little different from mine.

In Bain & Englehardt, Intro. to probability and mathematical statistics, 2e (1992), p300, the following are given as approximations to the MLE $\hat \alpha.$

For the ratio $M = \ln(\bar X/\tilde X),$ where $\tilde X = [\prod_{i-1}^n X_i]^{1/n}$ is the geometric mean:

For $0 \le M \le 0.5772,$ use $\hat \alpha = (0.5000876 + 0.1648852M - 0.0544247M^2)/M;$

for $0.5772 < M \le 17.$ use $\hat \alpha = \frac{8.898919 + 9.059950M + 0.9775373M^2} {M(17.79728 + 11.968477M + M^2};$

and for $M> 17,$ use $\hat \alpha = 1/M.$ Then $\hat \lambda = \bar X/ \hat \alpha.$

For our data above, the resulting approximate MLEs are $\hat\alpha = 2.96$ and $\hat \lambda = 0.097.$

a; g; M
[1] 30.53406
[1] 25.54886
[1] 0.1782502
alp.mle = (0.5000876 + 0.1648852*M - 0.0544247*M^2)/M; alp.mle
[1] 2.960722
lam.mle = alp.mle/a;  lam.mle
[1] 0.09696457

Notes: (1) R code for the plot is provided below:

par(mfrow = c(1,2))
 lbls = "CDF of GAMMA(3, .1) with ECDF of Sample of 100"
 curve(pgamma(x, 3, .1), 0, 100, lwd=2, ylab="CDF", main=lbls)
  abline(v = 0, col="green2");  abline(h=0:1, col="green2")
  lines(ecdf(x), pch=".", col="red")
 lble = "CDF of GAMMA(3, .1) with CDF Using MMEs"
 curve(pgamma(x, 3, .1), 0, 100, lwd=2, ylab="CDF", main=lble)
  curve(pgamma(x, 2.77, .091), add=T, lwd=1, lty="dashed", col="red")
  abline(v = 0, col="green2"); abline(h=0:1, col="green2")
par(mfrow = c(1,1))

(2) In the column to the right, under 'Related', you will find links to somewhat similar Q & A's. Please look at them to see if they contain information that is useful to you.

(3) The approximate MLE method quoted from Bain & Englehardt references Greenwood & Durand (1960) in Technomerics. Modern MLE software may be better or more convenient.

(4) Minitab 'Quality tools > Individual Dist'n ID' identifies this sample as fitting a gamma distribution (among others). Assuming it is gamma, here are Minitab's MLEs, which are consistent with those obtained above with the Greenwood-Durand approximation. I suspect many other statistical software packages have similar capabilities. (Perhaps we'll see Comments to that effect.)

ML Estimates of Distribution Parameters

Distribution  Location    Shape     Scale
Gamma                   2.96121  10.31134
$\endgroup$
2
  • 1
    $\begingroup$ +1 I could not imagine a better answer. $\endgroup$
    – Noah
    Commented Jul 12, 2019 at 0:28
  • 1
    $\begingroup$ This was incredibly helpful, thank you. $\endgroup$ Commented Jul 12, 2019 at 17:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.