[R] Help with optim() to maximize log-likelihood

Ben Bolker bbolker at gmail.com
Mon Mar 9 18:17:27 CET 2015


Sophia Kyriakou <sophia.kyriakou17 <at> gmail.com> writes:

> 
> hello, I am using the optim function to maximize the log likelihood of a
> generalized linear mixed model and I am trying to replicate glmer's
> estimated components. If I set both the sample and subject size to q=m=100
> I replicate glmer's results for the random intercept model with parameters
>  beta=-1 and sigma^2=1. But if I change beta to 2 glmer works and optim
> gives me the error message "function cannot be evaluated at initial
> parameters".
> 
> If anyone could please help?
> Thanks

 snip to make gmane happy.

It looks like you're getting floating-point under/overflow.  If you do
all the computations on the log scale first and then exponentiate,
it seems to work, i.e.:

        piYc_ir[i,] <- lchoose(m,Y[i]) + Y[i]*(z+beta) + (-z^2/(2*exp(psi))) - 
            m*(log1p(exp(z+beta))) - 0.5*(log(2*pi)+psi)
        piYc_ir[i,] <- exp(piYc_ir[i,])

follow-ups should probably go to r-sig-mixed-models at r-project.org
instead ...

  Ben Bolker



More information about the R-help mailing list