Four Faces of Leibniz

Starting identity

First Leibniz — pass $\frac{\partial}{\partial\theta}$ through $\int$

Second Leibniz — pass $\frac{\partial}{\partial\theta}$ through $\int$ again

Starting from $\int \frac{\partial \log f}{\partial\theta} \cdot f\,dx = 0$, differentiate the integrand via the product rule:

Integrate both sides (passing $\frac{\partial}{\partial\theta}$ through $\int$):

The definition is a variance — always non-negative, requires only $\mathbb{E}[U]=0$. The second form is a theorem: it requires the second interchange to be valid. When it fails — $\text{Uniform}(0,\theta)$, $\text{Laplace}(\theta)$ — the two expressions disagree, and the definition governs.

Starting identity

where $A(\eta) = \log \int h(x)\exp\{\eta^\top T(x)\}\,dx$ is the log-partition function. (For discrete $X$, replace $\int\cdots\,dx$ with $\sum_x$ throughout; the parameter $\eta$ is always continuous, so all derivatives with respect to $\eta$ proceed unchanged.)

First Leibniz — pass $\frac{\partial}{\partial\eta_j}$ through $\int$

Left side: $A'_j(\eta)\,e^{A(\eta)}$. Right side: $e^{A(\eta)}\,\mathbb{E}_\eta[T_j(X)]$. Cancel $e^{A(\eta)}$:

Second Leibniz — pass $\frac{\partial}{\partial\eta_k}$ through $\int$ again

Differentiate $\frac{\partial A}{\partial\eta_j} = \mathbb{E}_\eta[T_j] = \int T_j(x)\,f(x|\eta)\,dx$ with respect to $\eta_k$. By Leibniz:

We need $\frac{\partial f}{\partial\eta_k}$. Since $f(x|\eta) = h(x)\exp\{\eta^\top T(x) - A(\eta)\}$, the chain rule gives:

where the last step uses the first Leibniz result $\frac{\partial A}{\partial\eta_k} = \mathbb{E}_\eta[T_k]$. Substituting back:

$A$ is convex. Fisher information in natural parameters: $I(\eta) = A''(\eta) = \mathrm{Var}_\eta(T)$. Moments of any order by further differentiation.

Starting identity

Here $f(x)$ is a reward that does not depend on $\theta$, and $p(x|\theta)$ is a policy or sampling distribution parameterized by $\theta$. The goal is to find $\theta$ that maximizes $V(\theta)$ — we want the policy that produces, on average, the highest reward. The natural approach is gradient ascent: $\theta \leftarrow \theta + \alpha\,\nabla_\theta V(\theta)$. So we need $\nabla_\theta V(\theta)$.

The problem — gradient of an expectation you can only sample from

We want $\nabla_\theta V(\theta)$, but the expectation is taken over $p(x|\theta)$ — the distribution we're differentiating with respect to. We can draw samples from $p(x|\theta)$, but we cannot differentiate through the sampling process itself. If $x$ is discrete, or if $p$ is a black-box simulator, a computer looking at the sampling procedure is completely blind — there is no path to push $\nabla_\theta$ inside naively. The Leibniz interchange is what makes this computable at all: without it, there is no gradient signal, and gradient ascent cannot even start.

The solution — same idea as Chapter 2

Recall the fundamental Monte Carlo move from Chapter 2.1: to estimate an integral, rewrite it as an expectation under a distribution you can sample from, then average over samples. Algorithm 2.1 computes $\hat{I}_n = \frac{1}{n}\sum_i h(X_i)$ where $X_i \sim f$. REINFORCE does the same thing for a gradient: it rewrites $\nabla_\theta V$ as $\mathbb{E}_{p(x|\theta)}[\,\cdot\,]$, giving a Monte Carlo estimator

Draw samples, evaluate $f(x_i)\,\nabla_\theta \log p(x_i|\theta)$ at each one, average. The only new ingredient is the Leibniz step that makes the rewrite possible — the estimation itself is Algorithm 2.1 applied to a vector-valued integrand.

Leibniz — pass $\nabla_\theta$ through $\int$

Apply the same identity as Case 1: $\nabla_\theta p = p\,\nabla_\theta \log p$:

Variance reduction via Case 1

The estimator above is unbiased, but just as Chapter 2 warned that naive Monte Carlo has high variance for rare events, the REINFORCE estimator can have enormous variance — $f(x)$ can be large and $\nabla_\theta \log p(x|\theta)$ fluctuates across samples. The baseline trick exploits Case 1 to reduce this variance without introducing bias.

Since $\mathbb{E}_{p(x|\theta)}[\nabla_\theta \log p(x|\theta)] = 0$ (from Case 1), for any constant $b$ that does not depend on $x$:

So we can subtract $b$ from $f(x)$ without changing the expected gradient:

This helps because the variance of the estimator scales with $\mathbb{E}[(f(x) - b)^2 \|\nabla_\theta \log p(x|\theta)\|^2]$. By centering $f(x)$ around $b$, we shrink the magnitude of the terms being averaged. The optimal scalar baseline minimizes this variance:

In practice, $b^*$ is estimated from samples, and simpler choices work well: the running mean of rewards in RL, or a learned value function $V_\phi(s)$ in actor-critic methods. The key point is that any $b$ independent of $x$ is valid — the mean-zero property of the score guarantees unbiasedness.

The problem — flexible models with intractable normalizing constants

Suppose you want to model a complex distribution by specifying how “unlikely” each configuration $x$ is through an energy function $E_\theta(x)$ — a scalar-valued function where low energy means high probability and high energy means low probability. The model takes the Boltzmann form $p_\theta(x) = \exp(-E_\theta(x))\,/\,Z(\theta)$, where the normalizing constant $Z(\theta) = \int \exp(-E_\theta(x))\,dx$ ensures the density integrates to 1. The energy can be as flexible as you like — a neural network, a pairwise interaction model, any function that captures the structure you believe is in the data.

The catch: MLE requires evaluating the likelihood, which requires the normalizing constant $Z(\theta) = \int \exp(-E_\theta(x))\,dx$. For any model rich enough to be interesting, this integral is intractable — it is a sum or integral over all possible configurations of $x$. This is not a rare corner case. It arises in:

• Undirected graphical models (Markov random fields, Boltzmann machines) — the partition function sums over exponentially many configurations
• Spatial and network models — autologistic models for lattice data, exponential random graph models (ERGMs) for network data
• Energy-based models in machine learning — any model of the form $p_\theta(x) \propto \exp(-E_\theta(x))$ with a flexible $E_\theta$

Without $Z(\theta)$, you cannot compute the likelihood, so you cannot do MLE, you cannot do Bayesian updates, you cannot compare models by AIC or BIC. Score matching offers an escape: fit the model without ever evaluating $Z(\theta)$.

The key observation — $\nabla_x$ kills the normalizing constant

We want to fit a model $p_\theta(x)$ to data from unknown $p_{\text{data}}(x)$. The data-side score is $\mathbf{s}(x) = \nabla_x \log p(x)$ — gradient w.r.t. $x$, not $\theta$. Unlike the classical score from Case 1 (a scalar measuring sensitivity to the parameter), this is a vector field over data space: at every point $x$, it points in the direction where the density increases fastest. Think of the density $p(x)$ as a landscape — the data-side score is the slope at every point, telling you which way is “uphill” toward higher probability. Score matching fits this vector field rather than the probability surface itself, and this is precisely what lets the normalizing constant drop out.

If $p_\theta(x) = q_\theta(x)/Z(\theta)$, then:

$Z(\theta)$ does not depend on $x$, so its $x$-gradient vanishes.

The Fisher divergence

Expanding the square: the first term is constant in $\theta$, the third is computable. The cross term requires the unknown $\nabla_x \log p_{\text{data}}$.

Integration by parts in $x$ — the key move

For each component $i$:

Integrate by parts in $x_i$ (assuming $p_{\text{data}}(x)\,\frac{\partial \log p_\theta(x)}{\partial x_i} \to 0$ as $\|x\| \to \infty$, so the boundary term vanishes):

The right side involves only derivatives of $\log p_\theta$ evaluated at data points — no $p_{\text{data}}$ needed, no normalizing constant needed.

Where is Leibniz? Integration by parts is Leibniz

In Cases 1–3, we used $\int \frac{\partial}{\partial\theta}[\cdots]\,dx = \frac{\partial}{\partial\theta}\int[\cdots]\,dx$ — moving a $\theta$-derivative past an $x$-integral. Here the derivative is in $x$, the same variable we integrate over, so the tool is integration by parts rather than Leibniz in its standard form. But it is the same structural move. Start from the product rule:

Integrate both sides over all $x$. The left side is a total derivative — it integrates to a boundary term $\big[p_{\text{data}}\cdot\frac{\partial \log p_\theta}{\partial x_i}\big]_{-\infty}^{+\infty}$ that vanishes when $p_{\text{data}}(x)\,\frac{\partial \log p_\theta}{\partial x_i} \to 0$ as $\|x\|\to\infty$:

This is exactly the pattern: an integral equals something known (zero, because the boundary terms vanish), and the constraint transfers a derivative from one factor ($p_{\text{data}}$, which we don't know) to another ($\log p_\theta$, which we do). In Cases 1–3, the known constant was $\int f\,dx = 1$. Here it is $\big[p_{\text{data}}\cdot\frac{\partial \log p_\theta}{\partial x_i}\big]_{-\infty}^{+\infty} = 0$. Different constant, same move.