Stein Variational Gradient Descent: Particles That Repel Their Way to the Posterior

Bayesian inference is simple to state and hard to compute. Given data, we want the posterior distribution over parameters. In almost every interesting model, that posterior has no closed form, its normalising constant is intractable, and direct sampling is expensive. Two classical families of tools attack the problem from opposite directions. Markov chain Monte Carlo (MCMC) gives asymptotically exact samples but can be slow to mix. Variational inference (VI) turns inference into optimisation, but usually forces the posterior into a chosen parametric family, such as a Gaussian, and then hopes that family is flexible enough.

Stein variational gradient descent (SVGD), introduced by Liu and Wang, offers a third way. It keeps a set of particles, moves them with a deterministic gradient-like rule, and lets the particles themselves represent the posterior. There is no parametric family to fit and no accept-or-reject step. The trick that makes it work is a repulsive force between particles. Without it, every particle would collapse to the same mode. With it, the particles spread out and trace the shape of the posterior.

The problem, stated precisely

Let $θ$ be a parameter vector and let $p (θ)$ be the target posterior, known only up to a constant. That is, we can evaluate $lo g p (θ)$ and its gradient, often called the score function, but we cannot evaluate the normalising constant and we cannot sample from $p$ directly.

Suppose we have $n$ particles $θ_{1}, \dots, θ_{n}$ drawn from some initial distribution $q$ . Our goal is to move these particles, step by step, so that their empirical distribution gets closer to $p$ . Closeness here is measured by Kullback-Leibler (KL) divergence:

KL (q ∥ p) = E_{θ \sim q} [lo g q (θ) - lo g p (θ)]

If we apply a smooth transformation $T (θ) = θ + ϵ ϕ (θ)$ to every particle, with a small step size $ϵ$ and a vector field $ϕ$ to be chosen, the distribution of the transformed particles changes. The question is: which $ϕ$ reduces the KL divergence fastest?

The Stein identity and the best direction

The answer comes from Stein's identity. For a smooth test function $ϕ$ that vanishes suitably at infinity, the following expectation under $p$ is zero:

E θ \sim p [\nabla θ lo g p (θ) \cdot ϕ (θ) + \nabla_{θ} \cdot ϕ (θ)] = 0

This identity is the reason SVGD needs only the score function $\nabla_{θ} lo g p (θ)$ , not the normalising constant. The constant disappears when we take the log gradient.

Stein variational gradient descent restricts $ϕ$ to lie in a reproducing kernel Hilbert space (RKHS) defined by a positive definite kernel $k$ . In that space, the direction of steepest descent of the KL divergence has a closed form. It is the expectation, over the current particle distribution, of a kernel-weighted score plus a kernel gradient term:

ϕ^{*} (θ) = E θ^{'} \sim q [k (θ^{'}, θ) \nabla θ^{'} lo g p (θ^{'}) + \nabla_{θ^{'}} k (θ^{'}, θ)]

That single equation contains the whole algorithm. Read its two terms separately, because they do opposite jobs.

The first term, $k (θ^{'}, θ) \nabla_{θ^{'}} lo g p (θ^{'})$ , is an attraction term. It is a kernel-weighted average of the score at the other particles. It pulls a particle uphill, toward regions where the posterior density is high, with nearby particles contributing more weight than distant ones.

The second term, $\nabla_{θ^{'}} k (θ^{'}, θ)$ , is a repulsion term. For the common radial basis function (RBF) kernel, this term points away from neighbouring particles. It stops particles from piling up on top of each other at a single mode. Without it, SVGD would be little more than parallel gradient ascent, and every particle would converge to the nearest mode. With it, the particles behave like a system of mutually repelling charges that still feels the pull of the posterior.

In practice we replace the expectation by the average over the finite particle set. The update for particle $i$ at iteration $t$ is:

θ_{i} \leftarrow θ_{i} + ϵ ϕ^^{(} θ_{i}), ϕ^^{(} θ_{i}) = \frac{1}{n} j = 1 \sum n [k (θ_{j}, θ_{i}) \nabla_{θ_{j}} lo g p (θ_{j}) + \nabla_{θ_{j}} k (θ_{j}, θ_{i})]

That is the entire algorithm. Initialise particles, repeatedly apply this update, and use the final particles as approximate posterior samples.

The kernel matters more than it looks

The usual choice is the RBF kernel:

k (θ, θ^{'}) = exp (- \frac{∣ θ - θ ^{'} ∣ ^{2}}{2 h ^{2}})

The bandwidth $h$ controls the range of interaction. If $h$ is too small, each particle only sees itself, the repulsion term vanishes, and the particles collapse. If $h$ is too large, the kernel is nearly constant, the repulsion is weak and indiscriminate, and the particles move almost in lockstep.

The standard practical fix is the median heuristic: set $h$ to the median pairwise distance between particles, recomputed at each iteration. This adapts the interaction range to the current spread of the particles, which is exactly what you want as the particles move from a diffuse initialisation to a concentrated posterior. It is a heuristic, not a theorem, but it works well in a wide range of problems and it removes one fragile tuning knob.

For the RBF kernel, the gradient term in the update has an explicit form:

\nabla_{θ_{j}} k (θ_{j}, θ_{i}) = - \frac{( θ _{j} - θ _{i} )}{h ^{2}} k (θ_{j}, θ_{i})

Notice the sign. When $θ_{j}$ is close to $θ_{i}$ , this term pushes $θ_{i}$ away from $θ_{j}$ , with a strength that decays with distance. That is the repulsion, written out in full.

A concrete example in code

Consider a simple target: a mixture of two Gaussians in two dimensions, with modes that are well separated. A single Gaussian variational approximation would be forced to pick one mode or sit uselessly between them. SVGD particles, by contrast, should split across both modes. The code below implements the update in NumPy for this toy target. It is short on purpose: the algorithm really is this short once the score function is available.

import numpy as np def log_target_and_score(x): # Mixture of two 2D Gaussians, equal weights.  # Returns log density (up to a constant) and its gradient.  mu1 = np.array([2.0, 2.0]) mu2 = np.array([-2.0, -2.0]) d1 = x - mu1 d2 = x - mu2 log_p1 = -0.5 * np.sum(d1 ** 2, axis=-1) log_p2 = -0.5 * np.sum(d2 ** 2, axis=-1) m = np.maximum(log_p1, log_p2) log_p = m + np.log(0.5 * np.exp(log_p1 - m) + 0.5 * np.exp(log_p2 - m)) w1 = np.exp(log_p1 - log_p) w2 = np.exp(log_p2 - log_p) score = -(w1[:, None] * d1 + w2[:, None] * d2) return log_p, score def svgd(x, n_iter=500, eps=0.1): n = x.shape[0] for _ in range(n_iter): log_p, score = log_target_and_score(x) # Pairwise squared distances and the median heuristic.  diff = x[:, None, :] - x[None,Source：https://dev.to/umeshtharukaofficial/stein-variational-gradient-descent-particles-that-repel-their-way-to-the-posterior-1km8My favorites

PreviousTeaching by Building: How Sunmughan Swamy Turns What He Learns Into Guides Others Can UseNextArgon vs Wonder: Which AI Should Developers Build On?

`Related posts`


The Download: the agentic web arrives and the drive for safer roads
This is today’s edition of The Download , our weekday newsletter that provides a daily dose of what’s going on in the wo
0MIT Technology Review
Stop counting tokens, start counting what ships
When a developer gives an AI agent a coding task, the agent does more than answer once. It plans, generates code, runs t
0MIT Technology Review
OpenAI is finally merging ChatGPT’s chat and Work tabs
OpenAI said that ChatGPT and Codex will draw usage from the same allocation pool after this change.
0TechCrunch
10 Open-Source Robotics Repositories for 3D Perception, VLA Policies, and Sensor Calibration
Physical AI breaks down at the boundary where software assumptions meet sensor reality. Vision-Language-Action (VLA) mod
0dev.to
My self-modifying system tried to wire a gate-passing module into its own protection enforcer. It was refused three times.
Part 4 of a series on a self-modifying multi-agent LLM system. Report and data: doi.org/10.5281/zenodo.23622286 . This e
0dev.to
Last night an AI agent performed brain surgery on a mouse
Article URL: https://twitter.com/OcanaSantero/status/2115442964574908902 Comments URL: https://news.ycombinator.com/item
0Hacker News