I keep getting the above error when I try adding detection covariates on a beta-binomial N-mixture model in rags. According to Royle(2004). A binomial N mixture model can be used to model abundance data arising from repeat count surveys. The number of individuals on a site can be modeled by a Poisson [For simplicity I will stick to the Poisson model only] such that;
Ni ~ Poisson(λi)
yit ~ bin(pi't,Ni)
Ni - is the number of animals available in site i
yit - is the count of observed animals at site i visit t
λi - is the average number of animals at site i
pit - is the mean detection probability.
Covariate effects can be modeled as:
Abundance:
log(λi)= B0+ B1xi1 +...+Brxir where 1...r are covariates
detection:
logit(pit)= B0+ B1xi1 +...+Brxir where 1...r are covariates
Probability of detection pit is assumed to be constant for all present animals.
The *beta binomial model ease this assumption by letting detection probability follow a stochastic distribution instead, such that
pi't ~ Beta(pit(1-δ2)/δ2,(1-pit)(1-δ2)/δ2
for 0<δ<1
I tried implementing the model with a simulated data:
20 sites, 5 visits, site covariate = Location, and 2 observed covariates.
simulated data
library(modelr)
Location<-c("A","B","C","D")
Location<-data.frame(Location=rep(Location,5))
location=Location%>%model_matrix(~Location)%>%select(-1)
set.seed(100)
y<-matrix(rpois(100,0.5),ncol=5)
# Cov1
set.seed(100)
cov1<-matrix(rnorm(100,100,5),ncol=5)
# cov 2
set.seed(100)
cov2<-matrix(rnorm(100,50,2),ncol=5)
data<-list(y=y,
nSites=20,
nOcc=5,
nA=ncol(data$location),
location=location,
cov1=cov1,
cov2=cov2)
if i try estimating this model in rjags without the covariate effects for detection it works.
nx<-"
model{
for(i in 1:nSites) {
# Biological model
N[i] ~ dpois(lambda[i])
log(lambda[i])<-alphao+inprod(alpha, location[i,])
}
# Observation model
for(i in 1:nSites) {
for(t in 1:nOcc) {
y[i,t] ~ dbin(pit_h[i,t], N[i])
pit_h[i,t]~ dbeta((pit[i,t]*fac),(fac*(1-pit[i,t])))
}
}
# Priors
alphao ~ dnorm(5,1)
for(i in 1:nA){
alpha[i] ~ dnorm(2,.5)
}
for(i in 1:nSites) {
for(t in 1:nOcc) {
pit[i,t]~dunif(0,1)
}
}
sigma ~ dunif(0,1)
fac <-(1-sigma^2)/sigma^2
}"
writeLines(nx,con="mod.txt")
inits = function() list(N = apply(y,1,max,na.rm=T))
watch=c("alphao","alpha","lambda","pit")
set.seed(100)
mod<-jagsUI::jags(data,
parameters.to.save=watch,
model.file="mod.txt",n.iter=3,
n.chains=2)
mod
yielding
JAGS output for model 'mod.txt', generated by jagsUI.
Estimates based on 2 chains of 3 iterations,
adaptation = 100 iterations (sufficient),
burn-in = 0 iterations and thin rate = 1,
yielding 6 total samples from the joint posterior.
MCMC ran for 0.009 minutes at time 2023-01-11 16:52:03.
mean sd 2.5% 50% 97.5% overlap0 f Rhat n.eff
alphao 4.769 0.446 4.350 4.754 5.218 FALSE 1 27.188 2
alpha[1] 1.725 0.629 1.115 1.735 2.310 FALSE 1 48.915 2
alpha[2] 1.991 0.161 1.795 2.012 2.149 FALSE 1 7.456 2
alpha[3] 1.822 0.633 1.186 1.852 2.415 FALSE 1 28.072 2
lambda[1] 127.706 54.181 77.504 124.146 184.586 FALSE 1 19.433 2
lambda[2] 670.612 122.557 551.850 671.722 785.277 FALSE 1 48.321 2
lambda[3] 893.022 252.923 660.011 887.089 1138.013 FALSE 1 50.136 2
lambda[4] 739.326 137.772 604.421 735.512 877.068 FALSE 1 19.180 2
lambda[5] 127.706 54.181 77.504 124.146 184.586 FALSE 1 19.433 2
lambda[6] 670.612 122.557 551.850 671.722 785.277 FALSE 1 48.321 2
lambda[7] 893.022 252.923 660.011 887.089 1138.013 FALSE 1 50.136 2
lambda[8] 739.326 137.772 604.421 735.512 877.068 FALSE 1 19.180 2
lambda[9] 127.706 54.181 77.504 124.146 184.586 FALSE 1 19.433 2
lambda[10] 670.612 122.557 551.850 671.722 785.277 FALSE 1 48.321 2
lambda[11] 893.022 252.923 660.011 887.089 1138.013 FALSE 1 50.136 2
lambda[12] 739.326 137.772 604.421 735.512 877.068 FALSE 1 19.180 2
lambda[13] 127.706 54.181 77.504 124.146 184.586 FALSE 1 19.433 2
lambda[14] 670.612 122.557 551.850 671.722 785.277 FALSE 1 48.321 2
lambda[15] 893.022 252.923 660.011 887.089 1138.013 FALSE 1 50.136 2
lambda[16] 739.326 137.772 604.421 735.512 877.068 FALSE 1 19.180 2
lambda[17] 127.706 54.181 77.504 124.146 184.586 FALSE 1 19.433 2
lambda[18] 670.612 122.557 551.850 671.722 785.277 FALSE 1 48.321 2
lambda[19] 893.022 252.923 660.011 887.089 1138.013 FALSE 1 50.136 2
lambda[20] 739.326 137.772 604.421 735.512 877.068 FALSE 1 19.180 2
pit[1,1] 0.183 0.118 0.048 0.157 0.327 FALSE 1 0.946 6
pit[2,1] 0.267 0.238 0.036 0.218 0.576 FALSE 1 4.558 2
pit[3,1] 0.280 0.143 0.085 0.313 0.432 FALSE 1 0.929 6
pit[4,1] 0.354 0.236 0.045 0.396 0.622 FALSE 1 1.186 6
pit[5,1] 0.199 0.100 0.082 0.190 0.346 FALSE 1 1.128 6
pit[6,1] 0.130 0.076 0.037 0.118 0.233 FALSE 1 3.429 2
pit[7,1] 0.503 0.197 0.209 0.536 0.741 FALSE 1 1.120 6
pit[8,1] 0.369 0.199 0.076 0.414 0.597 FALSE 1 0.986 6
pit[9,1] 0.396 0.131 0.224 0.425 0.551 FALSE 1 1.751 4
pit[10,1] 0.281 0.141 0.122 0.271 0.468 FALSE 1 1.179 6
pit[11,1] 0.554 0.226 0.291 0.585 0.768 FALSE 1 1.031 6
pit[12,1] 0.304 0.165 0.139 0.296 0.558 FALSE 1 1.001 6
pit[13,1] 0.240 0.274 0.071 0.139 0.717 FALSE 1 1.648 4
pit[14,1] 0.199 0.111 0.074 0.197 0.332 FALSE 1 4.464 2
pit[15,1] 0.322 0.093 0.207 0.303 0.446 FALSE 1 0.849 6
pit[16,1] 0.380 0.226 0.060 0.381 0.690 FALSE 1 1.043 6
pit[17,1] 0.085 0.045 0.035 0.084 0.133 FALSE 1 1.130 6
pit[18,1] 0.193 0.109 0.066 0.217 0.335 FALSE 1 2.724 3
pit[19,1] 0.135 0.047 0.081 0.137 0.203 FALSE 1 1.110 6
pit[20,1] 0.370 0.219 0.099 0.368 0.702 FALSE 1 1.063 6
pit[1,2] 0.419 0.168 0.285 0.362 0.707 FALSE 1 1.275 6
pit[2,2] 0.535 0.287 0.250 0.502 0.869 FALSE 1 8.332 2
pit[3,2] 0.356 0.206 0.095 0.362 0.662 FALSE 1 1.030 6
pit[4,2] 0.330 0.123 0.177 0.335 0.486 FALSE 1 4.460 2
pit[5,2] 0.334 0.215 0.072 0.350 0.617 FALSE 1 1.160 6
pit[6,2] 0.032 0.026 0.005 0.023 0.072 FALSE 1 2.278 3
pit[7,2] 0.385 0.289 0.035 0.438 0.682 FALSE 1 1.426 5
pit[8,2] 0.537 0.162 0.378 0.488 0.778 FALSE 1 0.881 6
pit[9,2] 0.489 0.266 0.073 0.539 0.741 FALSE 1 1.344 5
pit[10,2] 0.194 0.182 0.024 0.180 0.414 FALSE 1 1.166 6
pit[11,2] 0.476 0.258 0.236 0.424 0.811 FALSE 1 6.065 2
pit[12,2] 0.536 0.225 0.232 0.610 0.744 FALSE 1 2.909 3
pit[13,2] 0.244 0.090 0.141 0.239 0.375 FALSE 1 3.298 2
pit[14,2] 0.432 0.175 0.257 0.403 0.660 FALSE 1 1.097 6
pit[15,2] 0.419 0.287 0.122 0.404 0.738 FALSE 1 8.309 2
pit[16,2] 0.522 0.146 0.378 0.502 0.744 FALSE 1 2.006 3
pit[17,2] 0.225 0.167 0.041 0.176 0.449 FALSE 1 3.682 2
pit[18,2] 0.264 0.079 0.164 0.265 0.356 FALSE 1 1.245 6
pit[19,2] 0.440 0.243 0.161 0.466 0.766 FALSE 1 3.425 2
pit[20,2] 0.238 0.139 0.099 0.216 0.446 FALSE 1 2.344 3
pit[1,3] 0.273 0.159 0.064 0.263 0.484 FALSE 1 1.031 6
pit[2,3] 0.332 0.115 0.200 0.327 0.497 FALSE 1 3.520 2
pit[3,3] 0.533 0.251 0.183 0.494 0.840 FALSE 1 1.101 6
pit[4,3] 0.324 0.250 0.117 0.205 0.685 FALSE 1 0.865 6
pit[5,3] 0.607 0.221 0.224 0.674 0.742 FALSE 1 1.470 5
pit[6,3] 0.298 0.113 0.160 0.293 0.461 FALSE 1 1.069 6
pit[7,3] 0.403 0.143 0.163 0.429 0.526 FALSE 1 1.613 4
pit[8,3] 0.415 0.170 0.261 0.363 0.682 FALSE 1 3.085 2
pit[9,3] 0.498 0.321 0.099 0.594 0.861 FALSE 1 3.386 2
pit[10,3] 0.258 0.222 0.055 0.185 0.611 FALSE 1 0.970 6
pit[11,3] 0.381 0.268 0.058 0.360 0.789 FALSE 1 1.756 4
pit[12,3] 0.162 0.072 0.084 0.159 0.268 FALSE 1 1.566 4
pit[13,3] 0.152 0.159 0.004 0.097 0.356 FALSE 1 2.475 3
pit[14,3] 0.057 0.042 0.010 0.060 0.099 FALSE 1 8.243 2
pit[15,3] 0.429 0.192 0.175 0.404 0.708 FALSE 1 1.080 6
pit[16,3] 0.099 0.044 0.045 0.108 0.143 FALSE 1 4.856 2
pit[17,3] 0.262 0.206 0.052 0.238 0.492 FALSE 1 9.405 2
pit[18,3] 0.400 0.153 0.177 0.416 0.573 FALSE 1 1.203 6
pit[19,3] 0.314 0.221 0.043 0.320 0.569 FALSE 1 0.955 6
pit[20,3] 0.150 0.108 0.045 0.114 0.325 FALSE 1 1.776 4
pit[1,4] 0.280 0.299 0.014 0.191 0.635 FALSE 1 5.639 2
pit[2,4] 0.329 0.317 0.113 0.225 0.877 FALSE 1 1.243 6
pit[3,4] 0.472 0.204 0.208 0.462 0.750 FALSE 1 1.074 6
pit[4,4] 0.457 0.293 0.146 0.460 0.806 FALSE 1 6.732 2
pit[5,4] 0.268 0.148 0.036 0.301 0.406 FALSE 1 0.943 6
pit[6,4] 0.251 0.213 0.047 0.193 0.561 FALSE 1 0.892 6
pit[7,4] 0.104 0.106 0.022 0.069 0.287 FALSE 1 1.347 6
pit[8,4] 0.220 0.079 0.140 0.194 0.346 FALSE 1 1.804 4
pit[9,4] 0.379 0.306 0.127 0.213 0.830 FALSE 1 2.284 3
pit[10,4] 0.482 0.209 0.164 0.582 0.660 FALSE 1 0.896 6
pit[11,4] 0.052 0.045 0.012 0.043 0.127 FALSE 1 1.629 4
pit[12,4] 0.465 0.136 0.341 0.418 0.671 FALSE 1 3.070 2
pit[13,4] 0.496 0.216 0.260 0.462 0.758 FALSE 1 2.080 3
pit[14,4] 0.365 0.194 0.176 0.301 0.642 FALSE 1 0.847 6
pit[15,4] 0.371 0.277 0.126 0.281 0.785 FALSE 1 1.523 4
pit[16,4] 0.371 0.254 0.088 0.375 0.699 FALSE 1 1.163 6
pit[17,4] 0.442 0.172 0.209 0.427 0.686 FALSE 1 1.040 6
pit[18,4] 0.527 0.300 0.141 0.629 0.804 FALSE 1 0.818 6
pit[19,4] 0.563 0.144 0.384 0.580 0.735 FALSE 1 3.056 2
pit[20,4] 0.198 0.096 0.104 0.187 0.343 FALSE 1 3.528 2
pit[1,5] 0.264 0.128 0.157 0.229 0.482 FALSE 1 1.119 6
pit[2,5] 0.397 0.102 0.229 0.419 0.496 FALSE 1 1.460 5
pit[3,5] 0.396 0.202 0.100 0.493 0.580 FALSE 1 0.888 6
pit[4,5] 0.423 0.092 0.312 0.403 0.532 FALSE 1 3.020 2
pit[5,5] 0.268 0.208 0.033 0.219 0.523 FALSE 1 3.875 2
pit[6,5] 0.330 0.158 0.172 0.326 0.514 FALSE 1 0.871 6
pit[7,5] 0.329 0.133 0.142 0.336 0.483 FALSE 1 0.905 6
[ reached 'max' / getOption("max.print") -- omitted 14 rows ]
**WARNING** Rhat values indicate convergence failure.
Rhat is the potential scale reduction factor (at convergence, Rhat=1).
For each parameter, n.eff is a crude measure of effective sample size.
overlap0 checks if 0 falls in the parameter's 95% credible interval.
f is the proportion of the posterior with the same sign as the mean;
i.e., our confidence that the parameter is positive or negative.
DIC info: (pD = var(deviance)/2)
pD = 71.8 and DIC = 239.53
DIC is an estimate of expected predictive error (lower is better).
>
if i ignore the stochastic distribution of pi't it works actually this is the constant pit i mentioned. everybody is doing it on every tutorial
The code
nx <- "
model{
# Abundance
for (i in 1:nSites) {
N[i]~ dpois(lambda[i])
log(lambda[i])<-alphao+inprod(alpha, location[i,])
}
for (i in 1:nSites) {
for(t in 1:nOcc){
y[i,t]~ dbin(pit[i,t],N[i])
#pit_h[i,t]~ dbeta(1,1)
logit(pit[i,t]) <- beta0+inprod(beta, location[i,])+inprod(beta1,c(cov1[i,t],cov2[i,t]))
}
}
# Priors
alphao ~ dnorm(1.2824,0.302)
for(i in 1:nA){
alpha[i] ~ dnorm(0.284,0.570)
}
beta0 ~ dunif(-1.67,0.61)
for(i in 1:3){
beta[i] ~ dnorm(-0.370,0.254)
}
for(i in 1:2){
beta1[i] ~ dnorm(-0.104,0.44)
}
# det
sigma ~ dunif(0,1)
fac <-(1-sigma^2)/sigma^2
# derived
}"
writeLines(nx,con="mod1.txt")
watch=c("alphao","alpha","lambda","beta0","beta","beta1","pit","sigma")
inits = function() list(N = apply(y,1,max,na.rm=T))
set.seed(100)
mod<-jagsUI::jags(data,
parameters.to.save=watch,inits=inits,
model.file="mod1.txt",n.iter=3,
n.chains=2,DIC=TRUE)
mod
yielding
JAGS output for model 'mod1.txt', generated by jagsUI.
Estimates based on 2 chains of 3 iterations,
adaptation = 100 iterations (sufficient),
burn-in = 0 iterations and thin rate = 1,
yielding 6 total samples from the joint posterior.
MCMC ran for 0.006 minutes at time 2023-01-11 17:02:18.
mean sd 2.5% 50% 97.5% overlap0 f Rhat n.eff
alphao 4.405 0.098 4.247 4.448 4.489 FALSE 1.0 2.264 3
alpha[1] 0.705 0.311 0.386 0.702 1.052 FALSE 1.0 8.735 2
alpha[2] 1.568 0.229 1.346 1.550 1.829 FALSE 1.0 10.655 2
alpha[3] 1.300 0.086 1.167 1.333 1.380 FALSE 1.0 2.526 3
lambda[1] 82.166 7.700 69.933 85.482 89.069 FALSE 1.0 2.337 3
lambda[2] 169.553 39.141 128.321 171.241 207.988 FALSE 1.0 12.838 2
lambda[3] 396.653 61.786 337.056 388.711 467.843 FALSE 1.0 9.972 2
lambda[4] 302.649 41.450 260.019 302.933 343.100 FALSE 1.0 19.472 2
lambda[5] 82.166 7.700 69.933 85.482 89.069 FALSE 1.0 2.337 3
lambda[6] 169.553 39.141 128.321 171.241 207.988 FALSE 1.0 12.838 2
lambda[7] 396.653 61.786 337.056 388.711 467.843 FALSE 1.0 9.972 2
lambda[8] 302.649 41.450 260.019 302.933 343.100 FALSE 1.0 19.472 2
lambda[9] 82.166 7.700 69.933 85.482 89.069 FALSE 1.0 2.337 3
lambda[10] 169.553 39.141 128.321 171.241 207.988 FALSE 1.0 12.838 2
lambda[11] 396.653 61.786 337.056 388.711 467.843 FALSE 1.0 9.972 2
lambda[12] 302.649 41.450 260.019 302.933 343.100 FALSE 1.0 19.472 2
lambda[13] 82.166 7.700 69.933 85.482 89.069 FALSE 1.0 2.337 3
lambda[14] 169.553 39.141 128.321 171.241 207.988 FALSE 1.0 12.838 2
lambda[15] 396.653 61.786 337.056 388.711 467.843 FALSE 1.0 9.972 2
lambda[16] 302.649 41.450 260.019 302.933 343.100 FALSE 1.0 19.472 2
lambda[17] 82.166 7.700 69.933 85.482 89.069 FALSE 1.0 2.337 3
lambda[18] 169.553 39.141 128.321 171.241 207.988 FALSE 1.0 12.838 2
lambda[19] 396.653 61.786 337.056 388.711 467.843 FALSE 1.0 9.972 2
lambda[20] 302.649 41.450 260.019 302.933 343.100 FALSE 1.0 19.472 2
beta0 0.598 0.010 0.583 0.598 0.609 FALSE 1.0 2.973 2
beta[1] -0.669 0.138 -0.825 -0.663 -0.526 FALSE 1.0 10.339 2
beta[2] 0.071 0.148 -0.087 0.075 0.217 TRUE 0.5 14.514 2
beta[3] -0.260 0.079 -0.348 -0.261 -0.178 FALSE 1.0 9.102 2
beta1[1] -0.027 0.136 -0.154 -0.027 0.101 TRUE 0.5 73.857 2
beta1[2] -0.241 0.272 -0.498 -0.241 0.014 TRUE 0.5 73.668 2
pit[1,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 14.106 1
pit[2,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 34.146 1
pit[3,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 48.475 1
pit[4,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 19.522 1
pit[5,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 4.367 1
pit[6,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 29.277 1
pit[7,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 51.736 1
pit[8,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 17.004 1
pit[9,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 19.399 1
pit[10,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 45.000 1
pit[11,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 46.585 1
pit[12,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 5.782 1
pit[13,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 9.274 1
pit[14,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 17.662 1
pit[15,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 46.144 1
pit[16,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 3.164 1
pit[17,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 12.268 1
pit[18,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 24.030 1
pit[19,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.667 1
pit[20,1] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 31.711 1
pit[1,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 13.064 1
pit[2,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 16.990 1
pit[3,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 44.056 1
pit[4,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 17.895 1
pit[5,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 19.220 1
pit[6,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 46.416 1
pit[7,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.210 1
pit[8,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 8.495 1
pit[9,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 24.812 1
pit[10,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 31.174 1
pit[11,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 48.592 1
pit[12,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 28.506 1
pit[13,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 8.271 1
pit[14,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 39.914 1
pit[15,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.118 1
pit[16,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 1.402 1
pit[17,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 3.396 1
pit[18,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 26.605 1
pit[19,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 20.285 1
pit[20,2] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 20.645 1
pit[1,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 7.704 1
pit[2,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 0.944 1
pit[3,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.881 1
pit[4,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 15.553 1
pit[5,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 14.434 1
pit[6,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 2.484 1
pit[7,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 50.668 1
pit[8,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 24.700 1
pit[9,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 5.467 1
pit[10,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 58.360 1
pit[11,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 51.129 1
pit[12,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 31.265 1
pit[13,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 3.455 1
pit[14,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 9.752 1
pit[15,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.242 1
pit[16,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 20.779 1
pit[17,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 28.627 1
pit[18,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 8.427 1
pit[19,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 5.913 1
pit[20,3] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 15.146 1
pit[1,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 10.232 1
pit[2,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 38.961 1
pit[3,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 50.758 1
pit[4,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 32.788 1
pit[5,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 4.176 1
pit[6,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 50.633 1
pit[7,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 35.684 1
pit[8,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 7.917 1
pit[9,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 7.211 1
pit[10,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 39.496 1
pit[11,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 40.445 1
pit[12,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 19.869 1
pit[13,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 24.888 1
pit[14,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 4.988 1
pit[15,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 52.557 1
pit[16,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 4.046 1
pit[17,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 23.680 1
pit[18,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 30.557 1
pit[19,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 22.686 1
pit[20,4] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 35.248 1
pit[1,5] 0.000 0.000 0.000 0.000 0.000 FALSE 1.0 6.609 1
[ reached 'max' / getOption("max.print") -- omitted 21 rows ]
**WARNING** Rhat values indicate convergence failure.
Rhat is the potential scale reduction factor (at convergence, Rhat=1).
For each parameter, n.eff is a crude measure of effective sample size.
overlap0 checks if 0 falls in the parameter's 95% credible interval.
f is the proportion of the posterior with the same sign as the mean;
i.e., our confidence that the parameter is positive or negative.
DIC info: (pD = var(deviance)/2)
pD = 19.7 and DIC = 1249.054
DIC is an estimate of expected predictive error (lower is better).
But if i include both the covariate effect for detection and the stochastic distribution for detection probability. things go south. see the code below
nx <- "
model{
# Abundance
for (i in 1:nSites) {
N[i]~ dpois(lambda[i])
log(lambda[i])<-alphao+inprod(alpha, location[i,])
}
for (i in 1:nSites) {
for(t in 1:nOcc){
y[i,t]~ dbin(pit[i,t],N[i])
pit[i,t]~ dbeta(1,1)
logit(pit[i,t]) <- beta0+inprod(beta, location[i,])+inprod(beta1,c(cov1[i,t],cov2[i,t]))
}
}
# Priors
alphao ~ dnorm(1.2824,0.302)
for(i in 1:nA){
alpha[i] ~ dnorm(0.284,0.570)
}
beta0 ~ dunif(-1.67,0.61)
for(i in 1:3){
beta[i] ~ dnorm(-0.370,0.254)
}
for(i in 1:2){
beta1[i] ~ dnorm(-0.104,0.44)
}
# det
sigma ~ dunif(0,1)
fac <-(1-sigma^2)/sigma^2
# derived
}"
writeLines(nx,con="mod1.txt")
watch=c("alphao","alpha","lambda","beta0","beta","beta1","pit","sigma")
inits = function() list(N = apply(y,1,max,na.rm=T))
set.seed(100)
mod<-jagsUI::jags(data,
parameters.to.save=watch,inits=inits,
model.file="mod1.txt",n.iter=3,
n.chains=2,DIC=TRUE)
mod
this is the error.
Error in jags.model(file = model.file, data = data, inits = inits, n.chains = n.chains, :
RUNTIME ERROR:
Compilation error on line 12.
Attempt to redefine node pit[1,1]
I understand it is telling me that; pit[i,t]~ dbeta(1,1) is being overwritten by logit(pit[i,t])<-beta0+inprod(beta1, location[i,])+inprod(beta2,c(cov1[i,t],cov2[i,t])) but how exactly should this model be implemented. here the model is implemented without detection covariates. Its not what I am looking for.
Note: I edited my question to match the parameters on the equations I wrote
I was getting this error because, my detection model is not properly specified. Two processes here are trying to determine pit[i,t] the first is
pit[i,t]~ dbeta(1,1)
then
logit(pit[i,t])<-beta0+inprod(beta1, location[i,])+inprod(beta2,c(cov1[i,j],cov2[i,j]))
so by the time the second process tries to set value for pit[i,t] it finds that another process has already done so hence this error.
what I didn't remember while writing my code was, pi't is a stochastic value, which follows a beta distribution, I should be working to estimate the parameters for the beta binomial distribution generating the stochastic pi't but here I went for the stochastic value.
here is the correct implementation
nx <- "
model{
# Abundance
for (i in 1:nSites) {
N[i]~ dpois(lambda[i])
log(lambda[i])<-alphao+inprod(alpha, location[i,])
}
for (i in 1:nSites) {
for(t in 1:nOcc){
y[i,t]~ dbin(pit_h[i,t],N[i])
pit_h[i,t]~ dbeta((pit[i,t]*fac),(fac*(1-pit[i,t])))
logit(pit[i,t]) <- beta0+inprod(beta, location[i,])+inprod(beta1,c(cov1[i,t],cov2[i,t]))
}
}
# Priors
alphao ~ dnorm(1.2824,0.302)
for(i in 1:nA){
alpha[i] ~ dnorm(0.284,0.570)
}
beta0 ~ dunif(-1.67,0.61)
for(i in 1:3){
beta[i] ~ dnorm(-0.370,0.254)
}
for(i in 1:2){
beta1[i] ~ dnorm(-0.104,0.44)
}
# det
sigma ~ dunif(0,1)
fac <-(1-sigma^2)/sigma^2
# derived
}"
writeLines(nx,con="mod1.txt")
watch=c("alphao","alpha","lambda","beta0","beta","beta1","pit","sigma")
inits = function() list(N = apply(y,1,max,na.rm=T))
set.seed(100)
mod<-jagsUI::jags(data,
parameters.to.save=watch,inits=inits,
model.file="mod1.txt",n.iter=3,
n.chains=2,DIC=TRUE)
mod
yielding
JAGS output for model 'mod1.txt', generated by jagsUI.
Estimates based on 2 chains of 3 iterations,
adaptation = 100 iterations (sufficient),
burn-in = 0 iterations and thin rate = 1,
yielding 6 total samples from the joint posterior.
MCMC ran for 0.051 minutes at time 2023-01-11 17:53:55.
mean sd 2.5% 50% 97.5% overlap0 f Rhat n.eff
alphao 2.114 0.887 1.205 2.127 2.970 FALSE 1.0 19.132 2
alpha[1] -0.181 0.774 -0.929 -0.232 0.674 TRUE 0.5 13.424 2
alpha[2] -0.792 0.437 -1.211 -0.821 -0.284 FALSE 1.0 10.708 2
alpha[3] -0.123 0.391 -0.537 -0.114 0.306 TRUE 0.5 10.409 2
lambda[1] 11.143 8.192 3.338 10.511 19.499 FALSE 1.0 16.476 2
lambda[2] 6.989 1.152 5.440 7.218 8.436 FALSE 1.0 2.041 3
lambda[3] 4.078 1.757 2.388 3.947 5.864 FALSE 1.0 18.172 2
lambda[4] 12.868 11.770 2.060 11.100 26.299 FALSE 1.0 10.437 2
lambda[5] 11.143 8.192 3.338 10.511 19.499 FALSE 1.0 16.476 2
lambda[6] 6.989 1.152 5.440 7.218 8.436 FALSE 1.0 2.041 3
lambda[7] 4.078 1.757 2.388 3.947 5.864 FALSE 1.0 18.172 2
lambda[8] 12.868 11.770 2.060 11.100 26.299 FALSE 1.0 10.437 2
lambda[9] 11.143 8.192 3.338 10.511 19.499 FALSE 1.0 16.476 2
lambda[10] 6.989 1.152 5.440 7.218 8.436 FALSE 1.0 2.041 3
lambda[11] 4.078 1.757 2.388 3.947 5.864 FALSE 1.0 18.172 2
lambda[12] 12.868 11.770 2.060 11.100 26.299 FALSE 1.0 10.437 2
lambda[13] 11.143 8.192 3.338 10.511 19.499 FALSE 1.0 16.476 2
lambda[14] 6.989 1.152 5.440 7.218 8.436 FALSE 1.0 2.041 3
lambda[15] 4.078 1.757 2.388 3.947 5.864 FALSE 1.0 18.172 2
lambda[16] 12.868 11.770 2.060 11.100 26.299 FALSE 1.0 10.437 2
lambda[17] 11.143 8.192 3.338 10.511 19.499 FALSE 1.0 16.476 2
lambda[18] 6.989 1.152 5.440 7.218 8.436 FALSE 1.0 2.041 3
lambda[19] 4.078 1.757 2.388 3.947 5.864 FALSE 1.0 18.172 2
lambda[20] 12.868 11.770 2.060 11.100 26.299 FALSE 1.0 10.437 2
beta0 -0.162 0.447 -0.600 -0.186 0.325 TRUE 0.5 14.944 2
beta[1] 0.924 0.693 0.243 0.892 1.644 FALSE 1.0 14.009 2
beta[2] 2.147 0.133 1.987 2.152 2.288 FALSE 1.0 7.489 2
beta[3] 1.047 0.911 0.072 1.073 1.920 FALSE 1.0 17.258 2
beta1[1] -0.088 0.008 -0.096 -0.088 -0.080 FALSE 1.0 42.539 2
beta1[2] 0.121 0.022 0.100 0.121 0.142 FALSE 1.0 55.542 2
pit[1,1] 0.067 0.041 0.029 0.061 0.114 FALSE 1.0 10.366 2
pit[2,1] 0.117 0.015 0.098 0.122 0.131 FALSE 1.0 0.923 6
pit[3,1] 0.341 0.181 0.169 0.324 0.538 FALSE 1.0 13.526 2
pit[4,1] 0.195 0.183 0.025 0.187 0.379 FALSE 1.0 29.881 2
pit[5,1] 0.060 0.037 0.026 0.054 0.102 FALSE 1.0 10.199 2
pit[6,1] 0.114 0.014 0.095 0.118 0.127 FALSE 1.0 0.926 6
pit[7,1] 0.360 0.186 0.184 0.345 0.563 FALSE 1.0 13.999 2
pit[8,1] 0.200 0.187 0.026 0.191 0.387 FALSE 1.0 30.322 2
pit[9,1] 0.071 0.043 0.031 0.064 0.120 FALSE 1.0 10.459 2
pit[10,1] 0.128 0.016 0.107 0.133 0.142 FALSE 1.0 0.915 6
pit[11,1] 0.334 0.179 0.164 0.318 0.530 FALSE 1.0 13.373 2
pit[12,1] 0.216 0.201 0.029 0.207 0.415 FALSE 1.0 31.999 2
pit[13,1] 0.063 0.039 0.028 0.058 0.108 FALSE 1.0 10.284 2
pit[14,1] 0.105 0.013 0.088 0.110 0.118 FALSE 1.0 0.933 6
pit[15,1] 0.333 0.178 0.163 0.316 0.529 FALSE 1.0 13.343 2
pit[16,1] 0.219 0.203 0.030 0.211 0.421 FALSE 1.0 32.359 2
pit[17,1] 0.065 0.040 0.029 0.060 0.112 FALSE 1.0 10.335 2
pit[18,1] 0.110 0.014 0.091 0.114 0.123 FALSE 1.0 0.929 6
pit[19,1] 0.374 0.189 0.194 0.358 0.578 FALSE 1.0 14.324 2
pit[20,1] 0.161 0.153 0.019 0.153 0.316 FALSE 1.0 26.658 2
pit[1,2] 0.066 0.040 0.029 0.060 0.113 FALSE 1.0 10.349 2
pit[2,2] 0.105 0.013 0.087 0.109 0.118 FALSE 1.0 0.934 6
pit[3,2] 0.327 0.177 0.160 0.311 0.522 FALSE 1.0 13.219 2
pit[4,2] 0.198 0.186 0.026 0.190 0.384 FALSE 1.0 30.170 2
pit[5,2] 0.070 0.043 0.031 0.064 0.120 FALSE 1.0 10.455 2
pit[6,2] 0.130 0.016 0.109 0.135 0.144 FALSE 1.0 0.914 6
pit[7,2] 0.366 0.187 0.188 0.350 0.569 FALSE 1.0 14.133 2
pit[8,2] 0.212 0.198 0.028 0.204 0.409 FALSE 1.0 31.621 2
pit[9,2] 0.075 0.046 0.034 0.069 0.127 FALSE 1.0 10.557 2
pit[10,2] 0.115 0.014 0.096 0.120 0.128 FALSE 1.0 0.925 6
pit[11,2] 0.341 0.181 0.169 0.325 0.539 FALSE 1.0 13.537 2
pit[12,2] 0.174 0.165 0.021 0.165 0.340 FALSE 1.0 27.825 2
pit[13,2] 0.062 0.038 0.027 0.057 0.107 FALSE 1.0 10.267 2
pit[14,2] 0.122 0.015 0.102 0.128 0.136 FALSE 1.0 0.919 6
pit[15,2] 0.365 0.187 0.187 0.349 0.568 FALSE 1.0 14.103 2
pit[16,2] 0.225 0.208 0.031 0.216 0.431 FALSE 1.0 32.922 2
pit[17,2] 0.059 0.036 0.026 0.054 0.101 FALSE 1.0 10.182 2
pit[18,2] 0.112 0.014 0.093 0.116 0.125 FALSE 1.0 0.928 6
pit[19,2] 0.297 0.167 0.139 0.280 0.483 FALSE 1.0 12.542 2
pit[20,2] 0.193 0.182 0.025 0.185 0.375 FALSE 1.0 29.672 2
pit[1,3] 0.062 0.038 0.027 0.056 0.106 FALSE 1.0 10.257 2
pit[2,3] 0.094 0.012 0.077 0.098 0.105 FALSE 1.0 0.944 6
pit[3,3] 0.409 0.195 0.223 0.395 0.619 FALSE 1.0 15.215 2
pit[4,3] 0.202 0.189 0.026 0.194 0.391 FALSE 1.0 30.561 2
pit[5,3] 0.067 0.041 0.030 0.061 0.114 FALSE 1.0 10.372 2
pit[6,3] 0.095 0.012 0.079 0.099 0.107 FALSE 1.0 0.943 6
pit[7,3] 0.352 0.184 0.177 0.336 0.552 FALSE 1.0 13.790 2
pit[8,3] 0.185 0.174 0.023 0.176 0.359 FALSE 1.0 28.825 2
pit[9,3] 0.060 0.037 0.026 0.055 0.104 FALSE 1.0 10.218 2
pit[10,3] 0.165 0.019 0.140 0.172 0.182 FALSE 1.0 0.893 6
pit[11,3] 0.355 0.184 0.180 0.339 0.556 FALSE 1.0 13.870 2
pit[12,3] 0.267 0.241 0.042 0.260 0.504 FALSE 1.0 37.848 2
pit[13,3] 0.059 0.036 0.026 0.054 0.101 FALSE 1.0 10.183 2
pit[14,3] 0.086 0.011 0.071 0.089 0.097 FALSE 1.0 0.953 6
pit[15,3] 0.430 0.198 0.241 0.416 0.641 FALSE 1.0 15.755 2
pit[16,3] 0.193 0.181 0.025 0.184 0.374 FALSE 1.0 29.647 2
pit[17,3] 0.078 0.047 0.035 0.072 0.133 FALSE 1.0 10.632 2
pit[18,3] 0.087 0.012 0.072 0.090 0.098 FALSE 1.0 0.952 6
pit[19,3] 0.286 0.163 0.131 0.269 0.467 FALSE 1.0 12.293 2
pit[20,3] 0.242 0.221 0.035 0.233 0.460 FALSE 1.0 34.824 2
pit[1,4] 0.064 0.039 0.028 0.058 0.109 FALSE 1.0 10.300 2
pit[2,4] 0.122 0.015 0.102 0.127 0.135 FALSE 1.0 0.920 6
pit[3,4] 0.352 0.184 0.178 0.336 0.553 FALSE 1.0 13.805 2
pit[4,4] 0.155 0.148 0.018 0.147 0.305 FALSE 1.0 26.121 2
pit[5,4] 0.059 0.037 0.026 0.054 0.102 FALSE 1.0 10.196 2
pit[6,4] 0.136 0.016 0.114 0.142 0.151 FALSE 1.0 0.910 6
pit[7,4] 0.313 0.173 0.150 0.296 0.504 FALSE 1.0 12.894 2
pit[8,4] 0.213 0.198 0.029 0.205 0.410 FALSE 1.0 31.703 2
pit[9,4] 0.062 0.038 0.027 0.056 0.106 FALSE 1.0 10.249 2
pit[10,4] 0.122 0.015 0.102 0.127 0.136 FALSE 1.0 0.920 6
pit[11,4] 0.320 0.175 0.155 0.304 0.513 FALSE 1.0 13.056 2
pit[12,4] 0.248 0.226 0.037 0.240 0.471 FALSE 1.0 35.556 2
pit[13,4] 0.075 0.046 0.034 0.069 0.128 FALSE 1.0 10.559 2
pit[14,4] 0.090 0.012 0.074 0.093 0.101 FALSE 1.0 0.948 6
pit[15,4] 0.421 0.197 0.233 0.407 0.632 FALSE 1.0 15.524 2
pit[16,4] 0.218 0.203 0.030 0.210 0.419 FALSE 1.0 32.238 2
pit[17,4] 0.074 0.045 0.033 0.068 0.126 FALSE 1.0 10.536 2
pit[18,4] 0.115 0.014 0.095 0.119 0.128 FALSE 1.0 0.925 6
pit[19,4] 0.299 0.168 0.140 0.283 0.486 FALSE 1.0 12.588 2
pit[20,4] 0.277 0.249 0.044 0.270 0.520 FALSE 1.0 39.039 2
pit[1,5] 0.052 0.032 0.022 0.047 0.089 FALSE 1.0 10.005 2
[ reached 'max' / getOption("max.print") -- omitted 21 rows ]
**WARNING** Rhat values indicate convergence failure.
Rhat is the potential scale reduction factor (at convergence, Rhat=1).
For each parameter, n.eff is a crude measure of effective sample size.
overlap0 checks if 0 falls in the parameter's 95% credible interval.
f is the proportion of the posterior with the same sign as the mean;
i.e., our confidence that the parameter is positive or negative.
DIC info: (pD = var(deviance)/2)
pD = 32.5 and DIC = 215.949
DIC is an estimate of expected predictive error (lower is better).
>