I cannot make use of all summary statistics, as I am not sure what they mean. Could you maybe give some explanations in the documentary or in this tutorial? I would appreciate it.
Originally posted by **liesel-groupie** November 14, 2023
Dear Liesels,
The code for the tutorial "Parameter transformation" has some functionalities from older Liesel versions. Here is a small update:
```python
import jax
import jax.numpy as jnp
import numpy as np
import matplotlib.pyplot as plt
# We use distributions and bijectors from tensorflow probability
import tensorflow_probability.substrates.jax.distributions as tfd
import tensorflow_probability.substrates.jax.bijectors as tfb
import liesel.goose as gs
import liesel.model as lsl
rng = np.random.default_rng(42)
# data-generating process
n = 500
true_beta = np.array([1.0, 2.0])
true_sigma = 1.0
x0 = rng.uniform(size=n)
X_mat = np.c_[np.ones(n), x0]
y_vec = X_mat @ true_beta + rng.normal(scale=true_sigma, size=n)
# Model
# Part 1: Model for the mean
beta_prior = lsl.Dist(tfd.Normal, loc=0.0, scale=100.0)
beta = lsl.Param(value=np.array([0.0, 0.0]), distribution=beta_prior,name="beta")
X = lsl.Obs(X_mat, name="X")
mu = lsl.Var(lsl.Calc(jnp.dot, X, beta), name="mu")
# Part 2: Model for the standard deviation
a = lsl.Var(0.01, name="a")
b = lsl.Var(0.01, name="b")
sigma_sq_prior = lsl.Dist(tfd.InverseGamma, concentration=a, scale=b)
sigma_sq = lsl.Param(value=10.0, distribution=sigma_sq_prior, name="sigma_sq")
sigma = lsl.Var(lsl.Calc(jnp.sqrt, sigma_sq), name="sigma")
# Observation model
y_dist = lsl.Dist(tfd.Normal, loc=mu, scale=sigma)
y = lsl.Var(y_vec, distribution=y_dist, name="y")
gb = lsl.GraphBuilder().add(y)
gb.transform(sigma_sq, tfb.Exp)
model = gb.build_model()
builder = gs.EngineBuilder(seed=1339, num_chains=4)
builder.set_model(lsl.GooseModel(model))
builder.set_initial_values(model.state)
builder.add_kernel(gs.NUTSKernel(["beta", "sigma_sq_transformed"]))
builder.set_duration(warmup_duration=1000, posterior_duration=1000)
# by default, goose only stores the parameters specified in the kernels.
# let's also store the standard deviation on the original scale.
builder.positions_included = ["sigma_sq"]
engine = builder.build()
engine.sample_all_epochs()
results = engine.get_results()
# what is the output here?
summary = gs.Summary(results).to_dataframe().reset_index()
```
However, I cannot make use of all summary statistics, as I am not sure what they mean. Could you maybe give some explanations in the documentary or in this tutorial? I would appreciate it.
Main point:
Discussed in https://github.com/liesel-devs/liesel/discussions/154