I am trying to translate an occupancy model from jags to julia and arrived at the following code:
@model function model1(y, n_cells, n_years, n_traps, n_period, z, z_su, ::Type{T} = Float64) where {T}
#
# -- define / preallocate
psi = Array{T}(undef, (n_years, n_cells))
mu_su = Array{T}(undef, (n_years, n_cells, n_traps))
theta = Array{T}(undef, (n_years, n_cells, n_traps))
mu_y = Array{T}(undef, (n_cells, n_years, n_traps, n_period))
p = Array{T}(undef, (n_cells, n_years, n_traps, n_period))
#
# -- priors
psi0 ~ filldist(Beta(1, 1), n_years)
theta0 ~ filldist(Beta(1, 1), n_years)
p0 ~ filldist(Beta(1, 1), n_years)
#
int_psi = logit.(psi0)
int_theta = logit.(theta0)
int_p = logit.(p0)
#
# -- Likelihood (for model structure)
for yy in 1:n_years
for i in 1:n_cells
# Occupancy in grid cell i for year yy
psi[yy, i] = logistic(int_psi[yy])
z[yy, i] ~ Bernoulli(psi[yy, i])
#
for j in 1:n_traps
# Site-use at camera location j for year yy
theta[yy, i, j] = logistic(int_theta[yy])
mu_su[yy, i, j] = z[yy, i] * theta[yy, i, j]
z_su[yy, i, j] ~ Bernoulli(mu_su[yy, i, j])
#
for k in 1:n_period
# detection probability at camera location j in period k and year yy
p[i, yy, j, k] = logistic(int_p[yy])
mu_y[i, yy, j, k] = z_su[yy, i, j] * p[i, yy, j, k]
y[i, yy, j, k] ~ Bernoulli(mu_y[i, yy, j, k])
#
end # End period
end # End camera trap
end # End grid cell
end # End year
end # End model
The MH sampler seems to work with dummy data (although I don’t get any standard errors), but I can’t get the NUTS sampler to work… It spits a lot of warnings about rejected proposals, but the
progress bar doesn’t change at all. So I am obviously doing some noob mistakes here…
I would really appreciate if a more experienced turing user would take a (quick) look a this code and point me to what is blatantly wrong in there… Many thanks!