NaNs produced by dbinom within CAR hierarchical model (reproducible minimal example attached)

[Felipe Montealegre]

Hi!

I've been trying to use an N-mixture model with a conditional autoregressive (CAR) component (essentially like Model 3 in: doi:10.1111/2041-210X.12983). I started running into error messages I couldn't quite diagnose so I wrote a minimal example model that reproduces the problem. I'll attach the full script for reproduction, but here is the model code:

```

code <- nimbleCode({
 
  d[1:nsites] ~ dcar_proper(
    mu=mu[1:nsites],
    adj=adj[1:len.adj],
    num=n.neigh[1:nsites],
    tau=1,
    gamma=0.0
    # between -1 and 1 since we omit M, C (see help(dcar_proper))
  )
  effect ~ dnorm(0,1)
 
  for(i in 1:nsites){
    counts[i] ~ dbinom(exp(d[i] + effect), 0.8)
  }
})

```

Here, counts is the only data node for the model.

When I attempt to run Markov chains for this model on a simple hand-made dataset, I run into the problem of getting `logProb_counts = -Inf`:

```

Warning messages:

1: In dbinom(model$counts[getNodeFunctionIndexedInfo(INDEXEDNODEINFO_,  : NaNs produced

[...]

```

Normally, I know that this can happen because some variable wasn't initialized and shouldn't give problems, but the problem is that the chains come out 'flat':

```

$chain1
Markov Chain Monte Carlo (MCMC) output:
Start = 1
End = 4000
Thinning interval = 1
        d[1] d[2] d[3] d[4] d[5] effect
   [1,]    2    2    2    2    2      0
   [2,]    2    2    2    2    2      0
   [3,]    2    2    2    2    2      0
   [4,]    2    2    2    2    2      0
   [5,]    2    2    2    2    2      0
   [6,]    2    2    2    2    2      0
   [7,]    2    2    2    2    2      0

```

So something has to be wrong with the way I specified the model, but I can't figure out what's going wrong. Does anyone happen to know what might be up? My guess is that it has to do with the CAR component, since it's my first time writing a CAR model. But I'm not sure, I might have just overlooked something very basic too haha.

Thank you so much in advance!!

[Felipe Montealegre]

I kinda realized it myself - this is embarassing haha but the 'size' parameter in the binomial distribution has to be an integer and in my model that wasn't the case. If I instead use lambda[i] <- exp(d[i] + effect), N[i] ~ dpois(lambda[i]), counts[i] ~ binom(N[i], 0.8) the problem goes away.

For completeness in case someone runs into the same issue in the future and googles this, I also tried rounding exp(d[i] + effect) to the nearest integer but I couldn't get that to work. I don't think that rounding would be a reasonable thing for a model in most situations anyway, but I just wanted to have my afternoon headache potentially help someone in the future haha!