CPNS LabUniversity of Exeter

CPNS / An intuitive introduction

Bayesian and
variational inference

How do we learn from uncertain data? Follow a prediction error from a simple Bayesian update to fitting a function—and then fitting a distribution.

Start with BayesTry the first experiment ↓

Three experiments · equations explained in words · no specialist background assumed

01 / Start with uncertainty

Bayes updates a distribution

A measurement does not tell us exactly what is true. It tells us which possibilities have become more plausible.

Suppose we want to infer an unknown quantity θ\theta. Before seeing the data, we describe what we believe with a prior, p(θ)p(\theta). A model describes how different values of θ\theta could generate data yy: this is the likelihood, p(y∣θ)p(y\mid\theta). Bayes’ rule combines them:

p(θ∣y)=p(y∣θ)p(θ)p(y).p(\theta\mid y)=\frac{p(y\mid\theta)\,p(\theta)}{p(y)}.

The posterior is our updated distribution over θ\theta. The evidence, p(y)p(y), normalises it so that its total probability is one. Evidence also matters when comparing models: it measures how well a model predicts the data after averaging over its prior uncertainty.

Prior

What was plausible before this measurement?

Likelihood

How compatible is this measurement with each possible parameter value?

Posterior

What is plausible after combining the two?

A likelihood is not automatically a probability distribution over parameters. It describes the probability density of the data given parameters. Multiplying it by a prior, then normalising, gives the posterior.

Why can the evidence be difficult to calculate?

We need to average over all parameter values:

p(y)=∫p(y∣θ)p(θ)dθ.p(y)=\int p(y\mid\theta)\,p(\theta)\,d\theta.

For some simple models this integral is available analytically. For a nonlinear model with many parameters, it can be expensive or intractable. This is one reason to approximate the posterior.

02 / Read the expression

What does e⊤Πee^\top\Pi e mean?

It is a prediction-error score that takes uncertainty into account.

Let e=y−f(θ)e=y-f(\theta) be the difference between the observed data and the model’s prediction. Write precision as Π\Pi (capital pi), rather than YY. Precision is inverse variance for a scalar, or inverse covariance for a vector:

π=1σ2,Π=Σ−1.\pi=\frac{1}{\sigma^2},\qquad \Pi=\Sigma^{-1}.

Small variance means high precision. Under the model, a precise measurement is expected to be close to its prediction, so the same error counts more strongly. Low precision makes a given mismatch less surprising.

For a single error, the expression is just πe2\pi e^2. If e=2e=2, a standard deviation of 11 gives a score of 44. A standard deviation of 44 gives a score of 0.250.25. The error has not changed; its meaning has changed because the expected uncertainty is different.

For a vector, e⊤e^\top is the same vector written as a row. Multiplication produces one scalar score. For independent errors, this reduces to a weighted sum:

e⊤Πe=∑iπiei2when Π is diagonal.e^\top\Pi e=\sum_i\pi_i e_i^2\qquad\text{when }\Pi\text{ is diagonal}.

Three expressions, three roles

  • ee: the prediction error.
  • Πe\Pi e: the precision-weighted error, still a vector.
  • e⊤Πee^\top\Pi e: the total weighted squared error, a scalar.

If errors are correlated, off-diagonal terms matter too. The score then asks whether the pattern of errors is surprising under the covariance model. Precision is an assumption about reliability; assigning high precision does not make a measurement trustworthy.

Where does the Gaussian distribution enter?

With y∣θ∼𝒩(f(θ),Σy)y\mid\theta\sim\mathcal N(f(\theta),\Sigma_y), the negative log likelihood is

−logp(y∣θ)=12e⊤Πye+12log|Σy|+n2log(2π).-\log p(y\mid\theta)=\frac12e^\top\Pi_y e+\frac12\log|\Sigma_y|+\frac n2\log(2\pi).

If the covariance is fixed, the last two terms do not affect the parameter gradient. If we estimate the noise covariance too, we must retain them: reducing precision cannot be treated as a cost-free way to make errors disappear.

03 / An exact Bayesian update

Move towards evidence, by the right amount

Suppose our prior is θ∼𝒩(μ0,σ02)\theta\sim\mathcal N(\mu_0,\sigma_0^2), and the measurement is y=θ+ϵy=\theta+\epsilon, where ϵ∼𝒩(0,σy2)\epsilon\sim\mathcal N(0,\sigma_y^2). The posterior is also Gaussian. Its precision is the sum of the two precisions, and its mean is their weighted average:

πpost=π0+πy,μpost=π0μ0+πyyπ0+πy.\pi_{\mathrm{post}}=\pi_0+\pi_y,\qquad \mu_{\mathrm{post}}=\frac{\pi_0\mu_0+\pi_y y}{\pi_0+\pi_y}.

The same result can be written as a correction to the old mean:

μpost=μ0+πyπ0+πy⏟gain(y−μ0)⏟prediction error.\mu_{\mathrm{post}}=\mu_0+\underbrace{\frac{\pi_y}{\pi_0+\pi_y}}_{\text{gain}}\underbrace{(y-\mu_0)}_{\text{prediction error}}.

This is the connection: a Bayesian update uses a prediction error, but the amount we move depends on how precise the new evidence is relative to our existing belief.

Experiment 1 / Exact inference

How much should one measurement change your mind?

Move the sliders. The posterior is calculated analytically, with no fitting algorithm.

Prior, likelihood and posterior Probability density curves update when the controls change.
Prior · dashedLikelihood · dottedPosterior · solid
Posterior mean13.20
Posterior SD0.89
Update gain80%

The observation moves the mean 80% of the way from 10 to 14. Posterior uncertainty is smaller than either starting uncertainty.

With the default values, the prior variance is 44 and the measurement variance is 11. The gain is 0.80.8: the mean moves from 1010 to 13.213.2, and posterior variance is 0.80.8. Make the measurement noisier and the correction shrinks. This direct-observation Gaussian example is exact; more complicated models need more work.

The likelihood curve is shown as a density over the horizontal axis for comparison. In this particular direct-observation model it also integrates to one over θ; that does not hold for likelihoods in general.

04 / From a score to a learning rule

Differentiate the mismatch

The quadratic error is a cost. Its derivative tells us how to change a parameter to reduce that cost.

Start with the Gaussian data cost ℰy=12e⊤Πye\mathcal E_y=\frac12e^\top\Pi_y e. Define the Jacobian J=∂f/∂θJ=\partial f/\partial\theta: it tells us how each prediction changes when we change each parameter. Because e=y−f(θ)e=y-f(\theta), increasing a prediction decreases its error. The chain rule gives

∇θℰy=−J⊤Πye.\nabla_\theta\mathcal E_y=-J^\top\Pi_y e.

Add a Gaussian prior with mean θ0\theta_0 and precision Π0\Pi_0. Ignoring constants that are fixed with respect to the parameters, the negative log posterior is

ℰ(θ)=12e⊤Πye+12(θ−θ0)⊤Π0(θ−θ0).\mathcal E(\theta)=\frac12e^\top\Pi_y e+\frac12(\theta-\theta_0)^\top\Pi_0(\theta-\theta_0).

A gradient-descent step is therefore

θnew=θ+η[J⊤Πye−Π0(θ−θ0)].\theta_{\mathrm{new}}=\theta+\eta\left[J^\top\Pi_y e-\Pi_0(\theta-\theta_0)\right].

The first term moves parameters to reduce the data mismatch. The second pulls them towards the prior. The step size η\eta controls the numerical optimisation: it is not, in general, the exact Bayesian gain from the previous example.

Why the factor of one half? Differentiating a square gives a factor of two. The half cancels it, leaving the clean precision-weighted error. The square and the update are related, but they are not the same expression.

Repeated updates seek the maximum a posteriori (MAP) estimate: the parameter value at the posterior’s highest point. We have used Bayesian assumptions to construct the objective, but tracking one value does not yet describe the whole posterior distribution.

05 / Fit a function

Would this reinvent variational Laplace?

It recovers part of the machinery. To reach variational Laplace, we also need an approximation to posterior uncertainty and a variational objective.

Imagine fitting a decaying response, f(t;A,k)=Aexp(−kt)f(t;A,k)=A\exp(-kt). The amplitude AA and decay rate kk are unknown. Each observation creates an error, and the Jacobian maps those errors back to parameter changes.

Curvature can make the updates more efficient. A Gauss–Newton approximation to the curvature of the negative log posterior is

HGN=J⊤ΠyJ+Π0.H_{\mathrm{GN}}=J^\top\Pi_yJ+\Pi_0.

Instead of choosing one scalar step size, it rescales and couples the parameter corrections:

Δθ=HGN−1[J⊤Πye−Π0(θ−θ0)].\Delta\theta=H_{\mathrm{GN}}^{-1}\left[J^\top\Pi_y e-\Pi_0(\theta-\theta_0)\right].

This uses the same information that helps describe local uncertainty: steep curvature implies a narrow range of plausible parameters; shallow curvature implies a wider range.

Experiment 2 / MAP + local curvature

Fit an exponential response

Ten fixed synthetic observations. Take a step, or run the fit. Compare the prediction with the local uncertainty ellipse.

Observed response and fitted prediction Dots show synthetic observations; the solid curve is the current model prediction.
● ObservationsPrediction · solidGenerating curve · dashed
Local uncertainty in amplitude and decay rate A local Gaussian 95 percent contour around the current parameter estimate.

95% contour of the local Gaussian. A tilted ellipse means the parameter uncertainties are coupled.

Amplitude A1.000
Decay rate k0.200
Objective 𝓔Initial parameters

Changing noise or prior precision restarts the fit. Priors: A = 1.2 ± 0.8, k = 0.25 ± 0.35 (mean ± SD, before the multiplier). The generating parameters are A = 2.5 and k = 0.65 s⁻¹. Steps use backtracking to reduce the objective.

Try this: run the fit, then increase the assumed noise and run it again. The observations have not changed, but they now exert less influence. Increase prior precision and the solution is pulled more strongly towards the prior. Switch to gradient descent to see how curvature affects the route to the solution.

This demonstration performs MAP optimisation and displays the inverse Gauss–Newton curvature. Before convergence, the ellipse is only a local curvature summary; near a valid optimum it approximates posterior uncertainty. It is not a full variational Laplace implementation.

06 / Fit a distribution

Variational inference asks a bigger question

Instead of asking only “which parameter value fits?”, ask “which tractable distribution best approximates the posterior?”

Choose a family of distributions q(θ)q(\theta), perhaps Gaussians described by a mean mm and covariance SS. We then optimise those distribution parameters. One common objective is the evidence lower bound (ELBO), also called variational free energy under the convention used here:

F(q)=𝔼q[logp(y∣θ)]−DKL(q(θ)∥p(θ)).F(q)=\mathbb E_q[\log p(y\mid\theta)]-D_{\mathrm{KL}}\!\left(q(\theta)\,\|\,p(\theta)\right).

The first term is expected log likelihood: how well the data are explained across the uncertainty in qq. The second is the divergence from the prior: how much the distribution changes its beliefs to achieve that fit. It is often called complexity, but it is not simply the number of parameters.

For a nonlinear model, optimising the mean of qq need not give exactly the MAP estimate. The objective considers a spread of parameter values, rather than evaluating the fit at just one point.

The key identity is

logp(y)=F(q)+DKL(q(θ)∥p(θ∣y)).\log p(y)=F(q)+D_{\mathrm{KL}}\!\left(q(\theta)\,\|\,p(\theta\mid y)\right).

KL divergence is non-negative. So FF is a lower bound on log evidence, and maximising FF brings qq closer to the posterior in this particular direction of KL divergence. Some texts define free energy with the opposite sign and minimise it; always check the convention.

Experiment 3 / A variational family

Make a Gaussian approximate the posterior

This uses the prior and observation from Experiment 1. Adjust the mean and width of q. Its score rewards both the fit and the right uncertainty.

Here the exact posterior is known, so we can show the gap. For difficult models, that exact comparison is usually unavailable.

Gaussian approximation and exact posterior The dashed curve is the adjustable approximation; the solid curve is the exact posterior.
Exact posterior · solidApproximation q · dashed
Expected log likelihood
KL from prior
ELBO F
Log evidence
Gap to posterior

Try this: match the exact posterior, then make qq very narrow. Concentrating on the best-looking parameter value does not keep improving the bound. We lose the posterior’s uncertainty, and the divergence grows. A distribution’s spread is part of what inference must get right.

What is being calculated in this experiment?

For q=𝒩(m,s2)q=\mathcal N(m,s^2), the expected log likelihood includes both the squared error of its mean and its variance:

𝔼q[logp(y∣θ)]=−12log(2πσy2)−(y−m)2+s22σy2.\mathbb E_q[\log p(y\mid\theta)]=-\frac12\log(2\pi\sigma_y^2)-\frac{(y-m)^2+s^2}{2\sigma_y^2}.

The prior divergence is

DKL(q∥p)=12[s2+(m−μ0)2σ02−1+logσ02s2].D_{\mathrm{KL}}(q\|p)=\frac12\left[\frac{s^2+(m-\mu_0)^2}{\sigma_0^2}-1+\log\frac{\sigma_0^2}{s^2}\right].

The exact evidence is p(y)=𝒩(y;μ0,σ02+σy2)p(y)=\mathcal N(y;\mu_0,\sigma_0^2+\sigma_y^2). All scores use natural logarithms (nats). “Match exact posterior” sets the known analytical solution; it does not simulate an optimisation algorithm.

07 / Put the pieces together

From Laplace approximation to variational Laplace

Near a well-behaved posterior mode, approximate the negative log posterior by a quadratic:

ℰ(θ)≈ℰ(θ̂)+12(θ−θ̂)⊤H(θ−θ̂).\mathcal E(\theta)\approx\mathcal E(\hat\theta)+\frac12(\theta-\hat\theta)^\top H(\theta-\hat\theta).

Exponentiating a negative quadratic gives a Gaussian. This is the Laplace approximation:

q(θ)≈𝒩(θ̂,H−1).q(\theta)\approx\mathcal N(\hat\theta,H^{-1}).

Here HH is the positive-definite Hessian at the mode. In nonlinear least squares, the Gauss–Newton matrix is a useful approximation to that Hessian; it omits terms involving residuals and second derivatives of the predictions.

Variational Laplace uses Gaussian approximations and local curvature within a variational free-energy scheme. In implementations such as SPM’s model inversion, mean updates, posterior covariance, and noise-precision estimates work together. The details go beyond taking a MAP step and appending an error bar.

Routine What it gives you
Minimise squared prediction errors A least-squares parameter estimate; with fixed Gaussian noise, a maximum-likelihood estimate.
Add a prior and optimise A MAP estimate: one best parameter vector under the posterior density.
Approximate curvature at a valid mode A local Gaussian posterior approximation: MAP plus approximate covariance.
Optimise a distribution using a variational bound Variational inference; the approximation can use many different families.
Use a Gaussian/local-curvature scheme for variational inversion Variational Laplace; typically with coordinated mean, covariance and precision updates.

So, would a function-fitting routine be reinventing variational Laplace? It could be heading there. Precision-weighted errors, derivatives and curvature are shared ingredients. The defining extra step is treating inference as an approximation to a distribution, with uncertainty included in the variational objective.

Where the approximation can fail

A single Gaussian can miss multiple posterior modes, asymmetry or curved dependencies. A local optimum may not be the global one. The covariance is conditional on the model and noise assumptions, and an approximate evidence is not exact evidence. Check whether the approximation is adequate for the scientific question.

08 / Bring it back to a generative model

The same logic scales to neural models

In a neural generative model, θ\theta might contain coupling strengths, time constants or other physiological parameters. A forward model predicts a signal from those parameters. Comparing the predicted and observed signals creates a vector of errors. The Jacobian connects changes in the signal to changes in the underlying parameters.

The logic remains: specify plausible parameters and a noise model; predict the observations; use precision to interpret the errors; update the parameter beliefs; and describe how much uncertainty remains. A good fit alone does not establish that a parameter is well identified or that the model is the right explanation.

Prediction

What signal should these parameters generate?

Update

Which parameter changes explain the mismatch?

Uncertainty

Which combinations remain plausible after seeing the data?

A small notation guide

Symbol Meaning
yy, f(θ)f(\theta), ee Observed data, model prediction, and their difference.
Σ\Sigma, Π\Pi Covariance and its inverse, precision. Lowercase π\pi is scalar precision.
JJ, HH Prediction Jacobian and local curvature of the negative log posterior.
θ0\theta_0, Π0\Pi_0 Prior mean and prior precision.
qq, mm, SS Approximate posterior, its mean and its covariance.
FF, DKLD_{\mathrm{KL}} Variational lower bound and Kullback–Leibler divergence.

Follow the ideas further

CPNS student guidePlace inference within the wider modelling workflow. Interactive conceptsExplore related ideas in the CPNS teaching collection. CPNS methodsConnect the introduction to the lab’s methodological work.

Primary references

The examples are deliberately small. The first and third have exact Gaussian answers; the nonlinear fit uses an approximation. Keeping those cases separate makes the connection between Bayes, optimisation and variational inference easier to see.