Betancourt 2015 - Hierarchical Geometry and the Funnel

Summary

A hierarchical model groups local parameters ${\theta }_i$ under shared global parameters $\varphi$. With sparse data the posterior is a funnel: a region of high density but low volume below a region of low density but high volume, with comparable probability mass and sharply varying curvature. Because the correlation between the log-scale $v$ and each ${\theta }_i$ flips strongly with position, no global rotation or rescaling removes it.

Overview

The exposition treats the one-level hierarchy as typical of all hierarchical pathologies. The footnote notes that not all parameters need share one hierarchy (multilevel models).

Main Content

One-level hierarchical model (eq. 1, Fig. 1)

Data $D=(D_i)$, local parameters $\theta =({\theta }_i)$, global parameters $\varphi$:

$$\pi (\theta ,\varphi \mid D)\propto \prod _{i=1}^n\pi (D_i\mid {\theta }_i)\pi ({\theta }_i\mid \varphi )\pi (\varphi ).$$

The ${\theta }_i$ interact only through $\varphi$, which lets data inform all of the ${\theta }_i$ (partial pooling).

One-way normal model (eq. 2)

$y_i\sim \mathcal{N}({\theta }_i,{\sigma }_i^2)$, ${\theta }_i\sim \mathcal{N}(\mu ,{\tau }^2)$, $i=\mathrm{1..}I$. In the general notation $D=(y_i,{\sigma }_i)$, $\varphi =(\mu ,\tau )$, $\theta =({\theta }_i)$. “Global” and “local” labels break down in multi-layer models.

Why it is pathological. All $n$ contributions at the bottom depend on $\varphi$, so a small change in $\varphi$ produces large changes in density. With sparse data the density has a funnel shape. Both regions carry comparable probability mass, so a successful sampler must cope with dramatic variation in curvature.

The funnel distribution (Neal 2003; Fig. 2)

With no data, $\mu =0$ and a log-normal prior on the variance ${\tau }^2=e^v$:

$$\pi ({\theta }_1,\dots ,{\theta }_n,v)\propto \prod _{i=1}^n\mathcal{N}({\theta }_i\mid 0,(e^{-v/2}{)}^2)\mathcal{N}(v\mid 0,3^2).$$

Notation note: the formula as printed in Sec. I uses $x_i$ for the local variables and scale $e^{-v/2}$, while the Appendix Stan model is theta ~ normal(0, exp(v/2)). The two differ by $v\to -v$ and describe the same funnel up to reflection; $v\sim \mathcal{N}(0,3^2)$ marginally either way. Footnote 2: the exponential link between latent $v$ and variance arises naturally whenever a positive parameter is transformed to an unconstrained one for sampling.

Position-dependent correlation

The hierarchy induces large correlations between $v$ and each ${\theta }_i$, with a strength that varies strongly with position (Fig. 2 visualises the eigenvectors of $\sqrt{|{\partial }^2\log p/\partial {\theta }_i\partial {\theta }_j|}$ scaled by eigenvalues). The funnel’s symmetry makes the correlations cancel globally, so no global correction (rotation or rescaling of parameters) can simplify the distribution.

Reproducing the funnel in numpy

import numpy as np
rng = np.random.default_rng(0)
n, N = 100, 100000
v = rng.normal(0, 3, N)                          # marginal of v
theta = rng.normal(0, np.exp(v/2)[:, None], (N, n))  # Stan model: theta ~ normal(0, exp(v/2))
# A fixed step size h must satisfy h*sqrt(curvature) < 2: curvature in theta_i is exp(-v), so
print(np.exp(-v).min(), np.exp(-v).max())        # curvature spans many orders of magnitude

With $v\sim \mathcal{N}(0,9)$ the curvature $e^{-v}$ varies over more than ten orders of magnitude across draws, which is why one step size cannot serve the whole funnel.

Examples (PyMC-style)

Centred (funnel-prone) vs non-centred parameterisation in PyMC 5 (untested sketch):

import pymc as pm
with pm.Model() as centred:
    mu = pm.Normal("mu", 0, 5); tau = pm.HalfCauchy("tau", 2.5)
    theta = pm.Normal("theta", mu, tau, shape=J)            # funnel: divergences when tau -> 0
    pm.Normal("y", theta, sigma, observed=y)
with pm.Model() as noncentred:
    mu = pm.Normal("mu", 0, 5); tau = pm.HalfCauchy("tau", 2.5)
    z = pm.Normal("z", 0, 1, shape=J); theta = pm.Deterministic("theta", mu + tau*z)
    pm.Normal("y", theta, sigma, observed=y)

Connections

Hierarchical Models and Hierarchical Linear Models introduce partial pooling; Betancourt 2015 - Gibbs, Random-Walk Metropolis and Reparameterisation explains why the funnel defeats naive samplers; Stochastic Volatility and the Inference Frontier shows the same geometry in a time-series state-space model.

See Also

Betancourt 2015 - HMC for Hierarchical Models - Overview, The Typical Set and the Log Posterior Density