if else statement errors

[Steven Forrest]

Hello,

I have a multi year model of count data and wrote code to condition on the data if there are NA's present or data and how the transitions occur. I wrote a model without conditionals for this same year with for loops in sequence and Nimble performed wonderfully. However I would prefer not to write unique coding for each year and this code should parse each year correctly. The issue is when I run this I get this specific error:

> constants <- list(Days = length(y), IS.NA = is.na(y)*1)
> Test_Model <- nimbleModel(Test_Code, inits = Test_Inits, data = list( y = y))
Defining model
Error in eval(code[[2]], constantsEnv) : object 'i' not found
Error in codeProcessIfThenElse(code[[i]], constants, envir) :
  Cannot evaluate condition of 'if' statement: IS.NA[i] == 1 & IS.NA[i - 1] == 1.
Condition must be able to be evaluated based on values in 'constants'.

This is the code with the conditionals. Im conditioning on  "IS.NA "which is a vector of 0's and 1's corresponding to NA's and data. I also conditioned direct on the data and receiced the same error.

Test_Code <- nimbleCode({
 
  n[1] <- z[1]
  if(is.na(y[1])){
     y[1] <- n[1]
  }else{  
     y[1] ~ dbin(p,n[1])
    }
 
  for(i in 2:Days){
    if(IS.NA[i]==1  & IS.NA[i-1]==1){
       n[i] <- z[i]
       y[i] <- n[i]
    }else if (IS.NA[i] != 1  & IS.NA[i-1] == 1) {
       n[i] <- z[i]
       y[i] ~ dbin(p,n[i])
    }else if(IS.NA[i] != 1 & IS.NA[i-1] != 1 ) {
       n[i] <- n[i-1]-y[i-1]+z[i]
       y[i]  ~ dbin(p,n[i])
    }else {
       n[i] <- n[i-1]-y[i-1]+z[i]
       y[i] <- n[i]
   }
  }
  for(i in 1:Days){
    z[i] ~ dnegbin(pz[i],r)
    pz[i] <- r/(r+mu[i])    
    mu[i] <- alpha * exp((m-i)/beta)*exp(-exp((m-i)/beta)) / beta
  }
   
 
  alpha <- exp(lnalpha)
  lnalpha ~ dnorm(0,sd=32)
 
 
  beta ~ dunif( 1,12 )
  m ~ dunif( 5,35 )
   
  r ~ dexp(0.001)
 
  p ~ dbeta(167,182)
  esc <- sum(y[1:48])
  })
 



y = c(NA,NA,NA,1,1,0,1,14,10,16,
      117,28,35,71,362,679,464,NA,NA,NA,
      NA,NA,NA,NA,725,754,NA,NA,NA,379,
      300,341,483,350,250,  357,139,53,54,96,
      142,NA,NA,NA,NA,54,34,29)
y

IS.NA <- is.na(y)*1.0

Test_Inits <- list(p=0.5,lnalpha=9.8, m=20.7 ,beta=7.3, r = 0.1,
                   z=c(1,1,1,2,1, 0,1,27,6,22,
                       140,17,42,107,653,  996,249,536,1000,1000,
                       1000,1000,1000,1000,1450,  783,246,1000,1000,758,
                       221,382,625,217,150,  226,66,60,55,138,
                       188,58,200,200,200,  108,14,24)
)
constants <- list(Days = length(y), IS.NA = is.na(y)*1.0)

Test_Model <- nimbleModel(Test_Code, inits = Test_Inits, data = list( y = y))

Thank you

Steve Forrest

[Perry de Valpine]

Hi Steven,

Thanks for the question. The way you are trying to write model code is not supported. The reason is that model code is declarative (order doesn't matter) and is processed to establish the directed acyclic graph of relationships among model variables. if-then-else statements are not supported because they imply potentially changing or branching DAG relationships. Given that your IS.NA vector is constant, there are several ways you could achieve what you want. One version would be with nested indexing. If I call your four cases A, B, C and D and you precompute their indices in advance and provide them in constants, you can do:

for(i in 1:n_caseA) {

  n[caseA_inds[i]] <- z[caseA_inds[i]]

  y[caseA_inds[i]] <- n[caseA_inds[i]]

}

for(i in 1:n_caseB) {

  n[caseB_inds[i]] <- z[caseB_inds[i]]

  y[caseB_inds[i]] ~ dbin(p,n[caseB_inds[i]])

}

etc.

Another version would be to define a nimbleFunction with some of the case splitting. However that would be trickier since sometimes you want y[i] to be stochastic and sometimes deterministic, so you'd still need multiple cases directly in the model code. Also the nested indexing approach above has the benefit that it tells the DAG specifically which relationships are where in the model.

A final relevant trick is to use logicals in arithmetic, like this:

n[i] <- z[i] + (IS.NA[i-1]!=1)*(n[i-1]-y[i-1])

which looks like it covers your cases for n[i]. A disadvantage of this approach is that in the DAG, when algorithms are determining which nodes depend on which, it will always look like n[i] depends on n[i-1] and y[i-1] even when they don't really.

HTH

Perry

--
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 visit https://groups.google.com/d/msgid/nimble-users/3754409c-8bf6-43b1-b077-d8595c371529n%40googlegroups.com.