# What is the canon way to randomly draw one (or more) parameter sets from a sampled distribution?
**URL:** https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935
**Category:** Probabilistic Programming
**Tags:** turing
**Created:** [March 1, 2026, 10:22am UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935 "2026-03-01T10:22:54Z")
**Posts on this page:** 12
**Page:** 1
### Author:  [@Torkel](https://discourse.julialang.org/u/Torkel)
#### Post date: [March 1, 2026, 10:22am UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/1 "2026-03-01T10:22:55Z")
Using turing I have done
```julia
n_steps = 1000
n_chains = 4
chain = sample(model, NUTS(), MCMCThreads(), n_steps, n_chains)
```
I can plot the output using
```julia
plot(chain)
```
However, if I want to draw a random parameter set to do something with, is there a canon way to do this? I.e. here I have done inference for parameters `γ, ν, σI` and want to draw a single value from that joint distribution. I practise, I think I could use
```julia
n_ps = 3
collect(chain.value[rand(1:n_steps), :, rand(1:n_chains)])[1:n_ps]
```
however, I also suspect that there might exist a more direct interface.
---
### Author:  [@eteppo](https://discourse.julialang.org/u/eteppo)
#### Post date: [March 1, 2026, 3:16pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/2 "2026-03-01T15:16:50Z")
With the `Array` constructor you can at least avoid accessing the `value` field, e.g. `Array(chain, append_chains=false)[rand(1:n_chains)][rand(1:n_steps), :]` ([Getting started · MCMCChains.jl](https://turinglang.org/MCMCChains.jl/stable/getting-started/#Exporting-Chains)). If the chains can be pooled together, it is `Array(chain)[rand(1:n_steps), :]`.
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [March 1, 2026, 5:11pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/3 "2026-03-01T17:11:57Z")
you can do:
```julia
sample(chain[:γ], 10)
```
---
### Author:  [@Torkel](https://discourse.julialang.org/u/Torkel)
#### Post date: [March 1, 2026, 7:06pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/4 "2026-03-01T19:06:02Z")
Thanks, works great! If I want to sample all three parameters, is there some version of
```julia
sample(chain[:γ, :ν, :σI], 1)[1] # errors
```
that I can do (which is a bit nicer than
```julia
[sample(chain[:γ], 1)[1], sample(chain[:ν], 1)[1] , sample(chain[:σI], 1)[1] ]
```
?
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [March 1, 2026, 7:16pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/5 "2026-03-01T19:16:36Z")
you could do `sample.((chain[:γ], chain[:ν]), 10)`
---
### Author:  [@Torkel](https://discourse.julialang.org/u/Torkel)
#### Post date: [March 1, 2026, 7:24pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/6 "2026-03-01T19:24:24Z")
That works, thanks a lot!
---
### Author:  [@Torkel](https://discourse.julialang.org/u/Torkel)
#### Post date: [March 1, 2026, 7:26pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/7 "2026-03-01T19:26:12Z")
> [@langestefan](#):
>
> .((chain[:γ], chain[:ν]), 10)
Just checking, when I do this, the different parameter values are correctly correlated, i.e. I am not sampling `γ` and `ν` values independently of each other, right?
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [March 1, 2026, 7:29pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/8 "2026-03-01T19:29:09Z")
You didn’t share your model, but I don’t see how `γ` and `ν` can be correlated?
---
### Author:  [@Torkel](https://discourse.julialang.org/u/Torkel)
#### Post date: [March 1, 2026, 7:32pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/9 "2026-03-01T19:32:39Z")
For reference, this is the full workflow:
```julia
using Catalyst
sir = @reaction_network begin
γ, S + I --> 2I
ν, I --> R
end
using OrdinaryDiffEqDefault, Plots
u0 = [:S => 999.0, :I => 1.0, :R => 0.0]
p_true = [:γ => 0.0005, :ν => 0.1]
oprob_true = ODEProblem(sir, u0, 100.0, p_true)
sol_true = solve(oprob_true)
plot(sol_true; label = ["S (true)" "I (true)" "R (true)"], lw = 4)
using Distributions
σI = 10.0
t_measurement = 1.0:2:100.0
I_observed = sol_true(t_measurement; idxs = :I)
I_observed = [rand(Normal(I, σI)) for I in I_observed]
plot!(t_measurement, I_observed; label = "I (measured)", color = 2)
using Turing, SymbolicIndexingInterface
setp_oop = SymbolicIndexingInterface.setp_oop(oprob_true, [:γ, :ν])
@model function sir_likelihood(I_observed, oprob_base, setp_oop, saveat)
# Defines the parameters we wish to estimate and their prior distributions.
γ ~ LogUniform(0.00001, 0.001)
ν ~ LogUniform(0.01, 1.0)
σI ~ LogUniform(0.1, 100.0)
# Simulate the model for parameter values γ, ν. Saves the solution at the measurement times.
p = setp_oop(oprob_base, [γ, ν])
oprob_base = remake(oprob_base; p)
sol = solve(oprob_base; saveat, verbose = false, maxiters = 10000)
# If simulation was unsuccessful, the likelihood is -Inf.
if !SciMLBase.successful_retcode(sol)
Turing.@addlogprob! -Inf
return nothing
end
# Computes the likelihood of the observations.
for idx in eachindex(I_observed)
I_observed[idx] ~ Normal(sol[:I][idx], σI)
end
end
n_steps = 1000
n_chains = 4
sir_model = sir_likelihood(I_observed, oprob_true, setp_oop, t_measurement)
chain = sample(sir_model, NUTS(), MCMCThreads(), n_steps, n_chains; progress = false)
```
The result should be a 3d posterior over the estimated parameters γ, ν, and σI, respectively(?). Imagine there is a correlation between γ and ν in this distribution. In this case, wouldn’t sampling a value of γ and one of ν independently yield a different result from sampling a single point in my posterior distribution and extracting its γ and ν values? I.e. there there are two different interpretations of the sampled point
```julia
first.(sample.((chain[:γ], chain[:ν]), 1))
```
where I am uncertain which of them it is.
---
### Author:  [@evan03](https://discourse.julialang.org/u/evan03)
#### Post date: [March 2, 2026, 6:39am UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/10 "2026-03-02T06:39:47Z")
This post was temporarily hidden by the community for possibly being off-topic, unfocused, inappropriate, or spammy.
---
### Author:  [@eteppo](https://discourse.julialang.org/u/eteppo)
#### Post date: [March 2, 2026, 2:05pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/11 "2026-03-02T14:05:07Z")
Yeah, just use the whole posterior, in general the posterior could be any shape and MCMC is run to infer it without too many assumptions. In any case, pairwise scatterplots are a handy check to go with the check on the mixing and convergence of chains.
```julia-auto
using PairPlots, CairoMakie
pairplot(chain)
```
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [March 2, 2026, 5:48pm UTC](https://discourse.julialang.org/t/what-is-the-canon-way-to-randomly-draw-one-or-more-parameter-sets-from-a-sampled-distribution/135935/12 "2026-03-02T17:48:15Z")
The method I shared samples directly from the marginals. So yeah in that case the correlation is lost. I think in your case you want to sample directly from the chain. But maybe @penelopeysm or someone else from turing can better help with that (there is also a turing channel on slack)