Accessing samples

FlexiChains stores data in a dict-of-array format (a column-oriented form, similar to dataframes). This means that it is very easy to access all the data for a given variable, but accessing all the data for a given sample (i.e., iteration / chain number) is a bit more involved.

Consider our now-familiar example chain:

using FlexiChainschn = FlexiChains.EXAMPLE_CHAIN

╭─FlexiChain (10 iterations, 2 chains) ────────────────────────────────────────╮
│ ↓ iter  = 101:110                                                            │
│ → chain = 1:2                                                                │
│                                                                              │
│ Parameters (3) ── VarName                                                    │
│  Float64          mu, tau                                                    │
│  Vector{Float64}  theta (8,)                                                 │
│                                                                              │
│ Extras (2)                                                                   │
│  Float64  logjoint                                                           │
│  Bool     is_divergent                                                       │
╰──────────────────────────────────────────────────────────────────────────────╯

Here, we have three parameters mu, tau, and theta. FlexiChains makes it very easy to access the value of a single parameter across all iterations and chains:

chn[@vn(mu)]

┌ 10×2 DimArray{Float64, 2} Parameter(mu) ┐
├─────────────────────────────────────────┴─────────────── dims ┐
  ↓ iter Sampled{Int64} 101:110 ForwardOrdered Regular Points,
  → chain Sampled{Int64} 1:2 ForwardOrdered Regular Points
└───────────────────────────────────────────────────────────────┘
   ↓ →   1         2
 101    11.3517    4.78008
 102    11.3517    4.59214
 103    11.3517    5.68752
 104    11.3517    7.10065
   ⋮               
 107     1.30113   5.48149
 108     0.915619  6.2456
 109     0.284842  6.25491
 110     0.516548  1.47832

However, suppose you want to obtain the values of all parameters at a given iteration and chain number. For example, you might want to know what happened at the fifth iteration.

One way is to subset a chain. The issue is that that gives you just another chain back:

subsetted = chn[iter=5, chain=1]

╭─FlexiChain (1 iteration, 1 chain) ───────────────────────────────────────────╮
│ ↓ iter  = 105:105                                                            │
│ → chain = 1:1                                                                │
│                                                                              │
│ Parameters (3) ── VarName                                                    │
│  Float64          mu, tau                                                    │
│  Vector{Float64}  theta (8,)                                                 │
│                                                                              │
│ Extras (2)                                                                   │
│  Float64  logjoint                                                           │
│  Bool     is_divergent                                                       │
╰──────────────────────────────────────────────────────────────────────────────╯

and to access the data you will still need to index into it, which is rather awkward:

(mu=only(subsetted[@vn(mu)]), tau=only(subsetted[@vn(tau)]))

(mu = 11.351650969356333, tau = 0.846032171561763)

We should definitely like to have a better way of accessing this collection of values at a given iteration/chain.

Generally speaking, the operation above can be viewed as a transformation from column-oriented data to row-oriented data, or from a dict-of-arrays to an array-of-dicts. FlexiChains provides a few convenient ways for you to do this transformation.

`values_at`

The function FlexiChains.values_at is at the core of this transformation. Given iteration and chain indices (as keyword arguments), it returns some container that holds all the values for a given iteration.

draw = FlexiChains.values_at(chn; iter=5, chain=1)

VarNames.Draw
 ├─ parameters
 │  VarNamedTuple
 │  ├─ mu => 11.351650969356333
 │  ├─ tau => 0.846032171561763
 │  └─ theta => [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468]
 └─ extras
    ├─ logjoint = -45.082557361617546
    └─ is_divergent = false

There is a subtlety here in that the exact type of the container returned depends on what is stored in the chain.

Structures

For the interested reader, or if you are developing FlexiChains, this is determined by the structures stored in the chain, which are meant to be a 'template' for the container. If you are just a user who has obtained a FlexiChain from some library, should not need to worry about this, as the structures should be set up correctly for you by that library.

For chains that have parameter type VarName, i.e. a VNChain, the returned container should usually be a VarNames.Draw. This is a struct that separately stores the parameters and the stats as a VarNames.VarNamedTuple and a NamedTuple respectively.

Not getting a `Draw`?

If you instead get an OrderedDict, you should probably report this as a bug to the library that produced the chain.

For more information about these types, the reader is referred to the VarNames.jl docs. The key point here is that this is a high-fidelity representation of the data. You can access the parameters and extras separately:

ps = FlexiChains.parameters(draw)

VarNamedTuple
├─ mu => 11.351650969356333
├─ tau => 0.846032171561763
└─ theta => [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468]

FlexiChains.extras(draw)

(logjoint = -45.082557361617546, is_divergent = false)

The parameter VarNamedTuple can be further indexed into with VarNames to obtain the values of individual parameters. Much like a FlexiChain itself, you have the ability to index into a VarNamedTuple with sub-variables, such as theta[1:3]:

ps[@vn(mu)], ps[@vn(theta[1:3])]

(11.351650969356333, [10.492901852254972, 11.46635882152557, 11.570533934319785])

Choosing a different container type

Sometimes, though, one might want different output formats. There is some support for converting to NamedTuple or AbstractDict, by passing an optional type argument.

Conversion to NamedTuple is not faithful

Note that conversion to NamedTuple is lossy if you have VarNames that contain indexing or field access syntax (e.g., x[1] or x.a).

FlexiChains.values_at(chn, NamedTuple; iter=5, chain=1)

(mu = 11.351650969356333, tau = 0.846032171561763, theta = [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468], logjoint = -45.082557361617546, is_divergent = false)

FlexiChains.values_at(chn, Dict; iter=5, chain=1)

Dict{Union{FlexiChains.Extra, FlexiChains.Parameter{var"#s91"}} where var"#s91"<:VarName, Any} with 5 entries:
  Parameter(tau)       => 0.846032
  Parameter(mu)        => 11.3517
  Extra(:is_divergent) => false
  Extra(:logjoint)     => -45.0826
  Parameter(theta)     => [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114,…

Obtaining only parameters

Often the stats are not of very much interest, and you just want the parameters. In that case, you can use parameters_at:

FlexiChains.parameters_at(chn; iter=5, chain=1)

VarNamedTuple
├─ mu => 11.351650969356333
├─ tau => 0.846032171561763
└─ theta => [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468]

This is exactly the same as values_at, except that it only returns the parameters.

Not getting a VarNamedTuple

Just as above, if you instead get an OrderedDict, you should probably report this as a bug to the library that produced the chain.

Arrays of samples

The functions above also accept multiple iteration or chain indices, and will return the corresponding array of containers.

FlexiChains.parameters_at(chn; iter=[5, 6], chain=1)

┌ 2-element DimArray{VarNames.VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, 1} ┐
├──────────────────────────────────────────────────────────────────────── dims ┤
  ↓ iter Sampled{Int64} [105, 106] ForwardOrdered Irregular Points
└──────────────────────────────────────────────────────────────────────────────┘
 105  …  VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 106     VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])

All of DimensionalData's selectors also work:

FlexiChains.parameters_at(chn; iter=Not(1..8), chain=1)

┌ 10-element DimArray{VarNames.VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, 1} ┐
├──────────────────────────────────────────────────────────────────────── dims ┤
  ↓ iter Sampled{Int64} [101, …, 110] ForwardOrdered Irregular Points
└──────────────────────────────────────────────────────────────────────────────┘
 101  …  VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 102     VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 103     VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 104     VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 105     VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 106  …  VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582])
 107     VarNamedTuple(mu = 1.30113, tau = 1.25758, theta = [-0.278603, 2.20435, 1.4987, 1.25086, 1.20394, 2.61129, 1.21561, 1.20296])
 108     VarNamedTuple(mu = 0.915619, tau = 0.973233, theta = [-0.312388, 1.68678, -0.409465, -0.343761, 0.0934012, 0.613312, 0.966749, 0.965976])
 109     VarNamedTuple(mu = 0.284842, tau = 1.82906, theta = [1.86213, -1.40004, 3.83221, -1.47719, 1.469, 1.02763, 0.180669, -0.839965])
 110     VarNamedTuple(mu = 0.516548, tau = 1.84247, theta = [2.57074, -2.26718, 3.4244, -0.980078, 0.672412, 0.74598, 1.65239, -1.42801])

This also means that you can convert a VNChain back into an Array of Draw objects by passing : for both the iteration and chain indices.

FlexiChains.values_at(chn; iter=:, chain=:)

┌ 10×2 DimArray{VarNames.Draw{VarNames.VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}, 2} ┐
├──────────────────────────────────────────────────────────────────────── dims ┤
  ↓ iter Sampled{Int64} 101:110 ForwardOrdered Regular Points,
  → chain Sampled{Int64} 1:2 ForwardOrdered Regular Points
└──────────────────────────────────────────────────────────────────────────────┘
   ↓ →  …  2
 101        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 4.78008, tau = 2.78233, theta = [5.06589, 7.0214, 5.71944, 6.71881, 2.90095, -2.12359, 3.62712, 5.59417]), (logjoint = -54.7921, is_divergent = false))
 102        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 4.59214, tau = 3.14813, theta = [5.82169, 6.69599, 5.35547, 6.51507, 2.97696, -1.63556, 3.65501, 6.79247]), (logjoint = -54.5201, is_divergent = false))
 103        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 5.68752, tau = 3.19854, theta = [9.76186, 8.79631, 6.02066, 8.1514, 2.35677, 4.94412, 1.58788, 6.26805]), (logjoint = -54.9841, is_divergent = false))
 104        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 7.10065, tau = 2.62221, theta = [12.5662, 6.24494, 7.5854, 5.73067, 4.84375, 7.32268, 5.19098, 7.05275]), (logjoint = -53.2478, is_divergent = false))
   ⋮    ⋱  
 107        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 5.48149, tau = 2.44535, theta = [2.74071, 5.78394, 8.56174, 8.99474, 5.77947, 7.40668, 5.41739, 5.32562]), (logjoint = -53.045, is_divergent = false))
 108        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 6.2456, tau = 4.63625, theta = [15.5577, 8.86566, 2.77102, 4.54069, 8.30045, 4.91825, 11.8373, 8.66637]), (logjoint = -57.7038, is_divergent = false))
 109        Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 6.25491, tau = 7.48281, theta = [8.97417, 8.43466, 2.55279, -8.15037, 8.45557, 0.756357, 18.4144, 3.61016]), (logjoint = -63.576, is_divergent = false))
 110    …   Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 1.47832, tau = 4.93959, theta = [6.06009, 4.07403, -0.727583, 6.03801, 1.85736, -4.66044, 2.52889, 12.7658]), (logjoint = -59.8085, is_divergent = false))

Drawing random samples

To get a random sample from the chain, you can use Base.rand. The container type returned is the same as for values_at:

rand(chn)

VarNames.Draw
 ├─ parameters
 │  VarNamedTuple
 │  ├─ mu => 11.351650969356333
 │  ├─ tau => 0.846032171561763
 │  └─ theta => [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468]
 └─ extras
    ├─ logjoint = -45.082557361617546
    └─ is_divergent = true

This method follows Julia's conventions closely, so you can have an optional rng first argument, and you can also specify the number of samples to draw:

using Random: Xoshirorand(Xoshiro(468), chn, 2, 2)

2×2 Matrix{VarNames.Draw{VarNames.VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}}:
 Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582]), (logjoint = -45.0826, is_divergent = false))  …  Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 4.78008, tau = 2.78233, theta = [5.06589, 7.0214, 5.71944, 6.71881, 2.90095, -2.12359, 3.62712, 5.59417]), (logjoint = -54.7921, is_divergent = false))
 Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 11.3517, tau = 0.846032, theta = [10.4929, 11.4664, 11.5705, 12.131, 10.9676, 12.0114, 10.794, 11.582]), (logjoint = -45.0826, is_divergent = false))     Draw{VarNamedTuple{(:mu, :tau, :theta), Tuple{Float64, Float64, Vector{Float64}}}, @NamedTuple{logjoint::Float64, is_divergent::Bool}}(VarNamedTuple(mu = 5.48149, tau = 2.44535, theta = [2.74071, 5.78394, 8.56174, 8.99474, 5.77947, 7.40668, 5.41739, 5.32562]), (logjoint = -53.045, is_divergent = false))

Weighted sampling is also supported (see the docstring below for more info).

If you only want parameters, pass the parameters_only=true keyword argument:

rand(chn, parameters_only=true)

VarNamedTuple
├─ mu => 11.351650969356333
├─ tau => 0.846032171561763
└─ theta => [10.492901852254972, 11.46635882152557, 11.570533934319785, 12.130958608708998, 10.967643507906876, 12.011427174304124, 10.793983473375265, 11.581970127184468]

Docstrings

FlexiChains.values_at — Function

FlexiChains.values_at(chn::FlexiChain; iter=:, chain=:)
FlexiChains.values_at(chn::FlexiChain, ::Type{Tout}; iter=:, chain=:)

Extract all values from the chain corresponding to a particular set of MCMC iterations(s).

The iter and chain keyword arguments can be anything used to index into the respective dimensions of a FlexiChain, such as an integer, a vector of integers, a range, or a DimensionalData.jl selector. Both default to : (i.e. all iterations / all chains).

In particular, you can convert the entire chain into a DimMatrix of the desired output type by calling this function with no keyword arguments.

To get only the parameter keys, use FlexiChains.parameters_at.

The output type can be specified with an (optional) positional argument. Possible options are:

unspecified: uses the structure stored in the chain to determine the output type. Specifically, for chains sampled with Turing.jl, the stored structure will be VarNamedTuple, and consequently the output of parameters_at will be a DynamicPPL.ParamsWithStats per iteration (which stores the parameters as a VarNamedTuple, and the stats/extras as a NamedTuple).
For chains that do not have a stored structure, the default output will be an OrderedDict, or NamedTuple for TKey = Symbol.
Tout <: AbstractDict: returns a dictionary mapping ParameterOrExtra{TKey} to their values. This is the most faithful representation of the data in the chain.
Tout = NamedTuple, or Tout = ComponentArray: attempts to convert every key name to a Symbol, which is used as the field name in the output NamedTuple or ComponentArray.

Using `NamedTuple` or `ComponentArray`

This will throw an error if any key cannot be converted to a Symbol, or if there are duplicate key names after conversion. If you have parameter names that convert to the same Symbol, you can either use OrderedDict, subset the chain before calling this function, or rename your parameters. Furthermore, please be aware that this is a lossy conversion as it does not retain information about whether a key is a parameter or an extra.

For order-sensitive output types, such as OrderedDict, the keys are returned in the same order as they are stored in the FlexiChain. This also corresponds to the order returned by keys(chn).

source

FlexiChains.parameters_at — Function

FlexiChains.parameters_at(chn::FlexiChain; iter=:, chain=:)
FlexiChains.parameters_at(chn::FlexiChain, ::Type{Tout}; iter=:, chain=:)

Extract all parameter values from the chain corresponding to a particular set of MCMC iteration(s), discarding non-parameter (i.e. Extra) keys.

In particular, you can convert the entire chain into a DimMatrix of the desired output type by calling this function with no keyword arguments.

To get all keys (not just parameters), use FlexiChains.values_at.

The output type for each iteration can be specified with the first (optional) positional argument. Possible options are:

unspecified: uses the structure stored in the chain to determine the output type. Specifically, for chains sampled with Turing.jl, the stored structure will be VarNamedTuple, and consequently the output of parameters_at will be a VarNamedTuple per iteration.
For chains that do not have a stored structure, the default output will be a NamedTuple if the key type is Symbol, and an OrderedDict otherwise.
Tout <: AbstractDict: returns a dictionary mapping TKey to their values.
Tout = NamedTuple or Tout <: ComponentArray: attempts to convert every parameter name to a Symbol, which is used as the field name in the output NamedTuple or ComponentArray.

Using `NamedTuple` or `ComponentArray`

For order-sensitive output types, such as OrderedDict, the parameters are returned in the same order as they are stored in the FlexiChain. This also corresponds to the order returned by FlexiChains.parameters(chn).

source

Base.rand — Function

Base.rand(rng::Random.AbstractRNG, chn::FlexiChain, dims::Int...; parameters_only=false, weights=nothing)

Sample uniformly from a FlexiChain with replacement.

Pass weights::StatsBase.AbstractWeights to sample with probabilities proportional to the weights instead. Supply one weight per sample, with all iterations of the first chain followed by all iterations of the second chain, and so on.

rand([rng,] chn) returns a single sample, whereas rand([rng,] chn, dims...) returns an Array of samples with dimensions dims.... The parameters_only keyword specifies whether the returned samples include only parameters, or both parameters and extras.

In general, the return type of rand is the same as the return type of FlexiChains.parameters_at or FlexiChains.values_at (but wrapped in an Array if dims is not empty). This means that the return type can depend on the 'structure' stored in the FlexiChain.

For example, if chn::FlexiChain{VarName} was constructed using Turing.jl, rand(chn) will return a DynamicPPL.ParamsWithStats, and rand(chn, parameters_only=true) will return a DynamicPPL.VarNamedTuple.

source