7
$\begingroup$

When we do maximum likelihood estimation of a (let's say scalar) parameter $\theta$, we get a point estimate wherever the likelihood function is maximized. However, if we want an interval estimate, the points giving nearly maximum likelihoods seem like good candidates for points to include in the interval estimate.

  1. (Main question) Is there a way to relate something like confidence level (I think I mean inverting the interval to give a hypothesis test) to give an interval estimate based on some criteria like, "I want all points that give likelihood within $\text{x}\%$ of the maximum likelihood"?

  2. If we form a confidence interval based on inverting an $\alpha$-level likelihood ratio test, is that basically what we're doing?

It would be nice to have decent conditions that will assure the interval estimate is indeed a connected interval, but that is of secondary importance.

EDIT

Here is a sketch of the situation. I mean the interval between the red vertical lines., where the red horizontal line is $\text{x}\%$ of the maximum value of the purple likelihood function.

enter image description here

$\endgroup$
2
  • 2
    $\begingroup$ I think you're asking for a procedure for calculating $1-\alpha$ confidence intervals that guarantees the parameter values within the interval have higher likelihood, given the data, than those without; & whether intervals obtained by inverting the LRT fit the bill (they do for some parametric families, but not in general). $\endgroup$ Mar 15, 2022 at 7:33
  • 1
    $\begingroup$ I'm pretty sure what you have in mind are what's called "profile-likelihood-based confidence intervals", see for instance Venzon & Moolgavkar. $\endgroup$
    – Durden
    Dec 30, 2022 at 23:06

4 Answers 4

6
$\begingroup$

The likelihood function $L(\theta)$ provides us with both, a point estimator and a confidence set estimator for $\theta$. In particular, to make things simple, let $X_1,\ldots, X_n$ be an i.i.d. random sample from a model having pdf $f_\theta$, with $\theta\in \mathbb{R}$.

The Maximum Likelihood Estimator is defined as

$$ \hat \theta = \underset{\theta\in\Theta}{\arg \max}\, L(\theta), $$

or, under appropriate smoothing conditions, as the solution to the likelihood equation $$ \frac{d\ell(\theta)}{d\theta} = 0, $$

where $\ell(\theta) = \log L(\theta)$. Under appropriate conditions, and if $\theta$ is the true parameter, we have

\begin{align}\label{1} \Lambda_n(\theta) = -2(\ell(\theta)-\ell(\hat\theta)) &\overset{d}\to \chi_1^2,\quad\quad\text{(*)} \end{align}

that is to say $$ \Lambda_n(\theta)\,\,\dot\sim\,\, \chi_1^2 $$

where the symbol "$\,\dot\sim\,$" means "asymptotically as $n\to\infty$ distributed as".

The point here is that $\Lambda_n(\theta)$ provides an asymptotic pivotal quantity, which we can use to do hypothesis testing and build confidence sets. To go straight on to the matter of your post, I'll focus here on the latter.

By definition, a confidence set of level $1-\alpha$ is a random set which traps the true value with a probability no lower than $1-\alpha$. Thus if $C_{1-\alpha}$ is some set s.t.

$$ P_{\theta}(\Lambda_n(\theta)\in C_{1-\alpha})\geq 1-\alpha,\quad\forall\theta, $$

then the set

$$ \{\theta:\Lambda_n(\theta)\in C_{1-\alpha}\} $$

forms a (random) confidence set with probability coverage $1-\alpha$. The set $C_{1-\alpha}$ can be determined in different ways, but the usual approach is to use the threshold $\chi_{1,1-\alpha}^2$, and look for values of $\Lambda_n(\theta)$ that are below $\chi_{1,1-\alpha}^2$. Indeed, for every fixed $\theta$, $$ P_\theta(\Lambda_{n}(\theta) \leq \chi_{1,1-\alpha}^2)\,\, \dot=\,\, 1-\alpha $$

where $\chi_{1,1-\alpha}^2$ is the upper $\alpha$th level quantile of the $\chi_1^2$ distributions. The usual likelihood-based confidence set is thus

$$ \{\theta:\Lambda_n(\theta)\leq \chi_{1,1-\alpha}^2\}. $$

The R-example below illustrates this in the case of $f_\theta$ being the Poisson distribution. In the figure, the horizontal dashed line represents the 0.95th quantile of $\chi_1^2$ distribution, and the dashed vertical lines mark the limits of the confidence set, which in this case happens to be an interval.

library(latex2exp)

# fix some observed data
y <- c(5, 4, 1, 0, 0, 1, 1, 2, 1, 1)

# the log-likelihood function
llik <- function(lambda, y){
  oo = dpois(y, lambda = lambda, log = TRUE)

  return(sum(oo))
}

ybar <- mean(y)

x <- seq(.1, 4, len=100)
ll <- sapply(x, function(x) llik(x,y=y))

llikv <- Vectorize(function(x) llik(x,y=y), "x")

lo <- uniroot(function(x) 2*(llik(ybar, y)-llikv(x))-qchisq(0.95, df=1),
        lower = 0.1, upper = ybar)
up <- uniroot(function(x) 2*(llik(ybar, y)-llikv(x))-qchisq(0.95, df=1),
              lower = ybar, upper=4)


plot(y = 2*(llik(ybar, y)-ll), x=x,
     type="l", ylab=NA,
     xlab = expression(theta),
     xlim=c(0.1,4))
abline(v = ybar, lwd=2, lty="dotted")
mtext(TeX("$\\widehat{\\theta}$"),side=1, at=1.6,padj = 1)
abline(h = qchisq(0.95, df=1), lwd=2, lty=2, col="gray")
segments(x0=lo$root, x1=lo$root,y0=-1.5,y1=8+qchisq(0.95, df=1), lwd=2, lty=2)
segments(x0=up$root, x1=up$root,y0=-1.5,y1=8+qchisq(0.95, df=1), lwd=2, lty=2)
mtext(TeX("$\\chi_{1,.05}^2$"),
      side=1,
      at=3,padj = -2.5,las=1,adj=0)
mtext("95% lik. conf. set",
      side=1,
      at=1,padj = -4.5,las=1,adj=0)
text(x = 0.5, y=50,TeX("$2(l(\\widehat{\\theta})-l(\\theta))$"))

enter image description here

$\endgroup$
2
$\begingroup$

Under some conditions you can show that the MLE of a parameter has the following distribution: $\sqrt n (\hat\theta-\theta)\sim\mathcal N(0,\hat\sigma^2)$, where $\hat\sigma^2$ is the corresponding estimate of the variance. The variance estimate $\hat\sigma^2/n$ can be used to estimate the confidence interval of the parameter estimate $\hat\theta$. The variance estimate is related to the Fisher information of the distribution (likelihood) function.

$\endgroup$
8
  • $\begingroup$ How does this help get the points within $\text{x}\%$ of the maximum likelihood? $\endgroup$
    – Dave
    Jun 2, 2020 at 19:25
  • $\begingroup$ If you agree that $(\hat\theta-\theta)\sim\mathcal N(0,\hat\sigma^2/n)$ then you can get all the percentile in usual way from the normal distribution $\endgroup$
    – Aksakal
    Jun 2, 2020 at 19:43
  • $\begingroup$ I don't think I'm talking about the probability of the MLE. I think I mean the y-axis values. If at the MLE the likelihood function has a value of 2, I am looking to get all inputs to the likelihood function that return values of at least 1.8, for instance. $\endgroup$
    – Dave
    Jun 2, 2020 at 20:04
  • $\begingroup$ Are you looking for inverse problem? The domain of x within likelihood interval? $\endgroup$
    – Aksakal
    Jun 2, 2020 at 20:07
  • $\begingroup$ Maybe. I've edited the original question with a sketch. $\endgroup$
    – Dave
    Jun 2, 2020 at 20:13
2
$\begingroup$

You have likelihood intervals/regions and confidence intervals/regions. They are different.

Confidence intervals stem from the fiducial distribution, which is a different thing from the likelihood function.


If we form a confidence interval based on inverting an $\alpha$-level likelihood ratio test, is that basically what we're doing?

The likelihood ratio test is using a likelihood ratio instead of a likelihood function. The goal is to make a hypothesis test that is most powerful for some alternative hypothesis.

  • The likelihood ratio is a function of the data (given a null hypothesis $\theta_0$ and alternative hypothesis $\theta_a$) $$\Lambda(x|\theta_0,\theta_a) = \frac{\mathcal{L}(\theta_0|x)}{\mathcal{L}(\theta_a|x)}$$
  • The likelihood or relative likelihood is a function of the distribution parameter (given data $x$). $$\mathcal{L}_{\text{relative}}(\theta|x) = \frac{\mathcal{L}(\theta|x)}{\mathcal{L}(\hat\theta_{ML}|x)}$$ it scales the likelihood function $\mathcal{L}(\theta|x)$ by the maximum $\mathcal{L}(\hat\theta_{ML}|x)$

The likelihood ratio $\Lambda(x)$ is a way to define a hypothesis test, and that hypothesis test can be used to define a confidence interval. But,

  • the bounds for the values $\Lambda(x)$ that cover $\alpha\%$ of the distribution $\Lambda(x)$

are different from

  • the bounds for the values $\theta$ that cover $\alpha\%$ of the distribution $\mathcal{L}_{\text{relative}}(\theta)$

In the question The basic logic of constructing a confidence interval the use of the likelihood function to construct a confidence interval is discussed. A condition for the two being the same is:

When are the two methods the same?

This horizontal vs vertical is giving the same result, when the boundaries $U$ and $L$, that bound the intervals in the plot $\theta$ vs $\hat \theta$ are iso-lines for $f(\hat \theta ; \theta)$. If boundaries are everywhere at the same height than in neither of the two directions you can make an improvement.

Or in other words when $\partial_{\hat\theta}f(\hat \theta ; \theta) = -\partial_{\theta}f(\hat \theta ; \theta)$ and the likelihood distribution and fiducial distribution coincidence.

This happens in the case described utobi. In the limit of sample size $n\to\infty$ the likelihood function becomes equivalent to a shifted chi-square distribution and the bounds of the likelihood and confidence interval (based on highest density) start to coincide and the likelihood interval can be interpreted gains (approximately) the properties a confidence interval.

$\endgroup$
2
  • $\begingroup$ I believe that one might even construct a contrived example where the maximum likelihood estimate is outside to the confidence region based on the $\alpha\%$ probability mass with the highest probability density. $\endgroup$ Nov 22, 2022 at 19:36
  • 1
    $\begingroup$ Is it correct to describe $\mathcal{L}_{\text{relative}}(\theta)$ as a distribution? For that matter, I'm not sure $\Lambda(x)$ is a distribution either. Isn't it more accurate to say it's a statistic with (approximately $\chi^2$) distribution? $\endgroup$
    – dipetkov
    Dec 30, 2022 at 0:09
1
$\begingroup$

TL;DR I use an example — inference for the rate $\theta$ in a Poisson model — to show how to construct a likelihood interval. I've borrowed this example from the book In All Likelihood: Statistical Modelling And Inference Using Likelihood by Yudi Pawitan.

Say we observe $x = 3$ from $\operatorname{Poisson}(\theta)$. The likelihood interval for the rate $\theta$ is: $$ \begin{aligned} \left\{\theta:\frac{L(\theta)}{L(\widehat{\theta})} > c\right\} \end{aligned} $$ where $\widehat{\theta}$ is the maximum likelihood estimate (MLE) and the likelihood is scaled so that the maximum is 1, which occurs at $\theta=\widehat{\theta}$; $c$ is a cutoff point (between 0 and 1; to be chosen next). In his answer @SextusEmpiricus refers to $L(\theta) / L(\widehat{\theta})$ as the relative likelihood.

This formalizes the concept of the likelihood "within some percentage of the maximum value" $L(\widehat{\theta})$: the likelihood interval is the set of parameter values $\theta$ with high likelihood. But how can we decide what is "high enough", ie. how can we calibrate the likelihood?

As @utobi explains, under regularity conditions, the likelihood ratio $2\log L(\widehat{\theta}) / L(\theta)$ has (approximately) Chi-squared distribution. Therefore, when the likelihood is reasonably regular:

$$ \begin{aligned} P\left\{\frac{L(\theta)}{L(\widehat{\theta})} > c\right\} = P\left\{2\log\frac{L(\widehat{\theta})}{L(\theta)} < - 2\log c\right\} = P\left\{\chi^2_1 < -2\log c\right\} \end{aligned} $$

If we want to choose $c$, so that the likelihood interval has approximate $(1-\alpha)$100% coverage probability, then $P(\chi^2_1 < -2\log c) = 1 - \alpha \Leftrightarrow c = \exp\left\{-\frac{1}{2}\chi^2_{1,1-\alpha}\right\}$.

Here is a plot of the (scaled) likelihood $L(\theta)/L(\widehat{\theta})$ as well as the likelihood interval calibrated to have 95% coverage probability.

# significance level
alpha <- 0.05
# corresponding cutoff
c <- exp(-qchisq(1 - alpha, 1) / 2)

x <- 3
theta <- seq(0, 10, len = 10000)
loglike <- dpois(x, theta, log = TRUE)

# Compute a likelihood interval with (1-alpha)100% coverage.
# See Example 5.14, "In All Likelihood"
likelihood_interval(theta, loglike, alpha)
#>   alpha    lower    upper
#>    0.05 0.746065 7.779286

R code to reproduce the figure and calculate the likelihood interval.

find_cutoff <- function(x, y, x0) {
  if (sum(x < x0) < 2) {
    return(min(x))
  } else {
    return(approx(y[x < x0], x[x < x0], xout = 0)$y)
  }
}

# Compute likelihood interval for a scalar theta.
# This implementation is based on the program `li.r` for computing likelihood
# intervals which accompanies the book "In All Likelihood" by Yudi Pawitan.
# https://www.meb.ki.se/sites/yudpaw/book/
likelihood_interval <- function(theta, llike, alpha = 0.05) {
  # Scale the negative log-likelihood so that the minimum is 0
  nllike <- -2 * (llike - max(llike))
  # Find the MLE
  theta.mle <- mean(theta[nllike == 0])
  # Shift by the Chi-squared critical value
  nllike <- nllike - qchisq(1 - alpha, 1)

  lower <- find_cutoff(theta, nllike, theta.mle)
  upper <- -find_cutoff(-theta, nllike, -theta.mle)

  data.frame(alpha, lower, upper)
}

# significance level
alpha <- 0.05
# corresponding cutoff
c <- exp(-qchisq(1 - alpha, 1) / 2)

# An observation from Poisson(theta). We will do inference on theta.
x <- 3
theta <- seq(0, qpois(0.999, x), len = 10000)

# Compute a likelihood interval with (1-alpha)100% coverage.
# See Example 5.14, "In All Likelihood"

loglike <- dpois(x, theta, log = TRUE)
like.ci <- likelihood_interval(theta, loglike, alpha)
like.ci # approximate confidence interval based on the likelihood ratio statistic

# Scale the Poisson likelihood, so that the maximum is 1 at the MLE theta = x.
like <- exp(loglike - dpois(x, x, log = TRUE))
plot(
  theta, like,
  xlab = "θ", ylab = "L(θ)", ylim = c(0, 1),
  type = "l", xaxs = "i", yaxs = "i",  xaxt = "n"
)
axis(1, at = c(0, 2, 3, 4, 6, 8, 10), labels = c(0, 2, "", 4, 6, 8, 10))
abline(h = c, col = "gray62")
mtext("θ = 3", side = 1, at = x, padj = 1)
segments(
  like.ci$lower, c, like.ci$upper, c,
  col = "#2297E6", lwd = 3
)

title(
  paste(
    "95% likelihood interval [blue segment] for Poisson rate θ on observing x =", x
  ),
  font.main = 1
)
$\endgroup$
7
  • 1
    $\begingroup$ Hi Dipetkov, two minor points. (i) calling the relative likelihood "normalized likelihood" is confusing since the normalization of a likelihood function (to a Bayesian) means that you make it integrate to 1; this definitely doesn't integrate do 1! (ii) I'm confused by your definition "probability-based confidence interval". The CI described in my post IS the likelihood-based CI; see, e.g., Severini (2000) "Likelihood Methods in Statistics", Pace and Salvan (1997) "Principles of Statistical Inference", Azzalini (1996) "Statistical Inference Based on the likelihood", among many many others. $\endgroup$
    – utobi
    Dec 30, 2022 at 18:02
  • 1
    $\begingroup$ Re-reading your post, we are talking about the same thing, the only difference is how we choose $c$. Choosing $c = \chi_{1,\alpha}^2$ leads to a calibrated interval. I do not get what do you mean by "probability-based CI" since leaving $c$ unspecified does not really define a CI, though it does define an interval of "plausible" values for $\theta$. $\endgroup$
    – utobi
    Dec 30, 2022 at 18:10
  • 1
    $\begingroup$ the LR statistic is related to the relative likelihood by a bijection, thus any LR-based CI is necessarily a likelihood interval and the other way around. $\endgroup$
    – utobi
    Dec 30, 2022 at 19:01
  • 1
    $\begingroup$ @utobi Okay, I see what you mean about the cutoff: setting alpha = 1 - pchisq(-2*log(c), 1) makes the two intervals the same. Obviously, I need to keep reading & thinking about the textbook. Thanks for pointing this out; I'll try to fix my answer. $\endgroup$
    – dipetkov
    Dec 30, 2022 at 19:47
  • 1
    $\begingroup$ @utobi Thank you again for taking the time to write comments & help me understand this theory better. (I've already upvoted your answer also.) $\endgroup$
    – dipetkov
    Dec 31, 2022 at 21:01

Your Answer

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