# Trouble translating a simple PyMC3 model to Turing

**URL:** <https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856>\
**Category:** New to Julia\
**Tags:** turing\
**Created:** [November 21, 2021, 10:27am UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856 "2021-11-21T10:27:12Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 10:27am UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/1 "2021-11-21T10:27:12Z")

</div>

I’ve been playing around with Turing solving this simple problem from Allen Downey’s _Think Bayes_:

23 grizzly bears visit a “trap” in the first season. In the 2nd season, there are 19 bears, 4 of which are repeat customers. We need to estimate the total number of bears.

This is the solution Downey gives using [PyMC3](https://github.com/AllenDowney/ThinkBayes2/blob/master/soln/chap19.ipynb):

```python
k = 23
n = 19
x = 4

with pm.Model() as model6:
    N = pm.DiscreteUniform('N', 50, 500)
    y = pm.HyperGeometric('y', N=N, k=k, n=n, observed=x)
    trace6 = pm.sample(1000, **options)

```

The mean of the posterior distribution for N is 182.

Now, this _should_ be the equivalent Julia code:

```julia
k = 23
n = 19
x = 4

@model function model6(y)
    N ~ DiscreteUniform(50, 500)
    y ~ Hypergeometric(k, N - k, n)
end
trace6 = sample(model6(x), IS(), 500);

```

(is IS the right sampler to use? NUTS crashes with an uninformative “Inexact” error)

the mean of this distribution is 284.

Where am I going wrong?

Thanks

DD

---

<div class="post-metadata">

**Author:** ![cjdoris](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cjdoris/32/213133_2.png) [@cjdoris](https://discourse.julialang.org/u/cjdoris)\
**Post date:** [November 21, 2021, 10:35am UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/2 "2021-11-21T10:35:39Z")

</div>

You sample 1000 in the Python case and 500 in the Julia case. Did you just divide by the wrong thing?

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 10:46am UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/3 "2021-11-21T10:46:23Z")

</div>

The code here is with a single chain instead of 2.  
The result doesn’t change if I increase the samples to 1000 or even 5000.

Neither switching to a different sampler - e.g. `SMC`.  
This is the result I get (with SMC & 1000 samples, but as I said, this hardly changes the result):

 ![Screen Shot 2021-11-21 at 12.43.22](https://global.discourse-cdn.com/julialang/original/3X/9/2/92f8139417e960e054af89bf385544c3926f2c33.png)

and this is Downey’s result:

 ![Screen Shot 2021-11-21 at 12.44.06](https://global.discourse-cdn.com/julialang/original/3X/f/d/fd98954f392b444d08e2c86abb431d24b50e013b.png)

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [November 21, 2021, 12:19pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/4 "2021-11-21T12:19:32Z")

</div>

Hi,

NUTS does not work on discrete parameters. Your full error message may have included `Int(some_value),` but I agree it is not the most informative.

I’m not sure what sampler would be best for this problem. What is used in PyMC3?

One unrelated point: you should pass variables to your model function rather than declare global variables. This is good practice in general, but specifically in Julia, global variables often reduce performance because the compiler cannot guarantee their types at the time of execution.

Here is a better approach:

```julia
using Turing, Random
Random.seed!(505)
k = 23
n = 19
y = 4

@model function model6(y, k, n)
    N ~ DiscreteUniform(50, 500)
    y ~ Hypergeometric(k, N - k, n)
end
trace6 = sample(model6(y, k, n), IS(), 500);

```

I can reproduce your results. Indeed, the posterior distribution is not changing from the prior.

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [November 21, 2021, 12:27pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/5 "2021-11-21T12:27:44Z")

</div>

`PG` seems to give reasonable results:

```julia
using Turing, Random
Random.seed!(505)
k = 23
n = 19
y = 4

@model function model6(y, k, n)
    N ~ DiscreteUniform(50, 500)
    y ~ Hypergeometric(k, N - k, n)
end
trace6 = sample(model6(y, k, n), PG(50), 1000)

```

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 12:28pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/6 "2021-11-21T12:28:19Z")

</div>

> [@Christopher\_Fisher](#):
>
> What is used in PyMC3?

Metropolis:

```julia
Multiprocess sampling (2 chains in 2 jobs)
Metropolis: [N]
100.00% [4000/4000 00:00<00:00 Sampling 2 chains, 0 divergences]
Sampling 2 chains for 1_000 tune and 1_000 draw iterations (2_000 + 2_000 draws total) took 1 seconds.
The number of effective samples is smaller than 25% for some parameters.

```

It would be nice if Turing chose sensible defaults automatically.

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [November 21, 2021, 12:31pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/7 "2021-11-21T12:31:25Z")

</div>

Does PYMC3 always default to Metropolis-Hastings or does it select a default based on your model properties?

Here is the Metropolis-Hastings version:

```julia
using Turing, Random
Random.seed!(505)
k = 23
n = 19
y = 4

@model function model6(y, k, n)
    N ~ DiscreteUniform(50, 500)
    y ~ Hypergeometric(k, N - k, n)
end
trace6 = sample(model6(y, k, n), MH(), 1000)

plot(trace6)

```

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 12:31pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/8 "2021-11-21T12:31:36Z")

</div>

> [@Christopher\_Fisher](#):
>
> `PG` seems to give reasonable results:

Thanks!  
But why aren’t `IS` and `SMC` working?

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [November 21, 2021, 12:32pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/9 "2021-11-21T12:32:46Z")

</div>

> [@dodoplus](#):
>
> Thanks!  
> But why aren’t `IS` and `SMC` working?

Good question. I am not exactly sure. It could be that they are not the best samplers for the problem or the hyper parameters need to be tweaked. Maybe one of the Turing devs will know.

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 12:35pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/10 "2021-11-21T12:35:30Z")

</div>

> [@Christopher\_Fisher](#):
>
> Does PYMC3 always default to Metropolis-Hastings or does it select a default based on your model properties?

It’s usually NUTS, I think.

So, apparently, it does sensible decisions regarding what algorithms to use for a specific problem.  
Why can’t we have a similar functionality in Turing? This is mainly a function of the distributions used in the model, right? Turing can definitely deduce whether e.g. a discrete random variable is required

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [November 21, 2021, 12:45pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/11 "2021-11-21T12:45:46Z")

</div>

Yeah. Thats a good question. I can’t speak for the Turing devs, but I suspect setting default samplers may not have a large benefit and may have the drawback of making the sampler more opaque. NUTS is reasonable continuous variables, but discrete parameters can be challenging and more problem-specific. In my experience, Metropolis-Hastings has poor scalability and can fail in many cases. I suspect it is a default that is overridden in many cases. However, most samplers have reasonable default settings when possible. For example, NUTS uses the same defaults in Stan. Sometimes you can change the target acceptance rate to improve performance. Even with NUTS, the most challenging part can be finding the right transformations to enable efficient sampling. It would be nice if Turing selected an optimal transformation, but I suspect that would be difficult to do 🙂

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 1:12pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/12 "2021-11-21T13:12:11Z")

</div>

Continuing this thread…

A similar issue occurs with this PyMC3 model (also from [Think Bayes](https://github.com/AllenDowney/ThinkBayes2/blob/master/soln/chap19.ipynb)):

```python
k10 = 20 - 3
k01 = 15 - 3
k11 = 3
num_seen = k01 + k10 + k11

with pm.Model() as model9:
    p0 = pm.Beta('p0', alpha=1, beta=1)
    p1 = pm.Beta('p1', alpha=1, beta=1)
    N = pm.DiscreteUniform('N', num_seen, 350)
    
    q0 = 1-p0
    q1 = 1-p1
    ps = [q0*q1, q0*p1, p0*q1, p0*p1]
    
    k00 = N - num_seen
    data = pm.math.stack((k00, k01, k10, k11))
    y = pm.Multinomial('y', n=N, p=ps, observed=data)
    trace9 = pm.sample(1000, **options)

```

PyMC3 chose Metropolis:

```julia
Multiprocess sampling (2 chains in 2 jobs)
CompoundStep
>NUTS: [p1, p0]
>Metropolis: [N]
100.00% [4000/4000 00:02<00:00 Sampling 2 chains, 0 divergences]
Sampling 2 chains for 1_000 tune and 1_000 draw iterations (2_000 + 2_000 draws total) took 3 seconds.
The acceptance probability does not match the target. It is 0.5430480274605854, but should be close to 0.8. Try to increase the number of tuning steps.
The rhat statistic is larger than 1.05 for some parameters. This indicates slight problems during sampling.
The estimated number of effective samples is smaller than 200 for some parameters.

```

Now, to Julia:

```julia
@model function model9(k01, k10, k11)
    p0 ~ Beta(1, 1)
    p1 ~ Beta(1, 1)
    num_seen = k01 + k10 + k11
    N ~ DiscreteUniform(num_seen, 350)
    
    q0 = 1-p0
    q1 = 1-p1
    ps = [q0*q1, q0*p1, p0*q1, p0*p1]
    
    y = [N - num_seen, k01, k10, k11]
    y ~ Multinomial(N, ps)
end
    
trace9 = sample(model9(k01, k10, k11), PG(50), 500)

```

We get:

```julia
Sampling: 100%|█████████████████████████████████████████| Time: 0:00:02
[112]:
Chains MCMC chain (500×9×1 Array{Float64, 3}):

Log evidence = 0.0
Iterations = 1:1:500
Number of chains = 1
Samples per chain = 500
Wall duration = 2.74 seconds
Compute duration = 2.74 seconds
parameters = p0, p1, N, y[1], y[2], y[3], y[4]
internals = lp, logevidence

Summary Statistics
  parameters mean std naive_se mcse ess rhat ⋯
      Symbol Float64 Float64 Float64 Float64 Float64 Float64 ⋯

          p0 0.5122 0.2965 0.0133 0.0125 449.8449 0.9981 ⋯
          p1 0.5232 0.2838 0.0127 0.0105 449.0627 0.9988 ⋯
           N 190.5040 90.7350 4.0578 4.4045 514.1959 1.0049 ⋯
        y[1] 44.0380 48.9169 2.1876 2.1886 572.9458 1.0044 ⋯
        y[2] 47.2080 48.1557 2.1536 1.8543 469.3461 0.9985 ⋯
        y[3] 47.9620 54.4876 2.4368 2.3048 503.4060 0.9983 ⋯
        y[4] 51.2960 54.6575 2.4444 2.4123 430.4700 1.0021 ⋯
                                                                1 column omitted

Quantiles
  parameters 2.5% 25.0% 50.0% 75.0% 97.5% 
      Symbol Float64 Float64 Float64 Float64 Float64 

          p0 0.0350 0.2440 0.5235 0.7708 0.9749
          p1 0.0447 0.2911 0.5222 0.7576 0.9844
           N 39.9500 108.0000 191.5000 267.0000 338.5250
        y[1] 0.0000 9.0000 25.0000 66.2500 180.7250
        y[2] 0.4750 10.0000 30.0000 68.0000 183.0500
        y[3] 0.0000 8.7500 27.0000 68.0000 199.2000
        y[4] 0.0000 11.0000 32.0000 72.0000 185.6750

```

and:

 ![Screen Shot 2021-11-21 at 15.08.14](https://global.discourse-cdn.com/julialang/original/3X/7/e/7e2823a77e6b9ab1346552682cfcd89235ef3277.jpeg)

This seems to be the same problem, but this time `PG(50)` doesn’t work…  
Neither does MH (which for some reason returns a normal distribution:

```julia
trace9 = sample(model9(k01, k10, k11), MH(), 500)

```

 ![Screen Shot 2021-11-21 at 15.11.56](https://global.discourse-cdn.com/julialang/original/3X/1/6/168b9b1c10ba7e0dd40966daa67570d9ca064d3d.png)

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [November 21, 2021, 4:04pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/13 "2021-11-21T16:04:31Z")

</div>

It doesn’t return a normal it returns a chain that didn’t move from its initial place and the density plot smoothed the delta function into a normal.

Notice how the normal is mean 146 and has no real mass by 145 or 147.

Sampling discrete distributions can be difficult but I’ll try to take a look at this later today to see if I can come up with some suggestions.

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [November 21, 2021, 4:11pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/14 "2021-11-21T16:11:41Z")

</div>

> [@dodoplus](#):
>
> ```julia
> NUTS: [p1, p0]
> >Metropolis: [N]
> 
> ```

PyMC chose mixed sampling, with p1 and p2 updated by NUTS and N updated by Metropolis (note Metropolis is not just one algorithm the question is what is the proposal it uses?)

You can do this kind of mixed sampling in Turing, I’m on my phone so hard to pull an exact link but check out the website docs and look for Compositional Gibbs sampling [https://turing.ml/dev/docs/using-turing/guide](https://turing.ml/dev/docs/using-turing/guide)

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [November 21, 2021, 5:05pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/15 "2021-11-21T17:05:19Z")

</div>

I think this _should_ work, but for me it segfaults ☹

```julia
using Turing,StatsPlots

@model function model9(k01, k10, k11)
    p0 ~ Beta(1, 1)
    p1 ~ Beta(1, 1)
    num_seen = k01 + k10 + k11
    N ~ DiscreteUniform(num_seen, 350)
    
    q0 = 1-p0
    q1 = 1-p1
    ps = [q0*q1, q0*p1, p0*q1, p0*p1]
    
    y = [N - num_seen, k01, k10, k11]
    y ~ Multinomial(N, ps)
end

k10 = 20-3
k01 = 15-3
k11 = 3

trace9 = sample(model9(k01, k10, k11), Gibbs(NUTS(1000,.8,:p0,:p1),PG(20,:N)),1000)

density(trace9[:N])

```

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 5:35pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/16 "2021-11-21T17:35:28Z")

</div>

> [@dlakelan](#):
>
> ```julia
> trace9 = sample(model9(k01, k10, k11), Gibbs(NUTS(1000,.8,:p0,:p1),PG(20,:N)),1000)
> 
> ```

I’m getting a gazillion warnings:

```julia
Warning: The current proposal will be rejected due to numerical error(s).
│ isfinite.((θ, r, ℓπ, ℓκ)) = (true, true, false, true)
└ @ AdvancedHMC ~/.julia/packages/AdvancedHMC/HQHnm/src/hamiltonian.jl:47

```

The sampled distributions are rubbish ☹

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [November 21, 2021, 5:45pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/17 "2021-11-21T17:45:56Z")

</div>

You may need to first optimize the initial position. Perhaps pick a value from the sample on python just as a way to test if you can initialize it and make it work

The keyword argument for the initial position is wrong in the documentation it should be `init_params=`

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 6:29pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/19 "2021-11-21T18:29:25Z")

</div>

> [@dlakelan](#):
>
> You may need to first optimize the initial position. Perhaps pick a value from the sample on python just as a way to test if you can initialize it and make it work

OK, tried to cheat, and pass initial parameters that are very close to their mean (according to PyMC3):

```julia
trace9 = sample(model9(k01, k10, k11), Gibbs(NUTS(1000,.8,:p0,:p1),PG(20,:N)), 1000, init_params=[0.2, 0.2, 100])

```

I’m still getting the warnings and results are even worse:

 ![Screen Shot 2021-11-21 at 20.26.17](https://global.discourse-cdn.com/julialang/original/3X/b/3/b350d2d578c6eb42ca7959db0ea9d280bae598e4.jpeg)

This does not strike me as an healthy chain (but at this point it should be clear that I’m not an expert in MCMC, Turing or PyMC3)

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [November 21, 2021, 6:33pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/20 "2021-11-21T18:33:44Z")

</div>

Definitely not a healthy chain but choosing the mean can be a bad idea it can actually be far from a typical value just try choosing the last sample used by pymc

meanwhile I’m updating some of my packages in hopes to eliminate the segfault so I can try it out. I’ve been meaning to try Gibbs mixed continuous / discrete sampling anyway.

---

<div class="post-metadata">

**Author:** ![dodoplus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dodoplus/32/30148_2.png) [@dodoplus](https://discourse.julialang.org/u/dodoplus)\
**Post date:** [November 21, 2021, 7:14pm UTC](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856/21 "2021-11-21T19:14:10Z")

</div>

OK, I’ve started the original on Colab:

> **[Google Colaboratory](https://colab.research.google.com/github/AllenDowney/ThinkBayes2/blob/master/soln/chap19.ipynb)**

Things aren’t that great here either apparently:

```julia
The rhat statistic is larger than 1.4 for some parameters. The sampler did not converge.
The estimated number of effective samples is smaller than 200 for some parameters.

```

(would be great if Turing gave similar commentary…)

Latest traces are:

```julia
{'p0_logodds __': -1.1201974288020848, 'p1_logodds__': -1.8108031775461457, 'N': 94, 'p0': 0.24597466437185683, 'p1': 0.1405410824792338}
{'p0_logodds __': -1.0759843501248898, 'p1_logodds__': -1.3370046772297803, 'N': 95, 'p0': 0.25426669304142646, 'p1': 0.208003069852844}
{'p0_logodds __': -1.309490295895163, 'p1_logodds__': -2.229973878615384, 'N': 97, 'p0': 0.21257214856225384, 'p1': 0.09709093087482638}
{'p0_logodds __': -1.2321949553810507, 'p1_logodds__': -2.4305296251833313, 'N': 96, 'p0': 0.22579748811056297, 'p1': 0.0808740895970291}

```

which aren’t far from the initial parameters I’ve tried (`p0=0.2, p1=0.2, N=100`).

I’ve tried:

```julia
trace9 = sample(model9(k01, k10, k11), Gibbs(NUTS(1000,.8,:p0,:p1),PG(20,:N)), 1000, init_params=[0.22, 0.08, 95])

```

But no joy:

 ![Screen Shot 2021-11-21 at 21.13.05](https://global.discourse-cdn.com/julialang/original/3X/3/4/34edb35bc207ead0fcab3499b4dbd76e96dab994.jpeg)

 ![Screen Shot 2021-11-21 at 21.13.33](https://global.discourse-cdn.com/julialang/original/3X/f/e/fe8c4e51a95aaf206112553c3470ff7318902d87.jpeg)

[Next page](https://discourse.julialang.org/t/trouble-translating-a-simple-pymc3-model-to-turing/71856.md?page=2)