Sergio Nolazco - 3 years ago 217

R Question

I am using

`dnbinom()`

`mle2()`

The problem is that I got 16 warnings for my negative binomial model, all of them NaNs produced like this one:

1: In dnbinom(y, mu = mu, size = k, log = TRUE) : NaNs produced

My code:

`# data`

x <- c(0.35,0.45,0.90,0.05,1.00,0.50,0.45,0.25,0.15,0.40,0.26,0.37,0.43,0.34,0.00,0.11,0.00,0.00,0.00,0.41,0.14,0.80,0.60,0.23,0.17,0.31,0.30,0.00,0.23,0.33,0.30,0.00,0.00)

y <- c(1,10,0,0,67,0,9,5,0,0,0,82,36,0,32,7,7,132,14,33,0,67,11,39,41,67,9,1,44,62,111,52,0)

# log-likelihood function

negbinglmLL = function(beta,gamma,k) {

mu= exp(beta+gamma*x)

-sum(dnbinom(y,mu=mu, size=k, log=TRUE))

}

# maximum likelihood estimator

model <- mle2(negbinglmLL, start=list(beta=mean(y), gamma= 0, k=mean(y)^2/(var(y)-mean(y))))

What do these warnings mean, and if this is a serious problem how can I avoid it?

Recommended for you: Get network issues from **WhatsUp Gold**. **Not end users.**

Answer Source

You're not restricting the negative log-likelihood function from trying negative values of `k`

. This *probably* doesn't mess up your final answer, but it's always best to avoid these kinds of warnings if you can. Two simple strategies:

- put a lower bound on
`k`

(switching to`method=L-BFGS-B`

) - fit the
`k`

parameter on the log scale, as follows:

```
negbinglmLL = function(beta,gamma,logk) {
mu= exp(beta+gamma*x)
-sum(dnbinom(y,mu=mu, size=exp(logk), log=TRUE))
}
model <- mle2(negbinglmLL,
start=list(beta=mean(y),
gamma= 0,
logk=log(mean(y)^2/(var(y)-mean(y)))))
```

By the way, for simple problems like this you can use a formula-based shortcut as follows:

```
mle2(y~dnbinom(mu=exp(logmu),size=exp(logk)),
parameters=list(logmu~x),
start=list(logmu=0,logk=0),
data=data.frame(x,y))
```

For this simple case `MASS::glm.nb`

should also work perfectly well (but perhaps this is the simplest version of something that will get more complicated/beyond the scope of `glm.nb`

).

Recommended from our users: **Dynamic Network Monitoring from WhatsUp Gold from IPSwitch**. ** Free Download**