So I’m making progress with implementing a bunch of models in Turing. But one I’ve had problems with is this.

```
using Turing, StatsPlots
sᵇ = [15,11,15,14,15,18,16,16,18,16,15,13,18,12,11,13,17,18,16,11,17,18,12,18,18,14,21,18,17,10,11,12,16,18,17,15,19,12,21,15,16,20,15,19,16,16,14,18,16,19,17,11,19,18,16,16,11,19,18,12,15,18,20, 8,12,19,16,16,16,12,18,17,11,20]
nᵇ = 21
sⁿ = [15,12,14,15,13,14,10,17,13,16,16,10,15,15,10,14,17,18,19,12,19,18,10,18,16,13,15,20,13,15,13,14,19,19,19,18,13,12,19,16,14,17,15,16,15,16,13,15,14,19,12,11,17,13,18,13,13,19,18,13,13,16,18,14,14,17,12,12,16,14,16,18,13,13]
nⁿ = 21
@assert length(sᵇ) == length(sⁿ)
Φ(x) = cdf(Normal(0, 1), x)
@model function model(sᵇ, nᵇ, sⁿ, nⁿ)
N = length(sᵇ)
μ ~ truncated(Normal(0, 1), 0, Inf)
σ ~ Uniform(0, 10)
δ ~ truncated(Normal(0, 1), 0, Inf)
σα ~ Uniform(0, 10)
μα = δ * σα
α = Vector(undef, N)
ϕⁿ = Vector(undef, N)
ϕᵇ = Vector(undef, N)
θⁿ = Vector(undef, N)
θᵇ = Vector(undef, N)
for i in eachindex(sᵇ)
α[i] ~ Normal(μα, σα)
ϕⁿ[i] ~ Normal(μ, σ)
ϕᵇ[i] = ϕⁿ[i] + α[i]
θⁿ[i] = Φ(ϕⁿ[i])
θᵇ[i] = Φ(ϕᵇ[i])
# likelihood
sⁿ[i] ~ Binomial(nⁿ, θⁿ[i])
sᵇ[i] ~ Binomial(nᵇ, θᵇ[i])
end
end
# sample from prior
prior_chain = sample(model(sᵇ, nᵇ, sⁿ, nⁿ), Prior(), 5000)
# sample from posterior
posterior_chain = sample(model(sᵇ, nᵇ, sⁿ, nⁿ), HMC(0.01, 10), 5000)
```

Inferences about the parameter of interest, `δ`

seem roughly correct compared to what I’m expecting.

The main problem is that sampling is extremely slow and sometimes seemingly stalls. I’m not sure if it’s because defining a function as `Φ(x) = cdf(Normal(0, 1), x)`

is a sub-optimal way of doing it, or because I’m not being specific enough about types of the deterministic variables. Any insights would be appreciated.