Profile Likelihood

Consider some simulated data

> set.seed(1)
> x=exp(rnorm(100))

Assume that those data are observed i.id. random variables with distribution$F_{\boldsymbol{\theta}}$, with $\boldsymbol{\theta}=(\alpha,\beta)$. The natural idea is to consider the maximum likelihood estimator

$\widehat{\boldsymbol{\theta}}=\underset{\boldsymbol{\theta}}{\text{argmax}}\left\lbrace\log \mathcal{L}(\boldsymbol{\theta})\right\rbrace$For instance, consider some maximum likelihood estimator,

> library(MASS)
> (F=fitdistr(x,"gamma"))
shape       rate
1.4214497   0.8619969
(0.1822570) (0.1320717)
> F$estimate[1]+c(-1,1)*1.96*F$sd[1]
[1] 1.064226 1.778673

Here, we have an approximated (since the maximum likelihood has an asymptotic Gaussian distribution) confidence interval for $\alpha$. We can use numerical optimization routine to get the maximum of the log-likelihood function

> log_lik=function(theta){
+   a=theta[1]
+   b=theta[2]
+   logL=sum(log(dgamma(x,a,b)))
+   return(-logL)
+ }

> optim(c(1,1),log_lik)
$par [1] 1.4214116 0.8620311$value
[1] 146.5909

And we have the same value.

Now, what if we care only about $\alpha$, and not $\beta$. The we can use profile likelihood. The idea is to solve

$\widehat{\alpha}=\underset{\alpha}{\text{argmax}}\left\lbrace \underset{\beta}{\text{max}}\left\lbrace\log \mathcal{L}(\alpha,\beta)\right\rbrace \right\rbrace$

i.e.

$\widehat{\alpha}=\underset{\alpha}{\text{argmax}}\left\lbrace \log \mathcal{L}_p(\alpha)\right\rbrace$or, equivalently,

$\widehat{\alpha}=\underset{\alpha}{\text{argmax}}\left\lbrace \log \mathcal{L}(\alpha,\widehat{\beta}_{\alpha})\right\rbrace$

> prof_log_lik=function(a){
+   b=(optim(1,function(z) -sum(log(dgamma(x,a,z)))))$par + return(-sum(log(dgamma(x,a,b)))) + } > vx=seq(.5,3,length=101) > vl=-Vectorize(prof_log_lik)(vx) > plot(vx,vl,type="l") > optim(1,prof_log_lik)$par
[1] 1.421094

$value [1] 146.5909 A few weeks ago, we have mentioned the likelihood ratio test, i.e. $2[\log\mathcal{L}(\widehat{\theta})-\log\mathcal{L}(\theta_0)]\sim n[\widehat{\theta}-{\theta}_0]^2I(\theta_0)\sim \chi^2(k)$ The analogous can be obtained here, since $2[\log\mathcal{L}_p(\widehat{\alpha})-\log\mathcal{L}_p(\alpha_0)]\sim\chi^2(1)$ (the 1 comes from the fact that $\alpha$ is a one-dimensional coefficient). The (technical) proof can be found in Suhasini Subba Rao’s notes (see also Section 4.5.2 in Antony Davison’s Statistical Models). From that property, we can easily obtain a confidence interval for $\alpha$ Hence, from our sample, we get the following 95% confidence interval, > abline(v=optim(1,prof_log_lik)$par,lty=2)
> abline(h=-optim(1,prof_log_lik)$value) > abline(h=-optim(1,prof_log_lik)$value-qchisq(.95,1)/2)

> segments(F$estimate[1]-1.96*F$sd[1],
-170,F$estimate[1]+1.96*F$sd[1],-170,lwd=3,col="blue")
> borne=-optim(1,prof_log_lik)$value-qchisq(.95,1)/2 > (b1=uniroot(function(z) Vectorize(prof_log_lik)(z)+borne,c(.5,1.5))$root)
[1] 1.095726
> (b2=uniroot(function(z) Vectorize(prof_log_lik)(z)+borne,c(1.25,2.5))\$root)
[1] 1.811809

that can be visualized below,

> segments(b1,-168,b2,-168,lwd=3,col="red")

In blue the obtained obtained using the asymptotic Gaussian property of the maximum likelihood estimator, and in red, the obtained obtained using the asymptotic chi-square distribution of the log (profile) likelihood ratio.

Cite this article as: Arthur Charpentier, "Profile Likelihood," in Freakonometrics, 16/11/2015, https://freakonometrics.hypotheses.org/20573.

2 thoughts on “Profile Likelihood”

1. Aminath says:

A good example. I am analyzing a system of nonlinear ode and need to compute profile likelihood. Do you have a similar example for ode. Or would you be able to recommend a reference with r code.
Thanks

This site uses Akismet to reduce spam. Learn how your comment data is processed.