Interactive Explainer

Probabilistic Programming, in Practice

PyMC, NumPyro, Stan — the productive recipe for Bayesian modelling. Declare the model, pick an inference engine, inspect diagnostics. Two live demos in the browser plus the checklist of what to look at after every run.

~20 min Bayesian · MCMC · VI

Step 1

What a PPL is — and isn't

A probabilistic programming language lets you write a generative model as code. You declare priors, sampling distributions, observations; the runtime takes care of posterior inference. The same model can be run with HMC, NUTS, ADVI, SMC, or Laplace approximation depending on a one-line switch.

Tool	Backend	Strengths	Quirks
Stan	Custom + C++	Best NUTS sampler; rock-solid for hierarchical models	Different language; slower compile
PyMC	PyTensor / JAX	Pythonic; large user base; rich tools (ArviZ)	Compile time on first run
NumPyro	JAX	Fast NUTS via JAX; GPU-friendly; cleanest deep-learning interop	Smaller community than PyMC
Turing.jl	Julia	Composes with Julia's ML stack; flexible	Julia-only
Pyro	PyTorch	Stochastic VI on big models; good for BNNs / latent-variable	HMC less mature than NumPyro's
Edward2 / TFP	TensorFlow	Scales; ecosystem hooks	API churn

The slogan. You write the math once. The PPL handles the inference. Switching between HMC and VI is a flag change, not a rewrite.

Step 2

Pseudo-code template

with model:
    theta = Beta("theta", alpha=1, beta=1)        # prior
    y = Binomial("y", n=20, p=theta, observed=14) # likelihood + obs

    trace = sample(2000, chains=4)                 # NUTS
    az.plot_trace(trace)
    az.summary(trace)                              # r-hat, ESS

Three things every PPL gives you:

Declarative model. Variables are named, distributions have priors, observations attach to data.
Inference engine. NUTS / HMC by default; VI / Laplace on demand.
Diagnostics. r-hat (between-chain convergence), ESS (effective sample size), divergent transitions, posterior predictive checks.

Step 3

Live: coin-bias model — HMC vs analytic

The simplest Bayesian model: $\theta \sim \mathrm{Beta}(1, 1)$, $y \sim \mathrm{Bernoulli}(\theta)$. Beta-Binomial conjugacy gives the analytic posterior; an HMC chain should match it. Add data with the buttons and run HMC live.

k / n 0/0 r-hat — ESS —

Posterior: HMC samples vs analytic Beta

HMC trace plot (chain 1)

Step 4

Live: 2-D Bayesian regression — HMC vs VI

A Bayesian linear regression with a small dataset. HMC gives honest posterior; mean-field VI fits a factorised Gaussian. Watch the credible bands differ: VI is usually too narrow, especially in directions where the parameters correlate.

n —

HMC posterior over (slope, intercept)

VI Gaussian fit (mean-field)

The VI trap. Mean-field VI factorises across parameters, so it can't represent correlated posteriors. The HMC contours often slant; the VI contours always align with the axes. For credible-interval-quality work, prefer HMC or full-rank VI.

Step 5

The diagnostic checklist

r-hat < 1.01. Between-chain consensus on each parameter. Above 1.05 = chains haven't converged; run longer or reparameterise.
ESS > 400. Effective sample size. Below this, your posterior summary is noisy.
Divergent transitions. NUTS-specific. Even 5 divergences ≥ 0% of the chain means the geometry is too curved — try a non-centred parameterisation or increase target_accept.
Posterior predictive checks. Simulate new data from the posterior and compare to observed. Bad fit here ≫ any single diagnostic.
Prior predictive checks. Sample from the prior alone. If those samples look implausible, your prior is uninformative-in-a-bad-way.
LOO / WAIC. Cross-validated predictive accuracy for model comparison. ArviZ ships both.

Step 6

What PPLs are for (and not for)

Good fits. Hierarchical models, small data, calibrated uncertainty, parameter interpretation, scientific models with mechanistic priors.
Reach for instead: deep learning on millions of examples (gradient descent is faster), pure point-estimation (sklearn / XGBoost), production serving (PPLs are great in research, awkward in production).
The deep BNN gap. Bayesian neural nets are hard. PyMC + NUTS chokes past a few thousand parameters; VI on Bayesian-deep is biased. SWAG, Laplace, deep ensembles, and SGLD are the practical compromises.
Cross-link. The Bayesian posterior playground covers MCMC vs VI from a different angle; the Neural Processes article shows the amortised-inference alternative.