1. Normal approximation around a posterior mode

Assume we have already:

  • Found (or approximated) a posterior mode $\hat{\theta}$
  • Maybe after using a boundary-avoiding prior as in the previous section.

We want to approximate the posterior density $p(\theta \mid y)$ with a multivariate normal centered at the mode.

Construction

Define:

  • $L(\theta) = \log p(\theta \mid y)$
  • Compute the Hessian (matrix of second derivatives) at the mode:
    • $H = \left.\frac{d^2 \log p(\theta \mid y)}{d\theta^2}\right|_{\theta=\hat{\theta}}$
  • Then the covariance of the normal approximation is:
    • $V_\theta = (-H)^{-1}$ (assuming the posterior is locally concave around $\hat{\theta}$, so $-H$ is positive definite).

So the normal approximation is:

$p_{\text{normal approx}}(\theta) = N(\theta \mid \hat{\theta}, V_\theta)$

To compute $H$ and $V_\theta$​:

  • You can derive the second derivative analytically, or
  • Approximate numerically with finite differences.

Parameter transformation

Before applying the normal approximation, it’s often a good idea to:

  • Transform parameters to be unconstrained and roughly symmetric:
    • positive → log,
    • (0,1) → logit, etc.
  • Then:
    • Work in the transformed space,
    • Add the Jacobian term to the log posterior,
    • Approximate that transformed posterior by a normal distribution.

This usually produces a better Gaussian approximation.

2. Laplace’s method for approximating expectations

Instead of approximating just $p(\theta \mid y)$ itself, we often want:

$E[h(\theta) \mid y]$

for some smooth function $h(\theta)$.

Laplace’s method approximates integrals of the form:

$\int h(\theta)\,p(\theta \mid y)\, d\theta$

Define:

  • $u(\theta) = \log\{h(\theta) p(\theta \mid y)\}$
  • Let $\theta_0$​ be the maximizer of $u(\theta)$:
    • $\theta_0 = \arg\max_\theta u(\theta)$
  • Let $d$ be the dimension of θ.

Then Laplace’s approximation yields:

$\int h(\theta)\, p(\theta \mid y)\, d\theta \approx h(\theta_0)p(\theta_0 \mid y)\,(2\pi)^{d/2} \left| -u”(\theta_0) \right|^{-1/2}$

and, after dividing by$\int p(\theta \mid y) d\theta$, we get an approximation to:

$E[h(\theta) \mid y]$

When is this reasonable?

  • When $p(\theta \mid y)$ is approximately normal around its main mass (large sample size, central limit behavior).
  • When $h(\theta)$ is smooth and doesn’t change too wildly near the high-density region.

Best use cases:

  • Unimodal posteriors,
  • Or applied separately to each mode of a multimodal posterior.

They mention: Laplace’s method is later used to compute relative masses of modes in a normal mixture approximation.

Laplace with unnormalized densities

Often you only know the posterior up to a constant:

  • You can compute $q(\theta \mid y) \propto p(\theta \mid y)$ (unnormalized).

We can still approximate:

$E[h(\theta)\mid y] = \frac{\int h(\theta) q(\theta \mid y)\, d\theta}{ \int q(\theta \mid y)\, d\theta }$

Apply Laplace’s method:

  • Once to the numerator $\int h(\theta) q(\theta \mid y)d\theta$,
  • Once to the denominator$\int q(\theta \mid y)d\theta$,

Then take the ratio of the two approximations. No need to know the normalizing constant.

3. Mixture of normals for multimodal posterior densities

If the posterior has multiple modes, say $K$ modes, we can approximate $p(\theta \mid y)$ as a mixture of K multivariate normals.

Assume we have:

  • Modes: $\hat{\theta}_1, \dots, \hat{\theta}_K$
  • Local covariance matrices: $V_{\theta 1}, \dots, V_{\theta K}$​ derived from Hessians at each mode.

We want:

$p_{\text{normal approx}}(\theta) \propto \sum_{k=1}^K \theta_k\, N(\theta \mid \hat{\theta}_k, V_{\theta k})$

where $\theta_k$​ are the relative weights (masses) of each component.

How to choose the weights $\theta_k$

Let $q(\theta \mid y)$ be the unnormalized posterior density. At each mode:

$q(\hat{\theta}_k \mid y)$

For a multivariate normal density:

  • The height at its mode is proportional to $|V_{\theta k}|^{-1/2}.$

So if we match the approximation to the true (unnormalized) posterior at the mode, we get:

$\hat{\theta}_k = q(\hat{\theta}_k \mid y) \, |V_{\theta k}|^{1/2}$

This leads to a mixture approximation:

$p_{\text{normal approx}}(\theta) \propto \sum_{k=1}^K q(\hat{\theta}_k \mid y)\, \exp\left\{ -\frac{1}{2} (\theta – \hat{\theta}_k)^T V_{\theta k}^{-1} (\theta – \hat{\theta}_k) \right\}$

After normalization, this defines a K-component Gaussian mixture approximating the posterior.

Assumptions:

  • Modes are well separated, so local normal approximations don’t overlap too badly.
  • Each mode is “locally normal”.

4. Using multivariate t instead of normal

Sometimes the posterior has heavier tails than a normal distribution or is not well captured by a Gaussian locally.

In that case, we can replace each normal component with a multivariate t distribution with small degrees of freedom ν:

$p_{t\text{ approx}}(\theta) \propto \sum_{k=1}^K q(\hat{\theta}_k \mid y)\, \left[ \nu + (\theta – \hat{\theta}_k)^T V_{\theta k}^{-1} (\theta – \hat{\theta}_k) \right]^{-(d+\nu)/2}$

where:

  • $d$ = dimension of θ,
  • ν is small (e.g. ν = 4), giving a distribution with heavier tails.

Properties:

  • t-mixture is more robust in tails,
  • Matches posterior better when it’s broader or heavy-tailed.

5. Improving approximations beyond local modes

They briefly mention other refinements:

  • Fit the approximation not only at modes, but also at:
    • saddle points,
    • tail regions,
  • Partially integrate out some components analytically or numerically,
  • Use more advanced approximate methods:
    • Variational Bayes,
    • Expectation Propagation (EP),

which are introduced later.

6. Sampling from the approximate posterior

Once we have a mixture approximation (normal or t), it’s easy to sample from it.

To draw a sample θ from the mixture approximation:

  1. Select a component k:
    • Draw k from {1,…,K} with probabilities proportional to $\theta_k$​ (or $\hat{\theta}_k$​ after normalization).
  2. Sample from that component:
    • If normal: $\theta \sim N(\hat{\theta}_k, V_{\theta k})$
    • If t: $\theta \sim t_\nu(\hat{\theta}_k, V_{\theta k})$

You can implement this with:

  • Multinomial random draw for the component index,
  • Then standard methods (e.g., via Cholesky factorization) to sample from the chosen multivariate normal or t.

Appendix A (of the book) describes:

  • How to sample from a multivariate normal or t using the Cholesky decomposition of the covariance / scale matrix.

7. Using the approximate samples for importance sampling / resampling

Once you can sample from the approximate posterior, you can use these samples to improve inference for the true posterior.

Two main techniques:

  1. Importance sampling
    • Treat the mixture approximation as the proposal distribution, say $g(\theta)$.
    • For each draw $\theta^{(s)} \sim g(\cdot)$, compute an importance weight:
      • $w^{(s)} \propto \frac{p(\theta^{(s)} \mid y)}{g(\theta^{(s)})}$
    • Approximate expectations:
      • $E[h(\theta) \mid y] \approx \frac{\sum_s w^{(s)} h(\theta^{(s)})}{\sum_s w^{(s)}}$
  2. Importance resampling (a.k.a. sampling-importance-resampling, SIR)
    • Draw many candidates from the approximation,
    • Weight them by importance weights,
    • Then resample according to these weights to get a set of approximately unweighted draws from the true posterior.

This can give better approximations than just using the mixture directly, especially if the mixture is only approximate.