Error in replaceConstantsRecurse(x, constEnv, constNames)

[Hanah]

Hi everyone,

I am fitting a spatial model as follows:

where \delta_{i} is the set of area i's neighbors, n_{\delta_{i}} is the number of neighbors.

I define the model in NIMBLE as:

####Model specification

model_code <- nimbleCode( {
  # Priors
  alpha~dflat()
  tau~dgamma(0.001,0.001)
  lamda[1:N] ~ dnorm(0,tau)
 
  #likelihood
  g <- rep(seq_along(vec1), vec1)
  sumld<--as.numeric(tapply(lamda[vec2], g, sum))
 
  for(i in 1:N) {
    u[i]<-1/nb[i]*sumld[i]
    log(theta[i]) <- alpha + u[i]
    y[i] ~ dpois(theta[i]*E[i])
  }
}
)

  ####  Data

vec1<-vec1  #this gives the number of neighbors for each area
vec2<-vec2  #this gives the list of neighbors for each area sequentially
L<-length(vec2)
N = nrow(data)
nb.data <- list(y=data$dthsim,E=data$E2020) 
nb.const <- list(N = N, L=L, vec1=vec1, vec2=vec2, nb=vec1)

set.seed(1)
#### Others
inits <- list(alpha=0,
              tau=1,
              u=rep(0,N))
params <- c("alpha", "tau", "u","lamda") # Params

#### Create the model
MR <- nimbleModel(
  code = model_code,
  name = "MR",
  constants = nb.const,
  data = nb.data,
  inits = inits)

I got the following error message while defining the model. Can someone please advise me what went wrong? Thank very much!

Best regards,

Hanah.

[Chris Paciorek]

Hi Hanah,

The specific issue is that you can't use arbitrary R code in defining your model (please see Section 5.2.5 of our manual). In particular, there's various R functions in these lines that won't work (rep, seq_along, as.numeric, tapply):

g <- rep(seq_along(vec1), vec1)
sumld<--as.numeric(tapply(lamda[vec2], g, sum))

However, you can do this with some additional work. Your case is similar to another discussion we had on the users group a while back. The basic idea is that you need to do the indexing as part of a nimbleFunction (because we don't support syntax like `sum(u[idx[1:4]])` directly in model code). 

Outside of your model code, write a nimbleFunction that does the mean calculation:

mymean <- nimbleFunction(
    run = function(lambda = double(1), nbrs = double(1), start = double(1), nb = double(1), i = double(0)) {
        return(mean(lambda[nbrs[start[i]:(start[i]+nb[i]-1)]])
        returnType(double(0))
    }
)

Then in model code:

u[i] <- mymean(lambda, nbrs, start, nb, i)

You'd need to provide nbrs, idx and nb as constants. 'nbrs' should be vector containing the deltas, all strung together in a single vector. 'start' should be a vector indicating the first delta for each i.  Suppose that the delta_1 is {1,3} and delta_2 is {1,7,11}, and delta_3 is {2, 4} and so forth. Then 'start' is {1,3,6,...}.

Alternatively, your spatial formulation is similar to a CAR model, so you may want to consider that instead. See Chapter 9 of our manual and this example.

-chris

--
You received this message because you are subscribed to the Google Groups "nimble-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nimble-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/1bd6bb70-fb67-458f-b0ba-d3fdfcf2134bn%40googlegroups.com.

[Hanah]

Hi Chris,

Thanks a lot for your detailed instructions! While it's clear to me how to proceed, my program still failed to run and the same errors persisted. I attach my script and a simulated dataset. I'd really appreciate if you could have a look.

Best regards,

Hanah. 

[Chris Paciorek]

Hi Hanah,

Two things:

1) In your declaration of 'u[i]' you need indexing on the vectors in the arguments to mymean:

    u[i]<-mymean(lamda[1:N],nbrs[1:Nnbrs],start[1:N],nb[1:N],i)

(and define Nnbrs as a constant equal to 13192)

2) Please provide nb, nbrs, and start as part of inits not as part of constants. That's what the (hard to interpret) error message is trying to indicate.

I think that should get you going.

-chris

To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/cf8c1f20-2731-451a-afac-8e6a620e629fn%40googlegroups.com.

[Hanah]

Hi Chris,

Thank you so much! I revised u[i] as you suggested and it works!

I would think that  nb, nbrs, and start are parts of constants. They are not parameters that need to be estimated, but are defined by the way I set up the neighborhood matrix. Maybe this is not clear from my previous communication, sorry for the confusion. 

Best regards,

Hanah

[Chris Paciorek]

Hi Hanah,

Yes, good point. Those values are conceptually constants in their role in the model.

As I was debugging your code, it looks like I got confused at some point and thought that for some reason that nimble needed those to be provided in inits for internal reasons, but I see now that it works if you provide them in constants (or in inits is fine too). For vector constants like these, we do create a model variable for them, but that's an internal detail not all that relevant here.

-chris

To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/948eea77-e56c-4a99-bb67-e0385fb19161n%40googlegroups.com.