The Mathematics of Diffusion — Generation Seen Through Scores and SDEs

A layer beneath 'add noise, then subtract it': diffusion models restated in the language of the score, the gradient of log-density. Why denoising and score estimation are literally the same job, what the forward SDE, reverse SDE and probability flow ODE actually assert, and how those equations turn into the knobs you set at inference time.

Generative Modeling by Estimating Gradients of the Data Distribution

Forget the altitude, remember the slope

Suppose you are somewhere in the mountains at night with no map and no GPS, trying to reach the summit. You have no way of knowing your altitude in metres. You can still climb, because knowing the slope under your feet is enough: step toward whichever direction goes up, repeat.

Generative modelling is in exactly that position. Try to learn the probability density $p(x)$ of "plausible images" directly, and a normalising constant always gets in the way. Write the model as an energy function $E_\theta(x)$ and the density looks like this:

$Z_\theta$ is the number you divide by so that everything sums to one, and the integral runs over every image there could be. For a $256\times256$ colour image that is an integral over roughly 200,000 dimensions. Nobody can compute it.

But take the logarithm of both sides and then differentiate with respect to $x$, and $Z_\theta$ vanishes — it does not depend on $x$ at all.

This is the score. The name is borrowed from the score function of classical statistics (the gradient of the log-likelihood), but the diffusion version differentiates with respect to the data $x$ rather than the parameters. In words: $\nabla_x$ means "differentiate along every component of $x$ and stack the results", so the score has exactly the same shape as its input. For an image, picture a field of arrows with one arrow per pixel, each saying nudge me this way and the picture becomes slightly more plausible. You give up on knowing the altitude and settle for knowing the slope.

The overall machinery of diffusion — the forward process that adds noise and the reverse process that removes it — is covered in Diffusion Models from the Ground Up. This article goes one layer beneath it and asks why that procedure is correct at all, in the language of scores.

Turning a slope into samples

Suppose you had the score. The recipe for extracting samples from it is Langevin dynamics, and it is one line:

Here $\varepsilon$ is the step size and $z_k$ is a fresh standard normal draw each iteration. In words: take half a step toward higher density, then stagger by however much the dice say, over and over. Make $\varepsilon$ small enough and run long enough, and the distribution of $x_k$ approaches $p_\theta$.

That second term is not decoration. Remove it and you are left with plain gradient ascent, which parks itself on the single most probable point and stops. That finds the mode; it does not draw from the distribution. The noise is what keeps the walker wandering the slope, lingering long where density is high and briefly where it is low — and that ratio of dwelling times is the probability.

A concrete case makes it tangible. For $p(x) = \mathcal{N}(0, \sigma^2)$ we have $\log p(x) = -x^2/(2\sigma^2) + \text{const}$, so differentiating gives

which says the score is a spring pulling back toward the origin: the further out you are the harder the pull, and the tighter the distribution (small $\sigma$) the stiffer the spring. So Langevin dynamics on a Gaussian is precisely a ball rolling in a parabolic valley with a random shove added at every step. And what happens when you crank the step size up is exactly what happens to that ball.

Why "denoising" and "score estimation" are the same job

There is an obvious problem. Nowhere do we have ground-truth labels for the score. If we knew the true $\nabla_x \log p(x)$, we would not need a generative model in the first place.

Plain score matching (Hyvärinen, 2005) solves this on paper. Integration by parts eliminates the unknown $p$ from the objective, leaving you to minimise the expectation of $\frac{1}{2}\|s_\theta(x)\|^2 + \operatorname{tr}(\nabla_x s_\theta(x))$. The theory is airtight, but that second term is the trace of a Jacobian, which costs one backward pass per dimension. On a 200,000-dimensional image it is a non-starter.

The way through was to deliberately corrupt the data. We do not know the distribution of clean data, but the conditional distribution after adding Gaussian noise is something we built ourselves, so we know it exactly. Corrupt with $\tilde{x} = x + \sigma z$ where $z \sim \mathcal{N}(0,I)$, and $p_\sigma(\tilde{x}\mid x) = \mathcal{N}(\tilde{x}; x, \sigma^2 I)$ — whose score is the spring we just derived: