Tag Archives: maximum likelihood

M1-statistique

Re-parametrization and Maximum Likelihood

3 Comments

The maximum likelihood estimator is invariant in the sense that for all bijective function  , if   is the maximum likelihood estimator of   then  . Let  , then   is equal to  , and the likelihood function in   is  . And since   is the maximum likelihood estimator of  ,

hence,   is the maximum likelihood estimator of  .

For instance, the Bernoulli distribution is   with   and

 

Given sample  , the likelihood is

 

The log-likelihood is then

 

with ICI

 https://latex.codecogs.com/gif.latex?\frac{\partial}{\partial%20p}\log\mathcal{L}(p)=\frac{\sum%20x_i}{p}-\frac{n-\sum%20x_i}{1-p}.

Thus, the first order condition

 https://latex.codecogs.com/gif.latex?\frac{\partial}{\partial%20p}\log\mathcal{L}(p)=0

is satisfied when  . In order to illustrate, consider the following data


> set.seed(1)
> X=sample(0:1,size=15,replace=TRUE)
> X
[1] 0 0 1 1 0 1 1 1 1 0 0 0 1 0 1

The (negative) log-likelihood is here


> loglik=function(p){
+ -sum(log(dbinom(X,size=1,prob=p)))
+ }

that we can visualize below


> u=seq(0,1,by=.025)
> v=-Vectorize(loglik)(u)
> plot(u,v,type="l",xlab="",ylab="")

From calculations above, we know that the maximum likelihood estimator for  is


> mean(X)
[1] 0.5333333

The numerical version is


> (opt=optim(.5,loglik))
$par
[1] 0.5333008

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL

Somehow, we were lucky here, because we did not say that the optimization was on the interval . Nevertheless, our estimator for the probability belongs to . In order to insure that the optimal value is in , we can consider some constrained optimization routine


> constrOptim(.5, loglik, grad=NULL,ui=matrix(c(1,-1),2,1), ci=c(0,-1))
$par
[1] 0.5333008

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL

$outer.iterations
[1] 2

$barrier.value
[1] 6.909277e-05

On the previous graph, we did – indeed – reach that maximum of the log-likelihood


> abline(v=opt$par,col="red")

An alternative is to consider   (as in the exponential family). The log-likelihood is then

 

since

 

Here

 

Thus, the first order condition

 

is satisfied when

i.e.

 

From a numerical perspective, we have the same optimal value


> loglik=function(theta){
+ -sum(log(dbinom(X,size=1,prob=exp(theta)/(1+exp(theta)))))
+ }
> (opt=optim(0,loglik))
$par
[1] 0.1335938

$value
[1] 10.36385

$counts
function gradient
20 NA

$convergence
[1] 0

$message
NULL
> exp(opt$par)/(1+exp(opt$par))
[1] 0.5333489