MLE Gamma

Author

Prof. Eric A. Suess

MME fitted gamma

\(X_1, X_2, ...,X_n\) iid \(Gamma(\alpha,\lambda)\), find the MMEs.

\(f(x|\alpha,\lambda)= \frac{\lambda^\alpha}{\Gamma(\alpha)} x^{\alpha-1} e^{-\lambda x}\)

n <- 100

alpha <- .5; a.min <- 0.1; a.max <- 2;
lambda <- 2; l.max <- 10


# generate random data

x <- rgamma(n, shape = alpha, rate = lambda)

x.orig <- x # save a copy of the original data

make a histogram

use the density estimator since this is a continuous random variable

hist(x.orig, probability=T)
lines(density(x.orig))

MMEs

sigma2.hat <- (1/n)*sum( (x - mean(x))^2 )

alpha.hat.mme <- mean(x)^2/sigma2.hat
alpha.hat.mme
[1] 0.3472268
lambda.hat.mme <- mean(x)/sigma2.hat
lambda.hat.mme
[1] 1.63452

The AIC and BIC are defined as follows, with \(n\) the sample size and \(k\) the number of parameters in the model.

\[AIC = -2\log(L(\hat{\theta})) + 2k\]

\[BIC = -2\log(L(\hat{\theta})) + k\log(n)\]

library(fitdistrplus)
Loading required package: MASS
Loading required package: survival
fit.gam <- fitdist(x.orig, "gamma", method = "mme")
summary(fit.gam)
Fitting of the distribution ' gamma ' by matching moments 
Parameters : 
       estimate Std. Error
shape 0.3472268  0.9672571
rate  1.6345201  5.3316147
Loglikelihood:  76.36084   AIC:  -148.7217   BIC:  -143.5113 
Correlation matrix:
          shape      rate
shape 1.0000000 0.8540047
rate  0.8540047 1.0000000

plot the fitted model to the simulated data

denscomp(list(fit.gam), legendtext = c("Gamma"))

cdfcomp(list(fit.gam), legendtext = c("Gamma"))

qqcomp(list(fit.gam), legendtext = c("Gamma"))

ppcomp(list(fit.gam), legendtext = c("Gamma"))

MLE fitted gamma

fit.gam <- fitdist(x.orig, "gamma", method = "mle")
summary(fit.gam)
Fitting of the distribution ' gamma ' by maximum likelihood 
Parameters : 
       estimate Std. Error
shape 0.4788249 0.05567706
rate  2.2542167 0.41812618
Loglikelihood:  79.86981   AIC:  -155.7396   BIC:  -150.5293 
Correlation matrix:
          shape      rate
shape 1.0000000 0.6268848
rate  0.6268848 1.0000000

plot the fitted model to the simulated data

denscomp(list(fit.gam), legendtext = c("Gamma"))

cdfcomp(list(fit.gam), legendtext = c("Gamma"))

qqcomp(list(fit.gam), legendtext = c("Gamma"))

ppcomp(list(fit.gam), legendtext = c("Gamma"))