Maximum Likelihood Estimate for Binomial Data

Question

I need to find the maximum likelihood estimate for a vector of binomial data.

one like this:

binvec <- rbinom(1000, 1, 0.5)

I tried to first create the function and then optimize it with optim().

llbinom <- function(theta, x) {return(sum(dbinom(x=x,
               size = theta[1],
               prob = theta[2],log=TRUE)))}
mybin <- optim(c(0.5,0.5),fn=llbinom,x=binvec)
mybin

I do get some result but also error messages that NaNs are being produced and the function cannot be evaluated at the initial parameters. I constructed it from an example that works for normally distributed data and believe I made a mistake in converting it.

Here, the original code I got:

ll <- function(theta,x) {return(-sum(dnorm(x=x,
   mean=theta[1],sd=theta[2],log=TRUE)))}
mle <- optim(c(5,3),fn=ll,x=binvec)

Ben Bolker Ben Bolker · Accepted Answer · 2018-12-04T23:31:56

Several issues here.

looks like you're missing a negative sign (optim() minimizes by default unless you set the control parameter fnscale=-1, so you need to define a negative log-likelihood function)
the size parameter must be an integer
it's unusual, and technically challenging, to to estimate the size parameter from data (this is often done using N-mixture models, if you want to read up on the technique: e.g. see the unmarked package); usually the number of trials is assumed to be known. So I would try

llbinom <- function(theta, x) {return(
   -sum(dbinom(x=x,
               size = 1,
               prob = theta[1],log=TRUE)))}
mybin <- optim(c(0.5),fn=llbinom,x=binvec)

There are plenty of reasons to do this the Hard Way (numerically); if you really only need to find the MLE of the probability of a single binomial sample x (independent observations with the same probability of success out of s trials), the analytical solution is sum(x)/sum(s) ...

Maximum Likelihood Estimate for Binomial Data

1 Answers