2.4 Method of Maximum Likelihood
The idea of maximum likelihood is to estimate \(\theta \) by the value in the parameter space that makes the observed data most probable.
Definition 2.4.1. Let \(X_1, \dots , X_n\) be a random sample from \(f(x;\theta )\), \(\theta \in \Omega \), with observed values \(\mathbf x = (x_1, \dots , x_n)\). The likelihood function of \(\theta \) is \[L(\theta ) = L(\theta ;\mathbf x) = \prod ^n_{i = 1}f(x_i;\theta ) ,\qquad \theta \in \Omega ,\] and the log likelihood is \(l(\theta ) = \ln L(\theta )\).
Definition 2.4.2. The value \(\hat \theta \) of \(\theta \) that maximises \(L(\theta )\), equivalently \(l(\theta )\), is the maximum likelihood (M.L.) estimate. As a function of the data, \(\hat {\theta } = \hat {\theta }(x_1, \dots , x_n)\); the corresponding random variable \(\hat {\theta }(X_1, \dots , X_n)\) is the M.L. estimator.
Definition 2.4.3. For a scalar \(\theta \), the score function and the information function are \[S(\theta ) = \frac {\partial l(\theta )}{\partial \theta } ,\qquad I(\theta ) = -\frac {\partial ^2 l(\theta )}{\partial \theta ^2} .\] For a vector \(\theta = (\theta _1, \dots , \theta _k)^{\mathsf T}\), the score vector has components \(\frac {\partial l}{\partial \theta _j}\) and the information matrix is the \(k\times k\) symmetric matrix with entries \[I(\theta )_{ij} = -\frac {\partial ^2 l(\theta )}{\partial \theta _i\partial \theta _j} .\]
These functions depend on the data as well as on \(\theta \); we write \(L(\theta ;\mathbf x)\), \(S(\theta ;\mathbf x)\) and \(I(\theta ;\mathbf x)\) when the dependence matters, and \(L(\theta ;\mathbf X)\) and so on when studying their sampling properties as functions of the random sample.
Finding M.L. Estimates
If \(l(\theta )\) is twice differentiable on an open interval, \(\hat \theta \) is usually found by solving the likelihood equation \(S(\theta ) = 0\), and for a scalar \(\theta \) the solution is confirmed to be a maximum by checking \(I(\hat {\theta }) > 0\). Two warnings: if the support of the distribution depends on \(\theta \), as for the uniform distribution on \((0,\theta )\), the maximum is usually at an end point and is not found from \(S(\theta ) = 0\); and often \(S(\theta ) = 0\) has no explicit solution and must be solved numerically.
Newton’s Method
Expanding the score about a current estimate \(\theta ^{(i)}\), \(S(\theta )\approx S\left (\theta ^{(i)}\right ) - I\left (\theta ^{(i)}\right )\left (\theta - \theta ^{(i)}\right )\), and setting this to zero gives the Newton–Raphson iteration \[\theta ^{(i + 1)} = \theta ^{(i)} + \left [I\left (\theta ^{(i)}\right )\right ]^{-1} S\left (\theta ^{(i)}\right ) ,\qquad i = 0, 1, 2, \dots ,\] started from an initial estimate \(\theta ^{(0)}\), for instance the method of moments estimate, and continued until successive values agree to the required accuracy.
Example 2.4.4. Let \(X_1, \dots , X_n\) be a random sample from \(f(x;\theta ) = \theta x^{\theta - 1}\), \(0< x < 1\), \(\theta > 0\). Find the score function, the information function and the M.L. estimator of \(\theta \).
Solution. The likelihood is \(L(\theta ) = \prod \theta x_i^{\theta - 1} = \theta ^n\left (\prod x_i\right )^{\theta - 1}\), so \[l(\theta ) = n\ln \theta + (\theta - 1) \sum ^n_{i = 1} \ln x_i ,\qquad S(\theta ) = \frac {n}{\theta } + \sum ^n_{i = 1} \ln x_i ,\qquad I(\theta ) = \frac {n}{\theta ^2} .\] Setting \(S(\theta ) = 0\) gives \[\hat {\theta } = -\frac {n}{\sum ^n_{i = 1}\ln X_i} ,\] which is positive because each \(\ln X_i<0\), and \(I(\hat \theta )>0\), so it is the maximum.
Example 2.4.5. Let \(X_1,\dots ,X_n\) be a random sample from the uniform distribution on \((0,\theta )\). Find the M.L. estimator of \(\theta \).
Solution. Each factor of the likelihood is \(\frac 1\theta \) if \(0<x_i<\theta \) and \(0\) otherwise, so \(L(\theta ) = \theta ^{-n}\) when \(\theta >\max x_i = y_n\), and \(L(\theta ) = 0\) when \(\theta <y_n\). On \(\theta >y_n\) the likelihood decreases, so it is largest at the smallest admissible value: \(\hat \theta = Y_n\), the largest observation. The likelihood equation \(S(\theta ) = -\frac {n}{\theta } = 0\) has no solution; the maximum is at the edge of the region where \(L>0\).
Example 2.4.6. Let \(X_1,\dots ,X_n\) be a random sample from the uniform distribution on \(\left (\theta - \frac {1}{2},\theta + \frac {1}{2}\right )\). Find an M.L. estimator of \(\theta \).
Solution. With \(I_A\) the indicator of the set \(A\), the likelihood is \[L(\theta ) = \prod ^n_{i=1}I_{\left (\theta - \frac 12,\theta + \frac 12\right )}(x_i) = I_{\left (\theta -\frac {1}{2},\infty \right )}(y_1)\cdot I_{\left (-\infty ,\theta + \frac {1}{2}\right )}(y_n),\] where \(y_1\) and \(y_n\) are the smallest and largest observations. So \(L(\theta ) = 1\) exactly when \(y_1>\theta - \frac 12\) and \(y_n<\theta + \frac 12\), that is \[y_n - \frac {1}{2} < \theta < y_1 + \frac {1}{2},\] and \(L(\theta ) = 0\) otherwise. Every \(\theta \) in this interval maximises the likelihood: the M.L. estimator is not unique. The usual choice is the midpoint, \(\hat {\theta } = \frac {Y_1 + Y_n}{2}\), the midrange.
Example 2.4.7. Let \(X_1, \dots , X_n\) be a random sample from the Weibull p.d.f. \(f(x;\beta ) = \beta x^{\beta - 1} e^{-x^{\beta }}\), \(x > 0\), \(\beta > 0\). Set up the Newton–Raphson iteration for the M.L. estimator of \(\beta \).
Solution. The log likelihood is \[l(\beta ) = n\ln \beta + (\beta - 1)\sum ^n_{i = 1}\ln x_i - \sum ^n_{i = 1}x_i^{\beta } ,\] so \[S(\beta ) = \frac {n}{\beta }+ \sum ^n_{i = 1}\ln x_i - \sum ^n_{i = 1} x_i^{\beta } \ln x_i ,\qquad I(\beta ) = \frac {n}{\beta ^2} + \sum ^n_{i = 1} x_i^{\beta }(\ln x_i)^2 .\] The equation \(S(\beta ) = 0\) has no explicit solution, but \(I(\beta )>0\) for every \(\beta \), so \(l\) is concave and has a single maximum, which the iteration \(\beta ^{(i + 1)} = \beta ^{(i)} + S\left (\beta ^{(i)}\right )/I\left (\beta ^{(i)}\right )\) finds from any reasonable starting value.
Theorem 2.4.8 (Invariance property). If \(\hat {\theta }\) is the M.L. estimator of \(\theta \) and \(\tau = h(\theta )\) is a one-to-one function, then \(\hat {\tau } = h(\hat {\theta })\) is the M.L. estimator of \(\tau \).
Proof. In terms of \(\tau \) the likelihood is \(L^*(\tau ) = L\left (h^{-1}(\tau )\right )\), which is maximised when \(h^{-1}(\tau ) = \hat \theta \), that is \(\tau = h(\hat \theta )\). □
Example 2.4.9. For a random sample from the Poisson distribution with mean \(\theta \), find the M.L. estimator of \(e^{-1/\theta }\).
Solution. The log likelihood is \[l(\theta ) = \left (\sum ^n_{i = 1} x_i\right )\ln \theta - n\theta - \sum ^n_{i = 1}\ln (x_i!) ,\qquad S(\theta ) = \frac {1}{\theta }\sum ^n_{i = 1} x_i - n ,\] so \(\hat {\theta } = \bar {X}\). By Theorem 2.4.8, the M.L. estimator of \(e^{-1/\theta }\) is \(e^{-1/\bar X}\). Similarly, the M.L. estimator of \(P(X = 0) = e^{-\theta }\) is \(e^{-\bar X}\).
Definition 2.4.10. For a scalar \(\theta \), the Fisher information is the expected information, \[J(\theta ) = E\left [I(\theta ;\mathbf X)\right ] = E\left [-\frac {\partial ^2 l(\theta ;\mathbf X)}{\partial \theta ^2}\right ].\] For a vector \(\theta \), the Fisher information matrix \(J(\theta )\) has entries \(E\left [-\frac {\partial ^2 l(\theta ;\mathbf X)}{\partial \theta _i\partial \theta _j}\right ]\).
For a random sample, the log likelihood is a sum of \(n\) independent terms \(\ln f(X_i;\theta )\), so \[J(\theta ) = n\,E\left [-\frac {\partial ^2 \ln f(X;\theta )}{\partial \theta ^2}\right ],\] \(n\) times the information in a single observation.
Example 2.4.11. Let \(X_1, \dots , X_n\) be a random sample from \(f(x;\theta ) = \frac {1}{\theta }e^{-x/\theta }\), \(x > 0\). Find the score function, the information function, the Fisher information and the M.L. estimator of \(\theta \).
Solution. The likelihood is \(L(\theta ) = \theta ^{-n}e^{-\sum x_i/\theta }\), so \[l(\theta ) = -n\ln \theta - \frac {1}{\theta }\sum ^n_{i = 1} x_i ,\qquad S(\theta ) = -\frac {n}{\theta } + \frac {1}{\theta ^2}\sum ^n_{i = 1} x_i ,\qquad I(\theta ) = -\frac {n}{\theta ^2} + \frac {2}{\theta ^3}\sum ^n_{i = 1} x_i .\] Since \(E(X_i) = \theta \), \[J(\theta ) = -\frac {n}{\theta ^2} + \frac {2}{\theta ^3}\cdot n\theta = \frac {n}{\theta ^2} .\] Setting \(S(\theta ) = 0\) gives \(\hat \theta = \frac 1n\sum X_i = \bar X\), and \(I(\hat \theta ) = \frac {n}{\bar x^2}>0\).
Questions on this section
Stuck on something here? Ask below and it stays attached to this topic.