leire
2014-04-01
Hi all,
I’d like to translate some fish population dynamics models from WinBUGS to JAGS. But I’ve a bounding issue and I’m not sure I’m dealing with it properly in JAGS.
Basically the model derives fish biomass (tonnes of fish) per year (deterministically) depending on the model parameters. Fish catches are assumed to be observed without error and are removed from the fish biomass. Therefore the biomasses need to be larger than the catches to avoid negative amount of fish. We used to implement this in WinBUGS using the zeros trick. Now in JAGS I have the following:
data{ for (y in 1:Y){ ones1[y] <- 1 } } model{ # rest of the code is omitted for (y in 1:Y){ logR[y] ~ dnorm(mur, psir) R[y] <- exp(logR[y]) B1[y] <- R[y]*exp(-g*f) - C1.per1[y]*exp(-g*(f-h1[y])) prob1[y] <- step(B1[y]) ones1[y] ~ dbern(prob1[y]) } }
I’ve tried to implement it in JAGS using dinterval() as follows:
data{ for (y in 1:Y){ ones1[y] <- 1 } } model{ # rest of the code is omitted for (y in 1:Y){ logR[y] ~ dnorm(mur, psir) R[y] <- exp(logR[y]) B1[y] <- R[y]*exp(-g*f) - C1.per1[y]*exp(-g*(f-h1[y])) ones1[y] ~ dinterval(B1[y], 0) } }
And I get quite similar results.
I’ve seen that in some implementations, they use max() to avoid negative values. Something like:
model{ for (y in 1:Y){ logR[y] ~ dnorm(mur, psir) R[y] <- exp(logR[y]) B1[y] <- max(R[y]*exp(-g*f) - C1.per1[y]*exp(-g*(f-h1[y])), 0.001) } }
But I’m not sure this is correct.
I’d like to know whether these approaches are valid or not. And if there is any preference for any of them or a better way of implementing this.
Thank you very much in advance,
Leire