Confidence regions of optimized parameters in maximum log likelihood fits

Question

I am using a numerical optimization algorithm to maximize a log-likelihood function, $\mathcal{L}$. The log-likelihood function has a fixed number of parameters, $\{\theta_i\}$. These parameters are optimized by the optimization algorithm to maximize the log-likelihood function, from which the best estimates of the values of the parameters $\theta_i$ are obtained.

I know that it is possible to derive, or calculate, some critical values of $\Delta \mathcal{L}$, which can be used to obtain confidence regions for the parameters $\theta_i$.

However, I do not know how to find this information, or where to find this information.

I looked through Statistical Inference by Casella and Berger, expecting to find a table of values for $\Delta \mathcal{L}$, or some method to calculate such a table of values. I expected to find this in the Chapter Interval Estimation, but did not. If the information is there, then I didn't understand what I was reading, and therefore it went over my head.

Can anyone point me in the right direction?

Here's some further information:

• Since $\mathcal{L}$ is a function of the parameters $\theta_i$, by varying the values of $\theta_i$ we also change the value of $\mathcal{L}$.
• The best estimates of the parameters, $\hat{\theta}_i$, are obtained by maximizing $\mathcal{L}$.
• If we want to estimate the confidence regions associated with $\theta_i$, we can change one of the values of $\theta_i$ until the value of the log-likelihood function $\mathcal{L}$ changes by some critical value $\Delta \mathcal{L}$.
• The value of $\Delta \mathcal{L}$ depends both on the number of parameters $i$ and the target confidence level.
• For example, it is larger for larger confidence bands.
• If we want a 99 % CL, then the value of $\Delta \mathcal{L}$ required will be larger than if we were to estimate a 95 % CL region.

Answer 1

The answer is pretty much implicit in this classic Cross Validated answer, which describes the 3 types of hypothesis tests used for models based on maximum likelihood. The following might help put that answer into the context of your question.

The first basic idea is that estimating a confidence interval (CI) is inverting a null hypothesis test at the corresponding type-I error. For example, a 95% confidence interval is constructed by finding parameter-value limits consistent with type-I error $\alpha < 0.05$.

The second basic idea is that inference with maximum likelihood estimates (MLE) typically relies on their asymptotic properties (given some regularity conditions) as the sample size grows large. See Chapter 10, "Asymptotic Evaluations," of Casella and Berger (second edition). Citations below are to that Chapter. Their presentation is on single-parameter estimation, but the principles extend to multiple parameters as in your question. The estimates are approximations in several senses, but they are generally what there is to work with.

Theorem 10.1.12 shows that MLE of parameters have an asymptotically normal distribution centered on the MLE. This Cross Validated answer shows how the estimate of the (in general, multivariate) normal distribution of parameters from MLE comes directly from the numerical procedures used to maximize the log-likelihood. The variance-covariance matrix of parameter estimates is the inverse of the matrix of second partial derivatives of the likelihood* evaluated at the MLE.

The variances of the individual parameter estimates are the diagonal elements of the variance-covariance matrix, so you can get corresponding CI from the properties of a normal distribution. Any interval that includes 95% of the distribution can be chosen in principle, but standard practice is to have them symmetric about the point estimate in some appropriate scale/transformation of the parameter. The multivariate normal distribution can be used to construct more general Wald tests on combinations of the parameter estimates. That allows for "chunk tests" on multiple parameter estimates at once, or on the linear combinations of parameter estimates used for making predictions from multiple-regression models.

Theorem 10.3.1 is the basis for a second way to construct confidence intervals, the likelihood-ratio test, considered more reliable with smaller sample sizes and unaffected by monotone transformations of parameters. Asymptotically, in distribution for a single parameter as the sample size increases:

$$-2 \log \lambda (X) \rightarrow \chi^2_1$$

where $\lambda(X)$ is the likelihood ratio at the MLE versus the null, given data $X$, and $\chi^2_1$ is a $\chi^2$ variable with one degree of freedom. With more parameters, the degrees of freedom for $\chi^2$ is the number of parameters estimated. The test of the null hypothesis is based on the corresponding $\chi^2$ distribution.

An important caution: this result, Wilks's Theorem, only holds if the MLE is in the interior of the parameter space. If the MLE is on a boundary (for example, a variance estimate of 0), then the result doesn't hold. Then this $\chi^2$ test (and its inversion to estimate CI) isn't valid.

Inverting a likelihood-ratio hypothesis test to get a confidence interval for a parameter requires re-fitting the rest of the model over a range of pre-specified values of the parameter in question, to find a range consistent with the desired type-I error based on the $\chi^2$ distribution. That's called profiling the likelihood, as illustrated in an answer linked above. That will also work in situations where the Wald test breaks down, as illustrated in this answer for the partial likelihood of a Cox survival model.

---

*or the negative of that matrix, depending on whether the algorithm was maximizing the likelihood or minimizing its negative

Answer 2

You need to use the Hessian matrix of the log-likelihood function

Note on notation: I am going to use slightly different notation to you, to conform with standard notation in the field. I will denote the log-likelihood function as $\ell$ (instead of $\mathcal{L}$) and the differentiation operator as $\nabla$ (instead of $\Delta$).

When you find the MLE from the log-likelihood function, typically its first derivative is going to be zero (or at least, approximately zero if you did things numerically). This occurs because the MLE typically falls at a critical point of the log-likelihood function. Consequently, the value $\nabla \ell (\hat{\boldsymbol{\theta}})$ is not particularly useful, other than to check if your numerical solution ended at or near a critical point.

To form a confidence interval for the parameter vector $\boldsymbol{\theta}$ (or the individual parameters that are elements of this vector) you need to use the second derivative of the log-likelihood function, which is the Hessian matrix (the matrix of second-order partial derivatives of the log-likelihood function):

$$\mathbf{H}(\boldsymbol{\theta}) \equiv \nabla^2 \ell(\boldsymbol{\theta}) = \bigg[ \frac{\partial^2 \ell}{\partial \theta_i \partial \theta_j} \bigg]_{i,j}.$$

Under certain regularity conditions (which basically require the MLE to fall at a critical point and require the log-likelihood to be sufficiently smooth) there is a classic theorem on the behaviour of the MLE which says that for large $n$ you have:$^\dagger$

$$\sqrt{n} (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}}) \overset{\text{Approx}}{\sim} \text{N}(\mathbf{0}, \mathbf{H}^{-1}(\boldsymbol{\theta})).$$

This approximate distributional result can be used to form a confidence region for the entire parameter vector or individual (marginal) confidence intervals for the individual parameter values.

---

CI for the parameter vector: Taking $\chi_{n, \alpha}^2$ to be the critical point of the chi-squared distribution with $n$ degrees-of-freedom and upper tail $\alpha$, you can form the following confidence region for the parameter vector:

$$\text{CI}(1-\alpha) \equiv \bigg\{ \boldsymbol{\theta} \bigg| (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}})^\text{T} \mathbf{H}^{-1}(\hat{\boldsymbol{\theta}}) (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}}) \leqslant \chi_{n, \alpha}^2 \bigg\}.$$

For large $n$ you have $\mathbb{P}(\boldsymbol{\theta} \in \text{CI}(1-\alpha) | \boldsymbol{\theta}) \approx 1-\alpha$, which means that the confidence interval has roughly the correct coverage probability for large $n$.

---

CI for individual parameters: Taking $z_{\alpha/2}$ to be the critical point of the standard normal distribution with upper tail $\alpha/2$, you can form the following confidence region for an individual parameter:

$$\text{CI}_k(1-\alpha) \equiv \bigg[ \hat{\theta}_k - z_{\alpha/2} \hat{\text{se}}_{k}, \hat{\theta}_k + z_{\alpha/2} \hat{\text{se}}_{k} \bigg]. \quad \quad \quad \hat{\text{se}}_{k} \equiv \sqrt{\mathbf{H}^{-1}(\hat{\boldsymbol{\theta}})_{k,k}}.$$

For large $n$ you have $\mathbb{P}(\theta_k \in \text{CI}_k(1-\alpha) | \boldsymbol{\theta}) \approx 1-\alpha$, which means that the confidence interval has roughly the correct coverage probability for large $n$.

---

$^\dagger$ For large $n$ you have $\hat{\boldsymbol{\theta}} \approx \boldsymbol{\theta}$ (due to consistency of the MLE under the regularity conditions) so the result is sometimes stated with the estimated parameter as the argument to the Hessian matrix. Both approximating results are valid for large $n$.