What's causing rjags "Node inconsistent with parents" error

44 Views Asked by At

I'm trying to run a bird occupancy model using rjags and I'm getting an error saying: Error in node ytb[1,2,9,1:5,1] Node inconsistent with parents. The vector that the error references is a vector of 5 0s because it represents a survey where no birds of species h were detected in any of the time bins. I've defined my initial values as follows:

Zin <- apply(jags.data$y, c(1,2,3), function(x) max(x,na.rm=T))
Zin[Zin == -Inf] <- NA # replace -inf (generated from max() function with NA)

inits<-function(){
  list(z=Zin,
       mean_alpha0.psi = runif(1, -3, 3),
       sd_alpha0 = runif(1),
       sd_eps = runif(1),
       mean_phi0=rnorm(nspeciesSub),
       sd_phi0=runif(nspeciesSub)
  )
} 

Zin is a matrix where Zin[h,i,j] = 1 if there was ever an observation of bird species h at property i and point j.

ytb is referenced in this line of my model: ytb[h, i, j, 1:ntbins, v] ~ dmulti(pi_pa_normalized[h, i, j, 1:ntbins, v], z[h,i,j])

Here is the full model:

model{
#-----------------------------------------------------------------------------------------------------
# 1) Priors
#-----------------------------------------------------------------------------------------------------
  
  ######## Abundance ########
  mean_alpha0.psi ~ dnorm(0, 0.01)   # Mean intercept for occupancy probability for all species
  sd_alpha0 ~ dunif(0, 10)      # Standard deviation on mean abundance for all species
  prec_alpha0 <- 1 / (sd_alpha0 ^ 2)
  
  for (h in 1:nspecies) {
     alpha0.psi[h] ~ dnorm(mean_alpha0.psi, prec_alpha0) # Mean log abundance per point for species h
    
    # Uninformative priors on covariates
    alpha1[h] ~ dnorm(mu_alpha1, tau_alpha1) # Effect of fire frequency for each species
    alpha2[h] ~ dnorm(mu_alpha2, tau_alpha2) # Effect of ground cover for each species
    alpha3[h] ~ dnorm(mu_alpha3, tau_alpha3) # Effect of vegetation height for each species
    alpha4[h] ~ dnorm(mu_alpha4, tau_alpha4) # Effect of total basal area for each species
    alpha5[h] ~ dnorm(mu_alpha5, tau_alpha5) # Effect of % deciduous area for each species
  } # close h
  
  sd_eps ~ dexp(1) # Standard deviation for random property effect
  prec_eps <- 1 / (sd_eps ^ 2) # Precision for random property effect
  
  for (i in 1:nprops) {
    eps[i] ~ dnorm(0, prec_eps) # Random property effect
  } # close i
  
  mu_alpha1 ~ dnorm(0, 0.001)
  tau_alpha1 ~ dgamma(0.01, 0.01)
  
  mu_alpha2 ~ dnorm(0, 0.001)
  tau_alpha2 ~ dgamma(0.01, 0.01)
  
  mu_alpha3 ~ dnorm(0, 0.001)
  tau_alpha3 ~ dgamma(0.01, 0.01)
  
  mu_alpha4 ~ dnorm(0, 0.001)
  tau_alpha4 ~ dgamma(0.01, 0.01)
  
  mu_alpha5 ~ dnorm(0, 0.001)
  tau_alpha5 ~ dgamma(0.01, 0.01)
  
  ######## Availability #########
  for (h in 1:nspecies) {
    mean_phi0[h] ~ dnorm(0, 0.001)    # Mean availability for each species
    sd_phi0[h] ~ dexp(1)              # Standard deviation of availability for each species
    prec_phi0[h] <- 1 / (sd_phi0[h] ^ 2)
    
    for (i in 1:nprops) {
      for (j in 1:npoints_site[i]) {
        phi0[h, i, j] ~ dnorm(mean_phi0[h], prec_phi0[h]) # Mean availability for species h at property i and point j
      } # close j
    } # close i
  } # close h
  
  phi1 ~ dnorm(0, 0.001) # coefficient of day of year in availiabity linear model
  phi2 ~ dnorm(0, 0.001) # coefficient of day of year^2 in availiabity linear model
  phi3 ~ dnorm(0, 0.001) # coefficient of time after sunrise in availiabity linear model
  phi4 ~ dnorm(0, 0.001) # coefficient of temp in availiabity linear model
  phi5 ~ dnorm(0, 0.001) # coefficient of wind in availiabity linear model
  
  
  
#-----------------------------------------------------------------------------------------------------
# 2) Ecological process model
#-----------------------------------------------------------------------------------------------------

  for (h in 1:nspecies) {
    for (i in 1:nprops) {
      for (j in 1:npoints_site[i]) {
        z[h,i,j] ~ dbern(psi[h,i,j]) # True occupancy based on occupancy probability
        psi[h,i,j] <- 1 / (1 + exp(-lpsi.lim[h,i,j]))
        lpsi.lim[h,i,j] <- min(999, max(-999, lpsi[h,i,j]))
        lpsi[h,i,j] <- alpha0.psi[h] + alpha1[h] * fire[i, j] + alpha2[h] * ground[i, j] + alpha3[h] *
          height[i, j] + alpha4[h] * tba[i, j] + alpha5[h] * decid[i, j] + eps[i]
       
        
        
        
#-----------------------------------------------------------------------------------------------------
# 3) Observational process model
#-----------------------------------------------------------------------------------------------------

        ###### Availability with time removal ######
        for (v in 1:vh[i, j]) {
          # loop over each visit (v) to each site
          
          # Availability is a function of intercept + covariates
          
          y[h,i,j,v] ~ dbern(mu.p[h,i,j,v]) # Detection/non-detection
          mu.p[h,i,j,v] <- z[h,i,j] * pa[h,i,j,v] # mu.p conditional on bird being present at site
          
          logit(p[h,i,j,v]) <- phi0[h, i, j] + phi1 * obs_date[i, j, v] + phi2 * obs_date[i, j, v] * obs_date[i, j, v] + phi3 * obs_time[i, j, v] + phi4 * obs_temp[i, j, v] + phi5 * obs_wind[i, j, v]
          
          for (z in 1:ntbins) {
            # Probability of being available in each time bin
            # pi_pa in time bin >1 conditional on not being detected before
            pi_pa[h, i, j, z, v] <-
              p[h, i, j, v] * pow(1 - p[h, i, j, v], (z - 1))
            pi_pa_normalized[h, i, j, z, v] <-
              pi_pa[h, i, j, z, v] / pa[h, i, j, v]
          } # close z
          
          # Probability ever available
          pa[h, i, j, v] <- sum(pi_pa[h, i, j, 1:ntbins, v])
          
          # Time data - when were detections mentally removed?
          ytb[h, i, j, 1:ntbins, v] ~ dmulti(pi_pa_normalized[h, i, j, 1:ntbins, v], z[h,i,j])
          
        } # close v
      } # close j
    } # close i
  } # close h
}

I've messed with the initial values for z, making them all 1s or all 0s, in addition to what I've tried written above, but I always run into an error saying either "invalid parent node" or "node inconsistent with parents."

0

There are 0 best solutions below