Fantastic. Thank you François.

Here is a MWE (well not that minimal I tried my best). Second model is almost always rejected.

```
using Turing
using DifferentialEquations
time=0.0:4.0
exposure=fill(0.0,5)
survival=[20,19,19,19,18]
loglogistic1(x, α, β) = x^β /(α^β + x^β)
loglogistic2(x,α, β) = 1/(1+(x/α)^(-β))
Turing.@model MWE_ll1(time, exposure, survival) = begin
# 1. SET PRIORS
kd_log10 ~ Normal(-1.38, 1.11)
hb_log10 ~ Normal(-1.88,0.86)
α_log10 ~ Normal(1.23, 0.16)
β_log10 ~ Uniform(-2.0, 2.0)
kd = 10 ^kd_log10
hb = 10 ^hb_log10
α = 10 ^α_log10
β = 10 ^β_log10
# 2. MECHANISTIC MODELS
prob = ODEProblem(odeTK, eltype(kd).([0.0]),(minimum(time), maximum(time)),[time, exposure, kd])
_saveat = time === nothing ? Float64[] : time
sol_tmp = solve(prob; saveat = _saveat)
Dmax_tmp = accumulate(max, sol_tmp.u)
N = length(sol_tmp)
pSurv = Vector{Any}(undef, N)
pSurv[1] = 1.0
for i in 2:N
pSurv[i] = exp(-hb*time[i])*(1-loglogistic1(Dmax_tmp[i][1], α,β))
end
# 3. INFERENCE
ratioSurv = Vector{Real}(undef, N)
ratioSurv[1] = 1.0
for i in 2:N
ratioSurv[i] = pSurv[i] / pSurv[i-1]
survival[i] ~ Binomial(survival[i-1], ratioSurv[i])
end
end
model_ll1 = MWE_ll1(time, exposure, survival)
chn_ll1 = sample(model_ll1, HMC(0.01, 10), 10)
Turing.@model MWE_ll2(time, exposure, survival) = begin
# 1. SET PRIORS
kd_log10 ~ Normal(-1.38, 1.11)
hb_log10 ~ Normal(-1.88,0.86)
α_log10 ~ Normal(1.23, 0.16)
β_log10 ~ Uniform(-2.0, 2.0)
kd = 10 ^kd_log10
hb = 10 ^hb_log10
α = 10 ^α_log10
β = 10 ^β_log10
# 2. MECHANISTIC MODELS
prob = ODEProblem(odeTK, eltype(kd).([0.0]),(minimum(time), maximum(time)),[time, exposure, kd])
_saveat = time === nothing ? Float64[] : time
sol_tmp = solve(prob; saveat = _saveat)
Dmax_tmp = accumulate(max, sol_tmp.u)
N = length(sol_tmp)
pSurv = Vector{Any}(undef, N)
pSurv[1] = 1.0
for i in 2:N
pSurv[i] = exp(-hb*time[i])*(1-loglogistic2(Dmax_tmp[i][1], α,β))
end
# 3. INFERENCE
ratioSurv = Vector{Real}(undef, N)
ratioSurv[1] = 1.0
for i in 2:N
ratioSurv[i] = pSurv[i] / pSurv[i-1]
survival[i] ~ Binomial(survival[i-1], ratioSurv[i])
end
end
model_ll2 = MWE_ll2(time, exposure, survival)
chn_ll2 = sample(model_ll2, HMC(0.01, 10), 10)
```