5.1 Gamma with unknown shape

This example deals with the Data Set B of the text which is supposed to represent twenty workers compensation medical benefit claims. The data is presented in the tables below - I demonstrate the use of two packages to nicely tabulate your data.

The gamma family of distributions is chosen as a model for this data, and the goal is to estimate the unknown shape and scale parameters. Just to fix the parametrization, the density is given by \[ f_{\alpha,\theta}(x)=\frac{1}{\Gamma(\alpha)\theta^\alpha} e^{-(x/\theta)} x^{\alpha-1 },\quad x\geq 0; \] the parameter space is given by \((0,\infty)^2\). We note that, the log-likelihood is given by, \[ l(\alpha,\theta):=-n\log\Gamma(\alpha) - n\alpha\log\theta - \frac{1}{\theta}\sum_{i=1}^n x_i +(\alpha-1)\sum_{i=1}^n \log x_i. \] Note that \[ \frac{\partial }{\partial \theta}l(\alpha,\theta) =-\frac{n\alpha}{\theta} +\frac{1}{\theta^2}\sum_{i=1}^n x_i, \] and equating it to zero yields, \[ \hat{\theta}_{MLE}= \frac{\bar{x}}{\hat{\alpha}_{MLE}}. \] This is so as \[ \frac{\partial^2 }{\partial \theta^2}l(\alpha,\theta)\Bigg\vert_{(\alpha=\hat{\alpha}_{MLE},\theta=\hat{\theta}_{MLE})} =\frac{n\hat{\alpha}_{MLE}}{\hat{\theta}_{MLE}^2} -\frac{2}{\hat{\theta}_{MLE}^3}\sum_{i=1}^n x_i=\frac{n\hat{\alpha}_{MLE}}{\hat{\theta}_{MLE}^3} \left(\hat{\theta}_{MLE}-\frac{2\bar{x}}{\hat{\alpha}_{MLE}}\right) = -\frac{n\hat{\alpha}_{MLE}}{\hat{\theta}_{MLE}^2} <0. \]

While we have an expression for \(\hat{\theta}_{MLE}\), it involves an unknown - \(\hat{\alpha}_{MLE}\). The nice thing about it is that it still helps reduce the optimization problem from one on \(\mathbb{R}_+^2\) to \(\mathbb{R}_+\). This is so because we reduce the optimization problem to that of maximizing \[ \begin{aligned} l(\alpha,\bar{x}/\alpha)&:=-n\log\Gamma(\alpha) - n\alpha\log\frac{\bar{x}}{\alpha} - \frac{\alpha}{\bar{x}}\sum_{i=1}^n x_i +(\alpha-1)\sum_{i=1}^n \log x_i\\ &=-n\log\Gamma(\alpha) - n\alpha\log{\bar{x}}+ n\alpha\log{\alpha} - n\alpha + (\alpha-1)\sum_{i=1}^n \log x_i. \end{aligned} \]

An alternate method is to find the root of \[ \begin{aligned} \frac{\partial}{\partial\alpha}l(\alpha,\bar{x}/\alpha)&= -n\frac{\Gamma'(\alpha)}{\Gamma(\alpha)} - n\log{\bar{x}}+ n\log{\alpha} + \sum_{i=1}^n \log x_i\\ &=n\left(-\frac{\Gamma'(\alpha)}{\Gamma(\alpha)} - \log{\bar{x}}+ \log{\alpha} + \frac{1}{n}\sum_{i=1}^n \log x_i\right) \end{aligned} \] Some comments on the expression within the brackets:

Note that the expression is scale invariant i.e. we can wolog assume that \(\theta=1\) (Why?). Now observe that for a \(Gamma(\alpha,1)\) random variable \(X\),

\[ \mathbb{E}\left({log(X)}\right)=\frac{\Gamma'(\alpha)}{\Gamma(\alpha)},\quad \hbox{and} \quad \mathbb{E}\left({X}\right)=\alpha. \]

The above comment implies that we seek a root of (with both sides being non-negative by Jensen’s inequality)

\[ \log{\mathbb{E}_\alpha\left({X}\right)} -\mathbb{E}_\alpha\left({log(X)}\right) = \left(\log{\bar{x}} - \frac{1}{n}\sum_{i=1}^n \log x_i \right) \]

Finally we note that by Jensen’s inequality (as \(log\) is a concave function) both sides of the above equation are non-negative. We below use these comments to implement another algorithm to find \(\hat{\alpha}_{MLE}\) - this uses the digamma function.

5.2 Grouped Data

In this example we will look at the analysis of Data Set C in the text which represents 227 payments from a general liability insurance policy block. The data set is given below (see Example 13.5 of the text).

This data set is an example of a grouped data set. We fit an exponential model to this. First, we realize that since the exponential is a scale family, we can change the unit to \(\$10,000\), estimate the scale parameter \(\theta\) and scale it back by a \(10,000\). To find the MLE we need the likelihood - we write a function and plot it below.