From Information Theory to Information Geometry

00

Outline: the path through this tutorial

By the end of this tutorial, you will have followed one continuous path from entropy and KL divergence to Fisher geometry, dual connections, and a natural-gradient variational-inference update. The five stages below preview that path and provide a practical sequence for revisiting the ideas afterward.

Learn it as a spiral: revisit the same few distributions with one new tool each time. Do not wait to “finish differential geometry.” The Bernoulli, categorical, and Gaussian families are enough to make nearly every new idea concrete.

1Rebuild the information layer

Entropy, cross-entropy, KL, JS, and coding. Derive each for a Bernoulli coin.

2Learn only the geometry you need

Coordinates, tangent vectors, metrics, connections, and geodesics—always tied to a distribution family.

3Compute Fisher information

Do Bernoulli first, then categorical and univariate Gaussian. Verify that coordinate changes preserve $ds^2$ .

4Meet exponential families

Natural $\theta$ and expectation $\eta$ coordinates make duality concrete instead of mystical.

5Apply it

Implement a natural-gradient update and compare it with an ordinary parameter gradient.

Sources and further reading

C. E. Shannon (1948), “A Mathematical Theory of Communication.” Bell System Technical Journal.
T. Cover and J. Thomas (2006), Elements of Information Theory, 2nd ed.
R. A. Fisher (1922), “On the Mathematical Foundations of Theoretical Statistics.” Philosophical Transactions of the Royal Society A.
C. R. Rao (1945), “Information and the Accuracy Attainable in the Estimation of Statistical Parameters.” Reprinted in Breakthroughs in Statistics.
S. Amari and H. Nagaoka (2000), Methods of Information Geometry.
F. Nielsen (2020), “An Elementary Introduction to Information Geometry”. Open-access survey.
D. J. C. MacKay (2003), Information Theory, Inference, and Learning Algorithms.

The bridge: distributions become points

Information theory begins with questions like: How surprising was this outcome? and How many bits does this message require? Information geometry asks a new question: What kind of space is formed by the probability distributions within a chosen statistical model?

Here, a statistical model is a particular family of probability distributions with the same mathematical form:

\mathcal M=\{p_\theta:\theta\in\Theta\}

The parameter $\theta$ selects one distribution $p_\theta$ from that family. That individual distribution becomes one point in the geometry, while its parameter values provide coordinates for locating it.

For example, the Bernoulli family contains every Bernoulli distribution obtained by varying $p$ , so its points occupy the interval $0<p<1$ . A three-outcome categorical family contains every probability vector $(p_1,p_2,p_3)$ satisfying $p_1+p_2+p_3=1$ , so its points occupy a triangle.

Ordinary geometryA pointa location

↔

Information geometryA pointan entire probability distribution

The central idea

Two distributions count as “nearby” when observations have difficulty distinguishing them—not merely when their parameter values look numerically close. Fisher information converts that operational notion of distinguishability into geometry.

Here is the whole bridge before we walk it slowly. Each step answers a question created by the previous one; none of the geometric language is being added merely for decoration.

1

Start with prediction

Reality is P; your model is Q

Reality draws an outcome $x$ according to $P(x)$ . Your model does not know that draw in advance; it assigns its own probability $Q(x)$ . The less probability Q assigned to what actually happened, the more surprised Q is—and the more information is needed to identify that outcome.

The surprisal, or pointwise information, of this one outcome is written with a lowercase $h_Q(x)$ . It is not yet the entropy of the whole distribution. Shannon entropy $H(Q)$ is obtained by averaging the surprisal of every possible outcome using Q's own probabilities.

h_Q(x)=-\log_2 Q(x)

H(Q)=\mathbb E_{X\sim Q}[h_Q(X)]=\mathbb E_{X\sim Q}[-\log_2 Q(X)]

Why define information like this?

In 1948, while working at Bell Labs, Claude Shannon formalized information to address the engineering problem of transmitting messages efficiently and reliably through a communication channel. He needed a numerical measure of how much information an observed outcome contributes. That measure had to satisfy three desiderata—three properties expected of any sensible definition:

Rarer outcomes carry more information. Learning that an unlikely event occurred tells us more than learning that an expected event occurred.
Information varies continuously with probability. A small change in probability should produce a small change in information.
Independent outcomes add their information. Learning two independent facts should provide the information from the first plus the information from the second.

The third requirement supplies the mathematical key. If two independent outcomes have probabilities $p$ and $q$ , their joint probability is $pq$ , while their combined information should be $h(pq)=h(p)+h(q)$ . The logarithm is precisely the function that converts multiplication into addition. Together, the three requirements determine $h(p)=-k\log p$ up to a positive scale factor. Choosing base two and $k=1$ measures information in bits.

This gives the same number two complementary meanings. $h_Q(x)$ measures how surprising $x$ is under Q, and it gives the optimal amount of information needed to encode that outcome under Q. More probable outcomes can be represented more economically; improbable outcomes require more information to distinguish them.

One outcome versus the whole random variableSurprisal measures one possible outcome; entropy averages those surprisals using the distribution's probabilities.

For a broader intuitive treatment, see Thomas Varley's Information Theory for Complex Systems Scientists (2025).

2

Average the prediction error

Cross-entropy becomes KL divergence

Averaging Q's code lengths over data from P gives the cross-entropy $H(P;Q)$ . But some bits are unavoidable even with a perfect model: that irreducible cost is $H(P)$ . Subtract it, and what remains is exactly the penalty for using Q when reality is P.

H(P;Q)=\mathbb E_{x\sim P}[-\log_2 Q(x)]

D_{\mathrm{KL}}(P\|Q)=H(P;Q)-H(P)

A concrete coin model

In the coin model introduced above, let P be the distribution that actually generates the flips and Q the model used to predict and encode them. Take $P=(0.5,0.5)$ while Q assigns $Q=(0.75,0.25)$ to heads and tails. Q assigns a probability—and therefore an expected code length—to every outcome P might generate.

Interactive toy model

Encode coin P using a code designed for coin Q

Reality P: heads probability 0.50Model Q: heads probability 0.75Message length: 100 flips

P — the coin generating the dataHeads 50% · Tails 50%

H

T

Q — your model and codebookHeads 75% · Tails 25%

H

T

Best possible1.000 bits/flipH(P): code matched to reality

Using model Q1.208 bits/flipH(P;Q): expected code length

Model penalty0.208 bits/flipDₖₗ(P‖Q): excess bits

Across 100 flips20.75 extra bitsthe mismatch accumulates

If the next result isHeads0.415 ideal bits

-\log_2 Q(H)

If the next result isTails2.000 ideal bits

-\log_2 Q(T)

Your model expects 75 heads per 100 flips, while reality produces about 50. Q rewards its expected outcomes with short codewords, but pays dearly whenever its “surprising” outcome occurs.

What does it mean for Q to model P?

We judge Q by averaging its surprises under P:

\mathbb{E}_{x\sim P}[-\log Q(x)]

. A good Q places probability where P actually produces data.

The price of the wrong codebook

If an outcome has model probability $Q(x)$ , its ideal information content is $-\log_2 Q(x)$ bits. Averaging that quantity over outcomes generated by P gives cross-entropy. Subtracting the best possible average code length $H(P)$ leaves KL divergence.

CROSS-ENTROPY

H(P;Q)=\mathbb{E}_{x\sim P}\!\left[-\log_2 Q(x)\right]

KL DIVERGENCE

D_{\mathrm{KL}}(P\|Q)=H(P;Q)-H(P)=\mathbb{E}_{x\sim P}\!\left[\log_2\frac{P(x)}{Q(x)}\right]

The coin model: 50/50 reality versus a 75/25 model

A Q-code gives heads an ideal length of $-\log_2(0.75)\approx0.415$ bits.
It gives tails a much longer length: $-\log_2(0.25)=2$ bits.
Reality produces each half the time, so $H(P;Q)\approx\tfrac12(0.415)+\tfrac12(2)=1.2075$ bits/flip.
A fair-coin code needs $H(P)=1$ bit/flip.
Therefore $D_{\mathrm{KL}}(P\|Q)\approx0.2075$ excess bits/flip.

Why the order matters

D_{\mathrm{KL}}(P\|Q)

asks how Q performs on data from P. Reversing the arguments asks a different question. If Q declares an event impossible while P can produce it, Q eventually encounters an outcome requiring an infinite ideal code length.

How can one coin-flip outcome require a fractional number of bits?

The value 0.415 bits is not a literal 0.415-bit string attached to one isolated flip. It is that outcome's ideal contribution to the average code length. Across long messages, arithmetic coding and related schemes can approach the expected $-\log_2 Q(x)$ bits per outcome. This is the source-coding interpretation developed from Shannon's theory and treated carefully by Cover and Thomas.

Coding interpretation: Shannon (1948); Cover & Thomas (2006).

Why this forces Step 3: KL compares two separated models, but geometry begins locally. We therefore ask what KL looks like when Q is only an infinitesimal change away from P.

3

Zoom in around one model

The local shape of KL is Fisher information

Begin with a parameterized family $\mathcal M=\{p_\vartheta:\vartheta\in\Theta\}$ —for example, every Bernoulli distribution obtained by changing its probability parameter. Fix one member as the current distribution $P=p_\theta$ . Changing its parameters by $\Delta\theta$ selects another member, $Q=p_{\theta+\Delta\theta}$ . We are asking how the probabilities assigned to all possible outcomes change as a function of that parameter displacement. Call the resulting KL divergence $F(\Delta\theta)=D_{\mathrm{KL}}(p_\theta\|p_{\theta+\Delta\theta})$ .

If the model has several parameters, $\nabla F(0)$ is the gradient: the column vector containing every first partial derivative of $F$ , evaluated at $\Delta\theta=0$ . The superscript ${\mathsf T}$ means transpose, which turns that column into a row. Multiplying it by the parameter-change column $\Delta\theta$ forms the dot product $\nabla F(0)^{\mathsf T}\Delta\theta$ : one number giving the first-order change in $F$ along that particular direction. The transpose has no additional calculus meaning; it is needed only because multidimensional gradients and parameter changes are conventionally written as column vectors. We could equivalently write the dot product as $\nabla F(0)\mathbin{\cdot}\Delta\theta$ . For a one-parameter model, it reduces to the familiar $F'(0)\Delta\theta$ .

Writing KL explicitly gives

F(\Delta\theta)=\mathbb E_\theta\!\left[\log p_\theta(X)-\log p_{\theta+\Delta\theta}(X)\right].

Only the nearby log-density depends on $\Delta\theta$ . Even though the nearby log-density can be evaluated exactly in simple models such as the Gaussian, Taylor-expanding it reveals how KL changes locally as the parameters move. It replaces the generally nonlinear dependence on $\Delta\theta$ with a quadratic approximation built from derivatives at the current model. In more complicated families, this local description can also be less expensive than repeatedly evaluating or optimizing the exact KL. To calculate the Taylor coefficients of $F$ , we therefore expand

\log p_{\theta+\Delta\theta}(X)

and then average its terms under $p_\theta$ . Under the regularity conditions, expanding inside the expectation produces the same Taylor expansion as differentiating the complete KL function directly; the step-by-step derivation later in the tutorial shows this calculation explicitly. At $\Delta\theta=0$ , the distributions match and $F$ is at its minimum value of zero. Its constant and linear Taylor terms therefore vanish, leaving the quadratic term as the first term that can survive.

F(\Delta\theta)=\underbrace{F(0)}_{\text{constant term }=0}+\underbrace{\nabla F(0)^{\mathsf T}\Delta\theta}_{\text{linear term }=0}+\underbrace{\frac12\Delta\theta^{\mathsf T}I(\theta)\Delta\theta}_{\text{quadratic term}}+\underbrace{O(\|\Delta\theta\|^3)}_{\text{cubic-and-higher remainder}}

D_{\mathrm{KL}}^{(\mathrm{nats})}(p_\theta\|p_{\theta+\Delta\theta})=\frac12\Delta\theta^{\mathsf T}I(\theta)\Delta\theta+O(\|\Delta\theta\|^3)

The two symbols in that equation play different roles. The KL divergence is one scalar assigned to the pair $(p_\theta,p_{\theta+\Delta\theta})$ . The Fisher information $I(\theta)$ is the Fisher information of the model family evaluated at the starting point $p_\theta$ . It describes how quickly KL grows when the parameters move away from that point in each possible direction. It is one number for a one-parameter family such as Bernoulli, and a matrix when $\theta$ contains several parameters.

One-dimensional cross-section through model space

The local picture: the bottom of a KL bowl

Scroll horizontally to compare the full curves →

Hold

p_\theta

fixed and treat

F(\Delta\theta)

as the height above each possible parameter displacement

\Delta\theta

. At the center,

\Delta\theta=0

, the models match and KL is zero. Because KL cannot be negative, this is the bottom of the bowl, with no first-order slope. The first visible local shape is therefore curvature:

I(\theta)

tells us whether KL rises gently or steeply in each direction. A steeper direction means that data can distinguish a small model change more easily. The solid curve is exact; the dashed quadratic curve is the close-up supplied by Fisher information.

Aside

What Fisher information was built to answer

Fisher information began as a tool of statistical estimation, before its later geometric interpretation. The geometry did not invent a new quantity. Mathematicians recognized that the same object used to quantify parameter precision also supplies a natural local metric on a family of probability distributions.

Suppose an observable random variable $Y$ is generated from a known family $p(y\mid\theta)$ , but the parameter value $\theta$ is unknown. The estimation problem is to infer that hidden value from measurements $y_1,\ldots,y_n$ . Fisher information asks: how well does an observation of $Y$ locate the parameter value that generated it?It assigns a larger value when nearby parameter choices produce observably different distributions, and a smaller value when many parameter choices would produce nearly indistinguishable data.

Shannon informationUncertainty about an outcome

How surprising or uncertain is the random state $Y$ ?

Fisher informationSensitivity to a parameter

How strongly does the distribution of $Y$ change when $\theta$ changes?

Fisher information is therefore not a bit count and not the mutual information between two random variables. In its classical use, $\theta$ is an unknown fixed value rather than a random variable with a prior. Fisher information is local: $I_Y(\theta)$ describes the measurement model's sensitivity near a particular parameter value. The subscript $Y$ reminds us that a different measurement procedure could carry a different amount of information about the same parameter.

A clean example: known Gaussian noise, unknown mean

Begin with a noise variable $\varepsilon\sim\mathcal N(0,\sigma^2)$ . The symbol $\sim$ means “is distributed according to”: the noise is Gaussian with known mean zero and known variance $\sigma^2$ . Now define the measurement by $Y=\mu+\varepsilon$ , where $\mu$ is the unknown signal. Adding $\mu$ shifts the noise distribution from being centered at zero to being centered at $\mu$ , without changing its variance. Therefore the observation has one Gaussian measurement model, $Y\mid\mu\sim\mathcal N(\mu,\sigma^2)$ ; the first Gaussian described its noise, while this one describes the resulting observation.

FROM ZERO-MEAN NOISE TO THE OBSERVATION MODEL

\varepsilon\sim\mathcal N(0,\sigma^2),\qquad Y=\mu+\varepsilon\qquad\Longrightarrow\qquad Y\mid\mu\sim\mathcal N(\mu,\sigma^2)

Capital $Y$ denotes the measurement before it is observed; lowercase $y_i$ denotes an actual observed value. Conditional on $\mu$ , the unknown mean is treated as one fixed constant; only $\varepsilon$ varies from one possible observation to another. Expectation therefore shifts by $\mu$ , while variance is unchanged because adding a constant moves every possible value equally without changing any value's displacement from the center.

The following derivation makes both relationships explicit, applying the rules of conditional expectation and variance one step at a time.

Y=\mu+\varepsilon,\qquad \mathbb E[\varepsilon]=0,\qquad \operatorname{Var}(\varepsilon)=\sigma^2

Begin with a fixed signal μ plus random zero-centered noise ε.

\mathbb E[Y\mid\mu]=\mathbb E[\mu+\varepsilon\mid\mu]

Substitute the measurement equation Y = μ + ε for Y inside the conditional expectation.

\mathbb E[\mu+\varepsilon\mid\mu]=\mathbb E[\mu\mid\mu]+\mathbb E[\varepsilon\mid\mu]

Use linearity of conditional expectation: the expectation of a sum equals the sum of the expectations.

\mathbb E[\mu\mid\mu]=\mu

Conditioning on μ means that μ is held fixed, and the expectation of a fixed value is that value itself.

\mathbb E[\varepsilon\mid\mu]=\mathbb E[\varepsilon]=0

The model assigns the same zero-mean noise distribution to every fixed μ, so knowing μ does not change the average noise.

\mathbb E[Y\mid\mu]=\mu+0=\mu

Combine the fixed signal μ with the zero expected contribution from the noise.

Y-\mathbb E[Y\mid\mu]=(\mu+\varepsilon)-\mu=\varepsilon

Subtracting the observation mean removes the constant shift and leaves exactly the noise.

\operatorname{Var}(Y\mid\mu)=\mathbb E[(Y-\mu)^2\mid\mu]=\mathbb E[\varepsilon^2]=\sigma^2

Variance averages squared displacement from the mean, so the observation inherits the noise variance.

The noise therefore does not move the distribution's mean; it moves individual observations away from that mean. Small $\sigma$ keeps likely observations near $\mu$ , so an observed value points back to a narrow range of plausible means. Large $\sigma$ permits likely observations much farther from $\mu$ . The same observed value can then be explained by many candidate means combined with different plausible noise values, making the true mean harder to locate.

From a noisy observation to evidence about the mean

The equations above describe the forward generative process by which observations of $Y$ are created: a fixed but unknown generating mean $\mu^\ast$ combines with fresh noise to produce each observation. Estimation reverses this by asking the question: once a particular value $Y=y$ has been observed, which candidate values of $\mu$ could plausibly have generated it? The same expression $p(y\mid\mu)$ becomes a likelihood over candidate means when $y$ is held fixed and $\mu$ is allowed to vary. Differentiating its logarithm with respect to those candidate means produces the score $s_\mu(y)=\partial_\mu\ln p(y\mid\mu)$ , a quantitative directional clue about the generating mean. If $y>\mu$ , the score is positive and recommends increasing the candidate mean; if $y<\mu$ , it is negative and recommends decreasing it. Its magnitude records the strength of that local recommendation. For known $\sigma$ , the likelihood is a Gaussian-shaped curve centered at $\mu=y$ : the observation points most strongly to itself as the possible mean, while the noise level determines how wide a range of other means remains plausible.

Y\mid\mu^\ast\sim\mathcal N(\mu^\ast,\sigma^2)

Forward view: the unknown generating mean μ* is fixed, while repeated observations vary because of noise.

L(\mu;y)=p(y\mid\mu)\propto\exp\!\left[-\frac{(y-\mu)^2}{2\sigma^2}\right]

Inverse view: hold the observed value y fixed and compare how plausible it is under each candidate mean μ.

\ell(\mu;y)=\ln p(y\mid\mu)=-\frac12\ln(2\pi\sigma^2)-\frac{(y-\mu)^2}{2\sigma^2}

Take the logarithm of the likelihood; this preserves its peak while making derivatives and independent evidence easier to combine.

s_\mu(y)=\frac{\partial}{\partial\mu}\ell(\mu;y)=\frac{y-\mu}{\sigma^2}

Differentiate with respect to the candidate mean. This log-likelihood slope is the score supplied by one observation.

One realized observation gives one score; Fisher information describes the expected strength of that clue across all observations the model could generate.

\mathbb E_\mu[s_\mu(Y)]=\frac{\mathbb E_\mu[Y]-\mu}{\sigma^2}=0

If μ is the generating mean, observations fall above and below it without a preferred direction, so their scores balance at zero.

I_1(\mu)=\mathbb E_\mu[s_\mu(Y)^2]=\frac{\mathbb E_\mu[(Y-\mu)^2]}{\sigma^4}=\frac{1}{\sigma^2}

Square the scores so opposite directions do not cancel, then average their strength. The Gaussian squared displacement averages to σ².

This is Fisher information from one observation. Small $\sigma$ produces a narrow likelihood: moving the candidate mean slightly changes the observation's likelihood sharply, so $I_1(\mu)=1/\sigma^2$ is large. Large $\sigma$ produces a broad likelihood and therefore less information about where the mean lies. Independent observations add information. With $n$ observations, $I_n(\mu)=nI_1(\mu)=n/\sigma^2$ , and their combined likelihood has standard deviation $\sigma/\sqrt n=1/\sqrt{I_n(\mu)}$ .

To turn those likelihoods into an estimate, multiply the likelihoods supplied by the observations—or, equivalently, add their log-likelihoods—and find the candidate mean at the resulting peak. At a smooth peak the total score, which is the slope of the combined log-likelihood, equals zero:

\ell_n(\mu)=\sum_{i=1}^{n}\ln p(y_i\mid\mu)

Independent observations multiply their likelihoods, so their log-likelihoods add.

\frac{\partial\ell_n}{\partial\mu}=\sum_{i=1}^{n}\frac{y_i-\mu}{\sigma^2}=\frac{n(\bar y-\mu)}{\sigma^2}

Add the scores from all observations and rewrite their sum using the sample mean ȳ.

\frac{n(\bar y-\hat\mu)}{\sigma^2}=0\qquad\Longrightarrow\qquad\hat\mu=\bar y

Set the combined slope to zero. This identifies the candidate mean at the stationary point.

\frac{\partial^2\ell_n}{\partial\mu^2}=-\frac{n}{\sigma^2}<0

The second derivative is negative everywhere, so the log-likelihood curves downward and that stationary point is its unique maximum.

The notation distinguishes the estimation rule from its realized answer. Before the data are observed, $\bar Y=n^{-1}\sum_iY_i$ is a random variable, and $\hat\mu=\bar Y$ names the estimator: the rule that will average whatever observations occur. After the particular values $y_1,\ldots,y_n$ are observed, the rule returns the number $\hat\mu=\bar y=n^{-1}\sum_i y_i$ . Averaged across repeated datasets, $\mathbb E_\mu[\bar Y]=\mu$ , so this estimator lands on the generating mean.

From observations to information about μ

One noisy observation gives a clue; independent observations sharpen it

Scroll horizontally to inspect the full calculation →

In the top row, the dashed line is the fixed value

\mu^\ast=0

before noise is added. Adding

\varepsilon\sim\mathcal N(0,1)

produces the teal Gaussian distribution of possible observations, and the dots are four values sampled from it. The bottom row asks the inverse question: which candidate mean could have produced those values? Each thin curve is the likelihood supplied by one observation. Their product is the bold combined likelihood, whose peak gives the estimate

\hat\mu=\bar y=0

for these data. One observation contributes

I_1(\mu)=1/\sigma^2

; four independent observations contribute

I_4(\mu)=4/\sigma^2

, reducing the likelihood's standard deviation from

\sigma

to

\sigma/2

. Fisher information is therefore already defined for one observation and accumulates as independent observations are added.

Why can we rely on the score being centered at zero?

The symmetry of the Gaussian model makes this centering plausible, but Fisher information is meant to work beyond Gaussians. The following argument establishes it for a broad class of models. The score is a random variable because it depends on the still-random observation $Y$ . Its average is zero whenever three ordinary assumptions hold: $p(y\mid\theta)$ changes smoothly with $\theta$ ; the set of possible observations does not change with $\theta$ ; and the parameter derivative may be moved through the integral. These assumptions are often bundled together under the name regularity conditions, and they hold for the Gaussian model above. They let us use the fact that the total probability $\int p(y\mid\theta)\,dy$ remains equal to one for every parameter value. Moving the derivative through the integral means exchanging the order of these two operations: $\int \partial_\theta p(y\mid\theta)\,dy=\partial_\theta\int p(y\mid\theta)\,dy$ . In other words, differentiating each infinitesimal contribution and then adding them must give the same result as adding the whole density first and then differentiating its total. Smoothness, fixed support, and a sufficiently well-behaved density make that exchange valid.

s_\theta(y)=\frac{\partial}{\partial\theta}\ln p(y\mid\theta)

Begin with the score: the parameter derivative of the log-density.

\mathbb E_\theta[s_\theta(Y)]=\int p(y\mid\theta)\,\frac{\partial}{\partial\theta}\ln p(y\mid\theta)\,dy

Write the expectation explicitly by weighting each possible score by its probability.

\frac{\partial}{\partial\theta}\ln p(y\mid\theta)=\frac{1}{p(y\mid\theta)}\frac{\partial p(y\mid\theta)}{\partial\theta}

Use the chain rule: the derivative of ln u is (1/u) times the derivative of u.

p(y\mid\theta)\,\frac{1}{p(y\mid\theta)}\frac{\partial p(y\mid\theta)}{\partial\theta}=\frac{\partial p(y\mid\theta)}{\partial\theta}

The probability multiplying the score cancels its reciprocal.

\mathbb E_\theta[s_\theta(Y)]=\int\frac{\partial p(y\mid\theta)}{\partial\theta}\,dy=\frac{\partial}{\partial\theta}\int p(y\mid\theta)\,dy

A regularity condition permits moving the parameter derivative outside the integral.

\mathbb E_\theta[s_\theta(Y)]=\frac{\partial}{\partial\theta}1=0

Every probability density integrates to one, and changing θ cannot change that total.

This result establishes the needed center point. Changing the parameter may increase the log-likelihood slope for some possible observations, but normalization forces compensating decreases elsewhere. When weighted by the model itself, those positive and negative scores balance to zero. The score spread is therefore centered on “no preferred direction,” allowing its spread to measure how strongly observations separate nearby parameter values rather than how strongly they collectively recommend moving somewhere else.

From balanced scores to Fisher information

The result above explains why the score can indicate a direction for one particular observation without favoring any direction across all possible observations generated at $\theta$ . To measure how strongly a typical observation distinguishes nearby parameter values, the opposing scores must be combined without allowing their signs to cancel. Squaring accomplishes this. Because $\operatorname{Var}(S)=\mathbb E[S^2]-\mathbb E[S]^2$ and the expected score is zero, the score variance becomes $\operatorname{Var}_\theta[s_\theta(Y)]=\mathbb E_\theta[s_\theta(Y)^2]$ . This is the scalar Fisher information $I_Y(\theta)$ already calculated for the Gaussian above. A large value means that observations typically give strong local evidence for distinguishing nearby parameter values; a small value means that nearby values remain difficult to tell apart.

Score view

I(\theta)=\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]

The covariance of the log-likelihood slopes.

Curvature view

I(\theta)=-\mathbb E_\theta[\nabla_\theta^2\ell_\theta]

The average sharpness of the log-likelihood peak.

Estimation view

\operatorname{Cov}(\hat\theta)\succeq I_n(\theta)^{-1}

A lower bound on unbiased-estimator uncertainty.

The Fisher matrix is not itself the gradient of a distribution. The score is one observation's log-likelihood gradient; the first formula above is its covariance across possible observations. In full, that covariance is $\mathbb E_\theta[(s_\theta-\mathbb E_\theta[s_\theta])(s_\theta-\mathbb E_\theta[s_\theta])^{\mathsf T}]$ . Because the expected score is zero, it reduces to $\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]$ . The second formula is not another covariance formula: it is the equivalent curvature view. It uses the Hessian, the matrix of second derivatives of the log-likelihood. In one dimension, the Hessian is a single ordinary second derivative $\partial_\theta^2\ell$ , or equivalently a $1\times1$ matrix. Fisher information is the negative expected value of that Hessian, not generally the Hessian produced by one particular observation. The equality between the score and curvature views follows from probability normalization under these same assumptions; Step 4 of the later KL Taylor derivation proves the identity explicitly.

The estimation view in the panel above now follows from the same zero-score identity. Suppose $\hat\theta$ is unbiased, meaning $\mathbb E_\theta[\hat\theta]=\theta$ across repeated datasets. Differentiating this equality, then applying the Cauchy–Schwarz inequality, connects estimator uncertainty to Fisher information:

\mathbb E_\theta[\hat\theta]=\theta\qquad\Longrightarrow\qquad\frac{\partial}{\partial\theta}\mathbb E_\theta[\hat\theta]=1

Unbiasedness says that the estimator averages to the true parameter; differentiate both sides with respect to that parameter.

1=\mathbb E_\theta[\hat\theta\,U_n(\theta)]=\mathbb E_\theta[(\hat\theta-\theta)U_n(\theta)]

Move the derivative through the expectation and use the zero expected score. Here Uₙ is the score of the complete dataset.

1\leq\mathbb E_\theta[(\hat\theta-\theta)^2]\,\mathbb E_\theta[U_n(\theta)^2]=\operatorname{Var}_\theta(\hat\theta)I_n(\theta)

Cauchy–Schwarz bounds the squared relationship between estimator error and score by the product of their squared spreads.

\operatorname{Var}_\theta(\hat\theta)\geq\frac{1}{I_n(\theta)}

Rearrange to obtain the general scalar Cramér–Rao bound for an unbiased estimator.

I_n(\mu)=\frac{n}{\sigma^2}\qquad\Longrightarrow\qquad\operatorname{Var}_\mu(\hat\mu)\geq\frac{\sigma^2}{n}

Only this final substitution is specific to the Gaussian mean model. Its sample mean attains the bound.

More observations accumulate more Fisher information, which reduces uncertainty about the value of the generating mean.

From one parameter to a parameter space

With several parameters, the score becomes a gradient vector and Fisher information becomes its covariance matrix. Large curvature in a direction means observations distinguish movements in that parameter combination well; small curvature means a poorly identified direction. A smoothly parameterized, identifiable family of distributions is a statistical manifold: each parameter choice names a point, and nearby parameters name nearby models. A categorical family is the special case that can be drawn as a probability simplex; a Gaussian family has a different parameter space. Attaching the Fisher matrix to every point equips the statistical manifold with a metric—a local rule for assigning statistical size to movements in each parameter direction.

Why Fisher information appears in the KL bowl

Both constructions ask how sensitively a probability model responds to a small parameter change. Estimation theory measures that sensitivity through score variance or expected log-likelihood curvature. Local KL measures how quickly a nearby model becomes distinguishable from the current one. When KL is Taylor-expanded, normalization removes the linear score term and the first surviving coefficient is the negative expected log-density Hessian—the same $I(\theta)$ defined above. The derivation that follows proves this identification explicitly.

Estimation perspective and equivalent formulas: measurement intuition, score and curvature visualization, and normalization-based derivation.

4

Follow the algebra

The full derivation of the local Fisher approximation to KL

We'll now show the full derivation of the local Fisher approximation to KL.

One convention before the algebra

Let $\ell_\theta(x)=\ln p_\theta(x)$ . Natural logarithms measure information in nats and give the standard Fisher matrix without an extra constant. KL in bits is the same quantity divided by $\ln 2$ .

1

Insert the nearby model into the cross-entropy definition

Set $P=p_\theta$ and $Q=p_{\theta+\Delta\theta}$ . Both expectations are over data drawn from the unchanged distribution $p_\theta$ .

D_{\mathrm{KL}}(P\|Q)=H(P;Q)-H(P)

Start with KL as the extra coding cost: cross-entropy minus entropy.

H(P;Q)=\mathbb E_{X\sim P}[-\ln Q(X)]

Write out cross-entropy: average Q’s negative log-probability over outcomes generated by P.

H(P)=\mathbb E_{X\sim P}[-\ln P(X)]

Write out entropy: average P’s own negative log-probability over the same outcomes.

H(p_\theta;p_{\theta+\Delta\theta})=\mathbb E_\theta[-\ln p_{\theta+\Delta\theta}(X)]

Replace P by pθ and Q by the nearby model. The averaging distribution remains pθ.

H(p_\theta)=\mathbb E_\theta[-\ln p_\theta(X)]

Make the same replacement in the entropy term.

D_{\mathrm{KL}}=\mathbb E_\theta[-\ln p_{\theta+\Delta\theta}(X)]-\mathbb E_\theta[-\ln p_\theta(X)]

Substitute those two expanded definitions back into cross-entropy minus entropy.

D_{\mathrm{KL}}=\mathbb E_\theta[-\ell_{\theta+\Delta\theta}(X)]-\mathbb E_\theta[-\ell_\theta(X)]

Use ℓϑ(x) = ln pϑ(x) as shorthand for each pointwise log-density.

D_{\mathrm{KL}}=\mathbb E_\theta[\ell_\theta(X)-\ell_{\theta+\Delta\theta}(X)]

Use linearity to combine the averages and simplify the minus signs.

The nearby log-density $\ell_{\theta+\Delta\theta}(X)$ is a value assigned to the random outcome $X$ ; it is not itself an entropy or cross-entropy. Only $\mathbb E_\theta[-\ell_{\theta+\Delta\theta}(X)]$ —the negative log-density averaged over data from $p_\theta$ —is the cross-entropy. Likewise, $\mathbb E_\theta[-\ell_\theta(X)]=H(p_\theta)$ is the entropy of the current model.

2

Taylor-expand the nearby model's log-density

The general problem a Taylor expansion addresses is this: when a smooth function is difficult to evaluate directly, but its value and derivatives are known at one input, how can we estimate the value of different inputs that are nearby? It replaces the difficult function locally with a polynomial built from its value, slope, curvature, and higher derivatives at that known input. Here the known input is $\theta$ and the nearby input is $\theta+\Delta\theta$ . In the exact KL expression above, all dependence on the nearby model sits inside $\ell_{\theta+\Delta\theta}(X)$ , so this is the term we expand.

Define the score $s_\theta=\nabla_\theta\ell_\theta$ , the vector of first derivatives, and the log-density Hessian $J_\theta=\nabla_\theta^2\ell_\theta$ . With one parameter, the Hessian is simply the second derivative. With several parameters, it is the matrix of every second partial derivative, recording curvature along—and interactions between—the parameter directions. Then

\ell_{\theta+\Delta\theta}(x)=\ell_\theta(x)+s_\theta(x)^{\mathsf T}\Delta\theta+\frac12\Delta\theta^{\mathsf T}J_\theta(x)\Delta\theta+O(\|\Delta\theta\|^3)

Apply the multivariable Taylor expansion to the log-density around θ.

D_{\mathrm{KL}}=\mathbb E_\theta[\ell_\theta(X)-\ell_{\theta+\Delta\theta}(X)]

Return to the expression for KL obtained in Step 1.

D_{\mathrm{KL}}=\mathbb E_\theta\!\left[\ell_\theta(X)-\left(\ell_\theta(X)+s_\theta(X)^{\mathsf T}\Delta\theta+\frac12\Delta\theta^{\mathsf T}J_\theta(X)\Delta\theta+O(\|\Delta\theta\|^3)\right)\right]

Replace the nearby log-density inside the exact KL formula by its Taylor expansion.

D_{\mathrm{KL}}=\mathbb E_\theta\!\left[-s_\theta(X)^{\mathsf T}\Delta\theta-\frac12\Delta\theta^{\mathsf T}J_\theta(X)\Delta\theta+O(\|\Delta\theta\|^3)\right]

The two copies of ℓθ(X) cancel inside the expectation.

D_{\mathrm{KL}}=-\mathbb E_\theta[s_\theta]^{\mathsf T}\Delta\theta-\frac12\Delta\theta^{\mathsf T}\mathbb E_\theta[J_\theta]\Delta\theta+O(\|\Delta\theta\|^3)

Take the expectation term by term; Δθ is fixed, so it passes through the expectation.

The complicated change in log-density has now become a local slope term and a local curvature term, plus an error that begins at third order. To learn the actual local shape of KL, we next determine what happens to each coefficient after averaging over $p_\theta$ .

3

The linear term vanishes because every distribution normalizes to one

Carry the result from Step 2 forward. Its first coefficient is the expected score $\mathbb E_\theta[s_\theta]$ :

D_{\mathrm{KL}}=-\mathbb E_\theta[s_\theta]^{\mathsf T}\Delta\theta-\frac12\Delta\theta^{\mathsf T}\mathbb E_\theta[J_\theta]\Delta\theta+O(\|\Delta\theta\|^3)

This is the KL approximation obtained at the end of Step 2; we now evaluate its linear coefficient.

\mathbb E_\theta[s_\theta]=\int p_\theta(x)\,\nabla_\theta\ln p_\theta(x)\,dx

Expand the expectation of the score as an integral over pθ.

\mathbb E_\theta[s_\theta]=\int \nabla_\theta p_\theta(x)\,dx

Use p∇ ln p = ∇p, which follows from the chain rule.

\mathbb E_\theta[s_\theta]=\nabla_\theta\int p_\theta(x)\,dx

Move the parameter derivative outside the integral under the regularity assumptions below.

\mathbb E_\theta[s_\theta]=\nabla_\theta 1

Every probability distribution is normalized: its total probability is one.

\mathbb E_\theta[s_\theta]=0

Differentiate the constant one. The score therefore has zero mean.

D_{\mathrm{KL}}=-\frac12\Delta\theta^{\mathsf T}\mathbb E_\theta[J_\theta]\Delta\theta+O(\|\Delta\theta\|^3)

Substitute Eθ[sθ] = 0 back into KL. The entire linear term disappears; the Hessian term remains.

This normalization identity is the algebraic reason “KL is minimal where the models match” removes the first-order term. We have not left the original KL formula; we have simplified its Taylor approximation one coefficient at a time.

4

The remaining Hessian is the Fisher information matrix

Step 3 has already completed the local KL approximation in Hessian form. If all we wanted were a quadratic coefficient, we could stop at $-\mathbb E_\theta[J_\theta]$ . The remaining task is to show that this coefficient is the same Fisher information introduced earlier as the expected score outer product $\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]$ . That identification connects the local curvature of KL to the variability of the score.

Step 3 differentiated normalization once because its unwanted coefficient contained a first derivative, the score. The coefficient that remains now contains a second derivative, the Hessian. Differentiating normalization twice is therefore how we produce an identity involving that second derivative:

D_{\mathrm{KL}}=-\frac12\Delta\theta^{\mathsf T}\mathbb E_\theta[J_\theta]\Delta\theta+O(\|\Delta\theta\|^3)

Carry forward the quadratic KL expression obtained at the end of Step 3.

-\mathbb E_\theta[J_\theta]\stackrel{?}{=}\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]

State what remains to be proved: the Hessian coefficient in KL equals the score-based definition of Fisher information.

\int p_\theta(x)\,dx=1\qquad\Longrightarrow\qquad\nabla_\theta^2\!\int p_\theta(x)\,dx=\nabla_\theta^2 1=0

Normalization holds for every θ. Differentiate it twice because the coefficient we need to identify contains second derivatives.

0=\int\nabla_\theta^2p_\theta(x)\,dx

Move the second parameter derivative through the integral under the regularity assumptions.

\nabla_\theta p_\theta(x)=p_\theta(x)\,\nabla_\theta\ln p_\theta(x)=p_\theta(x)s_\theta(x)

Rewrite the first derivative of the density using the score.

\partial_i\partial_jp_\theta(x)=\partial_i\!\bigl(p_\theta(x)s_{\theta,j}(x)\bigr)

To obtain one entry of the density Hessian, differentiate the first-derivative identity again.

\partial_i\partial_jp_\theta=(\partial_i p_\theta)s_{\theta,j}+p_\theta\,\partial_i s_{\theta,j}

Apply the product rule: one term differentiates pθ, while the other differentiates the score component.

\partial_i\partial_jp_\theta=p_\theta s_{\theta,i}s_{\theta,j}+p_\theta J_{\theta,ij}

Use ∂ᵢpθ = pθsθ,ᵢ and ∂ᵢsθ,ⱼ = Jθ,ᵢⱼ in the two product-rule terms.

\nabla_\theta^2p_\theta(x)=p_\theta(x)\bigl(s_\theta(x)s_\theta(x)^{\mathsf T}+J_\theta(x)\bigr)

Collect all component identities into matrix form. The density Hessian contains both the score outer product and the log-density Hessian.

0=\int p_\theta(x)\bigl(s_\theta s_\theta^{\mathsf T}+J_\theta\bigr)\,dx

Substitute the expanded density Hessian into the twice-differentiated normalization identity.

0=\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]+\mathbb E_\theta[J_\theta]

Recognize each probability-weighted integral as an expectation under pθ.

-\mathbb E_\theta[J_\theta]=\mathbb E_\theta[s_\theta s_\theta^{\mathsf T}]\equiv I(\theta)

Move the expected Hessian to the other side. The remaining score outer product is the Fisher information matrix.

In one dimension, this identity says $-\mathbb E_\theta[\ell_\theta''(X)]=\mathbb E_\theta[(\ell_\theta'(X))^2]$ . The average downward curvature of the log-likelihood and the average squared strength of its slope are two measurements of the same local sensitivity. We can now replace the Hessian coefficient in KL by the Fisher information matrix:

D_{\mathrm{KL}}=-\frac12\Delta\theta^{\mathsf T}\mathbb E_\theta[J_\theta]\Delta\theta+O(\|\Delta\theta\|^3)

Return to the Hessian form of the local KL approximation.

-\mathbb E_\theta[J_\theta]=I(\theta)

Use the identity just derived to rename the quadratic coefficient.

\boxed{D_{\mathrm{KL}}^{(\mathrm{nats})}(p_\theta\|p_{\theta+\Delta\theta})=\frac12\Delta\theta^{\mathsf T}I(\theta)\Delta\theta+O(\|\Delta\theta\|^3)}

Substitute I(θ) into KL. This is the local Fisher approximation.

D_{\mathrm{KL}}^{(\mathrm{bits})}=\frac{1}{\ln 2}\,D_{\mathrm{KL}}^{(\mathrm{nats})}

Convert nats to bits by dividing by ln 2.

This calculation assumes the model is smooth, its support does not change with $\theta$ , and differentiation may pass through the integral. For discrete outcomes, replace each integral by a sum.

From local KL to a geometric length

The derivation has produced a quadratic form. For any chosen parameter movement, multiplying the movement by the Fisher matrix in this way produces one scalar:

\Delta\theta^{\mathsf T}I(\theta)\Delta\theta

The three objects in this expression have separate roles:

$I(\theta)$ The Fisher information matrix at the current distribution; it supplies the local ruler.

$\Delta\theta$ The parameter movement being measured; it supplies a direction and coordinate magnitude.

$\Delta\theta^{\mathsf T}I(\theta)\Delta\theta$ The scalar produced by applying the ruler to that movement; it supplies a squared statistical length.

The Taylor derivation uses $\Delta\theta$ for a small, finite parameter change. Geometry takes the infinitesimal limit and writes that movement as $d\theta$ . The symbol $s$ denotes statistical path length, so information geometry defines its infinitesimal squared length by

THE FISHER LINE ELEMENT

ds^2\equiv d\theta^{\mathsf T}I(\theta)d\theta,\qquad D_{\mathrm{KL}}(p_\theta\|p_{\theta+d\theta})\approx\frac12ds^2

With one parameter, the matrix has one entry and the construction reduces to $ds^2=I(\theta)(d\theta)^2$ . With two parameters, the multiplication expands explicitly:

TWO-PARAMETER FORM

I(\theta)=\begin{bmatrix}a&b\\b&c\end{bmatrix},\quad d\theta=\begin{bmatrix}d\theta_1\\d\theta_2\end{bmatrix}\quad\Longrightarrow\quad ds^2=a(d\theta_1)^2+2b\,d\theta_1d\theta_2+c(d\theta_2)^2

The diagonal entries measure sensitivity along each coordinate direction, while the off-diagonal entry records how the two parameter directions interact. In every dimension, the complete multiplication returns one squared length.

KL has not become a global distance: it is still asymmetric for separated distributions. What loses orientation is only the second-order term, because the quadratic form $\Delta\theta^{\mathsf T}I(\theta)\Delta\theta$ is symmetric under $\Delta\theta\mapsto-\Delta\theta$ . The higher-order terms retain KL's directional information. This symmetric microscopic part is enough to supply a local metric; the directional higher-order part returns later when we reach dual connections.

What “zooming in” and “ruler” mean

Zooming in means holding one distribution $p_\theta$ fixed, restricting attention to nearby members of its parameterized family $p_{\theta+\Delta\theta}$ , and taking $\Delta\theta$ small enough that cubic and higher terms are negligible. We are studying how the entire probability output $p_\theta(x)$ changes with its parameters, not merely how the parameter numbers change on paper.

The Fisher line element $ds^2=d\theta^{\mathsf T}I(\theta)d\theta$ converts a small coordinate displacement into the leading-order KL distinguishability of the resulting distributions. Equal-sized parameter changes can therefore have different statistical sizes at different models or in different directions. This ruler is local: measuring a long route requires attaching it at every distribution along the path and adding the resulting infinitesimal lengths.

What the zoom accomplished: KL supplied a local ruler without becoming a distance itself. Repeating this construction at every

p_\theta

gives every distribution its own ruler. How that ruler changes from point to point is the geometry of the statistical manifold in Step 5.

5

Repeat at every possible model

The family becomes a statistical manifold

Collect every distribution in the model family into one space $\mathcal M$ . Each distribution is a point, its parameters $\theta$ are coordinates, and $d\theta$ is a tangent direction—an infinitesimal way to change the distribution. Attaching the Fisher matrix to every point gives the Fisher metric.

\mathcal M=\{p_\theta(x):\theta\in\Theta\},\qquad ds^2=d\theta^{\mathsf T}I(\theta)d\theta

A manifold is a space that looks like ordinary Euclidean space when viewed locally. The surface of Earth is curved globally, but a neighborhood admits an ordinary flat map. Likewise, a family $p(x\mid\theta)$ can be described locally by its parameter coordinates $\theta$ .

A family is a smoothly parameterized collection with a common mathematical form: all Bernoulli coins as $p$ varies, all Gaussians as $(\mu,\sigma)$ vary, or all Beta distributions as $(\alpha,\beta)$ vary. Paths and distances remain inside one chosen family, where the required derivatives are defined. Differently named model classes can be compared geometrically when they are embedded in a larger common statistical manifold, but unrelated parameterizations do not automatically provide a path between them.

PointOne complete distribution

p(x\mid\theta)

CoordinatesThe parameters

\theta

used to name it

Tangent vectorAn infinitesimal way to change the distribution

MetricA rule for measuring tiny displacements

ConnectionA rule for saying whether a direction changed while moving

GeodesicA path with zero turning according to a chosen connection

From points to paths

A continuously changing distribution traces a path through the manifold. If its parameters follow $\theta(t)$ , their derivative is a tangent direction: it records how the entire distribution is changing at that instant, just as the velocity of an ordinary curve records how a point is moving.

A DISTRIBUTION PATH AND ITS TANGENT

\dot p(x)=\frac{d}{dt}p_{\theta(t)}(x)=\nabla_\theta p_\theta(x)^{\mathsf T}\dot\theta

The Fisher metric assigns that tangent a statistical speed. Adding those infinitesimal lengths along the path gives its total Fisher–Rao length:

FISHER SPEED AND PATH LENGTH

\begin{aligned} \|\dot p\|_F^2&=\dot\theta^{\mathsf T}I(\theta)\dot\theta\\[4pt] L[\theta]&=\int\sqrt{\dot\theta^{\mathsf T}I(\theta)\dot\theta}\,dt \end{aligned}

Metric distance is obtained by comparing the lengths of paths joining two distributions within the chosen manifold. This is where the familiar shortest-path intuition enters—but it must be separated from the more general meaning of a geodesic.

Important correction about geodesics

“The shortest path” is a good first intuition for a Riemannian geodesic using the metric's Levi–Civita connection. More generally, a geodesic means “move without turning according to this connection.” An affine geodesic need not minimize metric distance. This distinction is exactly what we need for dual geometry.

Probability coordinate $p\in(0,1)$ The fair coin is

p=0.5

.

same coins, renamed

Log-odds coordinate $\theta=\log\frac{p}{1-p}$ The fair coin is

\theta=0

.

Changing coordinates changes the numbers in the Fisher matrix, but not the intrinsic length $ds$ . That coordinate independence is what earns Fisher information the word metricrather than merely “a useful Hessian.” Rao made this geometric viewpoint explicit.

Why this still is not the end: A metric tells us lengths and angles. It does not, by itself, tell us how to carry a direction from one distribution to another and decide whether a path kept going “straight.” That requires a connection.

6

Add a rule for preserving direction

Two statistical notions of straightness appear

Information geometry supplies two natural connections. The mixture connection calls a path straight when probabilities change at a constant rate. The exponential connection calls it straight when log-probability ratios change at a constant rate. They disagree about the route, but the Fisher metric pairs them as dual descriptions.

mixture-straight

p_t=(1-t)p_0+tp_1

exponential-straight

p_t(x)\propto p_0(x)^{1-t}p_1(x)^t

Interactive categorical manifold

Two routes between the same distributions

Travel from P₀ to P₁: t = 0.50

Mixture path: average probabilities Exponential path: average log-probabilities

m(t)A: 0.400B: 0.225C: 0.375

e(t)A: 0.356B: 0.301C: 0.343

Mixture route

m_i(t)=(1-t)P_{0i}+tP_{1i}

Constant velocity in probability coordinates.

Exponential route

e_i(t)\propto P_{0i}^{1-t}P_{1i}^{t}

Constant velocity in log-ratio coordinates.

m-connectionMixture-straight

Probabilities change linearly. In probability coordinates, this route looks visually straight.

e-connectionExponential-straight

Log ratios change linearly. In probability coordinates, this route usually looks curved.

Levi–Civita connectionMetric-straight

The Fisher–Rao geodesic is locally distance-minimizing. It is generally a third route.

So why call the first two dual?

The Fisher metric pairs the mixture and exponential connections so that what one connection does to a direction is balanced by what the other does to the direction it is paired with. In exponential families, natural parameters

\theta

are straight for the e-connection, expectation parameters

\eta

are straight for the m-connection, and

\theta

and

\eta

are related by Legendre duality.

Where KL asymmetry enters

The asymmetry is not the entire definition of dual geometry; it is the visible clue. Expanding a divergence around $P=Q$ produces one symmetric metric from its second-order term. Its third-order terms retain orientation and generate a pair of dual connections. Reversing KL swaps their roles.

THE LOCAL PICTURE

D\xrightarrow{\text{second derivatives}}g\qquad D\xrightarrow{\text{third derivatives}}(\nabla,\nabla^{*})

The simplest honest summary

A metric tells you how much an infinitesimal move costs. A connection tells you whether you turned during the move. Information geometry has one canonical metric and two statistically meaningful, mutually paired notions of not turning.

Dual connections and divergence geometry: Amari & Nagaoka (2000); Nielsen (2020).

This is the destination: “Dual geometry” means one Fisher metric together with two mutually paired rules for not turning. KL asymmetry helps generate this pair, and the probability triangle above makes their different routes visible.

07

Fisher information: zoom into KL

KL divergence compares separated distributions. Now place two models extremely close together: $p_\theta$ and $p_{\theta+d\theta}$ . As derived in the bridge above, normalization makes the first-order score average vanish; the leading term is quadratic. Its coefficient is Fisher information.

LOCAL KL EXPANSION (NATS)

D_{\mathrm{KL}}(p_\theta\|p_{\theta+d\theta})\approx\frac12 d\theta^{\mathsf T}I(\theta)d\theta

FISHER INFORMATION FOR A BERNOULLI COIN

I(p)=\mathbb{E}\!\left[\left(\frac{\partial}{\partial p}\log P(X\mid p)\right)^2\right]=\frac{1}{p(1-p)}

The formula for $I$ depends on the family being studied. A Bernoulli family has one parameter $p$ , so its Fisher information is one number. A Gaussian family parameterized by $(\mu,\sigma)$ has a $2\times2$ Fisher matrix; a Beta family parameterized by $(\alpha,\beta)$ has another. In every case, the entries describe how sensitive that particular family is to small movements in its own parameter directions.

Interactive local geometry

How distinguishable are two nearby coins?

Starting coin p = 0.50Small change dp = 0.020Independent observations N = 100

P

p = 0.500

local separation⟷

Nearby Q

p = 0.520

Fisher informationI(p) = 4.000sensitivity at this coin

Exact KL0.000801 natsDₖₗ(P‖Q)

Local approximation0.000800 nats½ I(p) dp²

Best-case estimator SD≥ 0.0500Cramér–Rao, N observations

Exact KL0.0008006

Fisher approximation0.0008000

Fisher–Rao distance0.04001

Try moving $p$ toward 0 or 1. The same numerical change $dp$ becomes more informative because a rare outcome changes by a large proportion. That is why $I(p)=1/[p(1-p)]$ grows near the edges.

The same numerical step need not be the same statistical step

For a Bernoulli coin,

I(0.5)=4

, while

I(0.01)\approx101

. Therefore, for the same sufficiently small change

dp

, the quadratic KL cost near

p=0.01

is about twenty-five times the cost near

p=0.5

. Near the center, a tiny probability change is difficult to detect from data. Near an edge, that same absolute change is a large relative change in a rare event and is much easier to distinguish. The Fisher ruler records this difference even though an ordinary number-line ruler sees the same

dp

.

Fisher information versus the Fisher information metric

At one modelFisher information

$I(p)$ is one local sensitivity value—or a matrix when there are multiple parameters.

Across all modelsFisher metric

Attach that matrix to every point. The smoothly varying field defines $ds^2=I_{ij}\,d\theta^i d\theta^j$ .

Same object, different viewpoint

“Fisher information matrix” emphasizes estimation at a particular distribution. “Fisher information metric” emphasizes how those matrices fit together over the whole family. The metric is what turns local distinguishability into lengths and angles.

Derive I(p) for the Bernoulli distribution

For $X\in\{0,1\}$ , the Bernoulli log-likelihood is:

\log P(X\mid p)=X\log p+(1-X)\log(1-p)

\frac{\partial}{\partial p}\log P(X\mid p)=\frac{X}{p}-\frac{1-X}{1-p}

The score is $1/p$ when $X=1$ and $-1/(1-p)$ when $X=0$ . Square it and average using probabilities $p$ and $1-p$ :

I(p)=p\left(\frac1{p^2}\right)+(1-p)\left(\frac1{(1-p)^2}\right)=\frac1{p(1-p)}

Statistical information and its geometry: Fisher (1922); Rao (1945).

08

Bonus: natural-gradient variational inference

Fisher geometry becomes operational when a probability distribution is itself being optimized. Variational inference does exactly that: it adjusts an approximate posterior $q_\phi(s)$ until it resembles the Bayesian posterior over the hidden state $s$ . The subscript $\phi$ names the adjustable coordinates used to represent the beliefs. If $\phi$ is the vector of probabilities, changing $\phi$ changes those probabilities directly. If $\phi$ is a vector of logits, the logits change first and $q_\phi$ changes through a softmax.

$s$ The hidden state whose value is uncertain

$\theta$ Parameters of the generative model

p_\theta(o,s)

; they may be fixed or learned

$\phi$ Coordinates of the approximate posterior

q_\phi(s)

that variational inference adjusts

These roles should not be collapsed. Inferring $s$ with $q_\phi(s)$ does not require learning the generative parameters $\theta$ . Those parameters may be supplied mechanically, estimated from previous data, or placed under their own prior and inferred in a larger model. Merely writing $p(s\mid\theta)$ says that the state distribution depends on $\theta$ ; it does not make $\theta$ an inferred quantity. Parameter inference additionally requires a prior $p(\theta)$ and data, producing a posterior such as $p(\theta\mid o)$ or a variational approximation $q(\theta)$ .

From a free-energy gradient to a belief update

A gradient is a vector: an arrow through the coordinate space. At the current coordinates $\phi$ , the free-energy gradient $\nabla_\phi F$ points in the direction where a small movement increases free energy fastest. Multiplying a vector by $-1$ reverses its arrow, so $-\nabla_\phi F$ points in the direction where free energy decreases fastest locally. The components of this reversed vector may be positive, negative, or zero, depending on the current location. Ordinary gradient descent scales the reversed vector by a positive learning rate $\rho$ :

ORDINARY COORDINATE-GRADIENT STEP

\Delta\phi=-\rho\nabla_\phi F

The gradient magnitude and $\rho$ together determine how far the coordinates move. The resulting belief update is $q_\phi(s)\rightarrow q_{\phi+\Delta\phi}(s)$ . But coordinate distance is not yet belief distance. A useful belief distance should describe how differently two models predict observations—how readily observations could distinguish them—not merely how far their parameter labels moved. Consider what happens if we try to define the distance between distributions using only the distance between their coordinates.

Write a binary distribution as $q_p=(p,1-p)$ , where the single number $p$ is its coordinate within the binary family. Consider $q_A=(0.50,0.50)$ and $q_B=(0.25,0.75)$ . Their coordinates are $p_A=0.50$ and $p_B=0.25$ , so the signed coordinate displacement from $q_A$ to $q_B$ is $\Delta p=p_B-p_A=-0.25$ and its coordinate distance is $|\Delta p|=0.25$ . Component by component, the probability vector changes by $q_B-q_A=(-0.25,+0.25)$ .

This coordinate distance does not yet say how distinguishable the two distributions are. To see the problem, compare two movements in the same probability coordinate with exactly the same numerical displacement, $\Delta p=-0.001$ :

q_{0.500}=(0.500,0.500)\longrightarrow q_{0.499}=(0.499,0.501)

Near the center of the family, the first probability decreases by 0.001.

D_{\mathrm{KL}}(q_{0.500}\|q_{0.499})\approx2.00\times10^{-6}\ \text{nats}

These two nearly balanced distributions are extremely difficult to distinguish.

q_{0.020}=(0.020,0.980)\longrightarrow q_{0.019}=(0.019,0.981)

Near the edge of the family, the first probability again decreases by exactly 0.001.

D_{\mathrm{KL}}(q_{0.020}\|q_{0.019})\approx2.64\times10^{-5}\ \text{nats}

The same parameter step causes more than thirteen times as much distributional change near the rare-event edge.

KL is useful here because it compares the probabilities that the two models assign to the same possible observations: $D_{\mathrm{KL}}(P\|Q)=\mathbb E_{X\sim P}[\ln(P(X)/Q(X))]$ . The expected log-ratio measures how much evidence an observation generated by $P$ provides for $P$ over $Q$ ; equivalently, it is the expected excess coding cost of using $Q$ for data generated by $P$ . Models that make nearly identical predictions have a small KL divergence, while models whose predictions are easier to distinguish have a larger one.

If we naively use coordinate distance as the distance between distributions, both pairs are exactly $|\Delta p|=0.001$ apart. Their KL divergences, however, are $2.00\times10^{-6}$ and $2.64\times10^{-5}$ nats. KL therefore distinguishes distributional changes that coordinate distance treats as identical.

Raw coordinate displacement has another limitation: its numerical size depends on how the family is parameterized. The same binary family can instead be parameterized by the log-odds $\phi=\ln\!\bigl(p/(1-p)\bigr)$ . For example, here is how we calculate the log-odds coordinates for distributions $q_A$ and $q_B$ :

\phi_A=\ln\frac{p_A}{1-p_A}=\ln\frac{0.50}{0.50}=\ln1=0

The 50/50 distribution has odds of one to one, so its log-odds coordinate is zero.

\phi_B=\ln\frac{p_B}{1-p_B}=\ln\frac{0.25}{0.75}=\ln\frac13\approx-1.10

The 25/75 distribution has odds of one to three, which gives the negative log-odds coordinate.

\Delta\phi=\phi_B-\phi_A=\ln\frac13-0\approx-1.10

The same two distributions are now separated by about 1.10 units of log-odds rather than 0.25 units of probability.

The identical change from $q_A$ to $q_B$ therefore has displacement $\Delta p=-0.25$ in probability coordinates but $\Delta\phi\approx-1.10$ in log-odds coordinates. The distributions did not change; only the coordinates used to describe them changed. Together, the two examples expose both problems with raw coordinate distance: different distributional changes can produce equivalent coordinate distances, and these coordinate distances change depending on the particular parameterization.

KL divergence avoids both problems by comparing the probability distributions themselves rather than their coordinate labels. For a sufficiently small movement from $\phi$ to $\phi+\Delta\phi$ , the local KL expansion approximates that distributional change. The Fisher matrix supplies the location-dependent conversion from coordinate displacement to local statistical change: it records how sensitively the current distribution $q_\phi$ responds to movement in each parameter direction.

KL is not the only useful global comparison between distributions; Jensen–Shannon and Hellinger divergences, for example, answer related questions. Information geometry supplies the deeper local result: KL's quadratic term is the Fisher metric, and a broad class of other smooth statistical divergences produce that same metric locally, up to an overall constant.

Canonical Fisher geometry: Amari & Nagaoka (2000); Nielsen (2020).

LOCAL BELIEF CHANGE

D_{\mathrm{KL}}\!\left(q_\phi\middle\|q_{\phi+\Delta\phi}\right)\approx\frac12\Delta\phi^{\mathsf T}I(\phi)\Delta\phi

For the binary family in probability coordinates, $I(p)=1/[p(1-p)]$ . Substituting the two starting locations and the shared step $\Delta p=-0.001$ shows how Fisher information produces the different local KL changes:

I(0.500)=\frac1{0.500(0.500)}=4

At the balanced distribution, Fisher information assigns relatively low sensitivity to movement in p.

\frac12 I(0.500)(\Delta p)^2=\frac12(4)(0.001)^2=2.00\times10^{-6}

The local approximation reproduces the tiny KL change near the center.

I(0.020)=\frac1{0.020(0.980)}\approx51.02

At the rare-event distribution, the same probability coordinate is much more sensitive.

\frac12 I(0.020)(\Delta p)^2\approx\frac12(51.02)(0.001)^2\approx2.55\times10^{-5}

The parameter step is unchanged; the larger Fisher information makes its statistical effect much larger.

Large Fisher information in a directionThe beliefs are sensitive

A small coordinate movement already causes a large statistical change, so the inverse Fisher suppresses movement in that direction.

Small Fisher information in a directionThe beliefs are insensitive

A larger coordinate movement is needed to cause the same statistical change, so the inverse Fisher permits more movement.

Imagine every possible small parameter movement of the same length and ask which direction decreases $F$ most. That is precisely the question answered by the ordinary gradient, and it can be visualized as equally long arrows extending from the current point on the statistical manifold. The negative gradient selects the arrow that produces the greatest predicted decrease in free energy.

Equal coordinate lengths, however, need not produce equal changes in the represented distribution. A short movement in a Fisher-sensitive direction may change $q_\phi$ substantially, while a longer movement in an insensitive direction may change it very little. To make “equal coordinate length” precise, suppose $d=(d_1,\ldots,d_k)$ describes a candidate arrow through the parameter coordinates. Its ordinary Euclidean length is $\sqrt{d_1^2+\cdots+d_k^2}$ . Giving every candidate arrow the same length $c$ therefore means $d_1^2+\cdots+d_k^2=c^2$ , written more compactly as $\lVert d\rVert_2^2=c^2$ . The subscript $2$ names the Euclidean, or $L^2$ , norm. These are unit arrows only when $c=1$ .

Fisher information lets us replace “equal coordinate length” with “equal local belief change.” Instead of placing every candidate on that Euclidean circle, we give every candidate the same locally predicted KL budget:

\frac12d^{\mathsf T}I(\phi)d=\varepsilon

The earlier local-KL derivation explains where this new constraint comes from. Taylor-expanding the nearby model's log-density $\ln q_{\phi+d}(s)$ introduces the score as its linear term and the log-density Hessian as its quadratic term. After averaging inside KL, the expected score vanishes and the quadratic coefficient becomes the Fisher information matrix:

D_{\mathrm{KL}}\!\left(q_\phi\middle\|q_{\phi+d}\right)=\frac12d^{\mathsf T}I(\phi)d+O(\lVert d\rVert^3)

The Fisher ellipse is therefore the local quadratic approximation to an exact equal-KL boundary. Choosing the point where the lowest reachable linearized free-energy contour touches that ellipse produces the natural-gradient update.

Interactive three-state variational inference

The distributions determine the gradient, Fisher ellipse, and selected update

Current logit

\phi_1

0.70Current logit

\phi_2

-0.40Local KL budget

\varepsilon

0.018 nats

Observed outcome o₂

Prior p(s)fixed generative-model prior

s1s2s3

[0.460, 0.340, 0.200]

Likelihood p(o₂|s)selected row of the observation model

s1s2s3

[0.180, 0.640, 0.180]

Target p(s|o)normalized likelihood × prior

s1s2s3

[0.246, 0.647, 0.107]

Current qφ(s)softmax(φ₁, φ₂, 0)

s1s2s3

[0.547, 0.182, 0.271]

F(q_\phi)=\sum_s q_\phi(s)\ln\frac{q_\phi(s)}{p(o\mid s)p(s)}=1.547\;\text{nats}

\nabla_\phi F=(0.186,\,-0.314)

I(\phi)=\begin{bmatrix}0.248&-0.099\\-0.099&0.149\end{bmatrix}

Scroll horizontally to compare the two rulers →

Each row below shows the new approximate posterior obtained by adding the displayed step to the current logits. Both steps use the same Fisher-predicted KL budget $\varepsilon$ . The “actual KL” is then recalculated directly between the current and updated distributions; its departure from $\varepsilon$ shows the error in the local quadratic approximation.

Ordinary-updated q(s)qφ → qφ+d actual KL = 0.0192 nats

s1s2s3

[0.464, 0.257, 0.279]

Natural-gradient-updated q(s)qφ → qφ+d actual KL = 0.0199 nats

s1s2s3

[0.495, 0.266, 0.238]

Nothing in the geometry is chosen by eye. The logits produce the current three-state distribution

q_\phi=\operatorname{softmax}(\phi_1,\phi_2,0)

. The selected observation combines with the prior to produce the target posterior. Here

F

is the variational free energy for that observation; minimizing it moves

q_\phi

toward the target posterior. These distributions determine both

\nabla_\phi F

and

I(\phi)

. The circle contains parameter steps with the same coordinate length as the displayed ordinary step. The Fisher ellipse contains parameter steps with the same predicted belief change

\varepsilon

. The ordinary and natural gradients select the best downhill direction under those two definitions of equal-sized movement. The solid gold boundary shows where the KL change is exactly

\varepsilon

; its agreement with the dashed Fisher boundary shows the local approximation at work. The parallel lines are equal predicted-free-energy contours, and the highlighted endpoint is where the lowest reachable contour touches each boundary.

Variational inference contains two nested questions. The outer problem is to find the coordinates $\phi$ that minimize $F(\phi)$ . At the current coordinates, the inner step-selection problem asks which nearby displacement $d=\Delta\phi$ should be taken next. Write the current free-energy gradient as $g=\nabla_\phi F$ . For a small candidate step, the predicted change in free energy is $g^{\mathsf T}d$ .

The question becomes: among the steps allowed to change the belief distribution by this amount, which one produces the greatest local decrease in free energy? This is written compactly as:

MOST FREE-ENERGY DECREASE PER UNIT OF BELIEF CHANGE

\min_{\Delta\phi}\;\bigl(\nabla_\phi F\bigr)^{\mathsf T}\Delta\phi \quad\text{subject to}\quad \frac12\Delta\phi^{\mathsf T}I(\phi)\Delta\phi\leq\varepsilon \qquad\Longrightarrow\qquad \Delta\phi\propto-I(\phi)^{-1}\nabla_\phi F

Read the expression from left to right. The symbol $\min_{\Delta\phi}$ means “search over the possible candidate steps $\Delta\phi$ and choose the one that makes the following quantity as small as possible.” Each part then has a direct interpretation:

F(\phi+\Delta\phi)-F(\phi)\approx(\nabla_\phi F)^{\mathsf T}\Delta\phi

This dot product predicts the first-order change in free energy. Making it as negative as possible gives the greatest local decrease.

\frac12\Delta\phi^{\mathsf T}I(\phi)\Delta\phi\approx D_{\mathrm{KL}}(q_\phi\|q_{\phi+\Delta\phi})

This quadratic form predicts how much the entire belief distribution changes under the candidate step.

\frac12\Delta\phi^{\mathsf T}I(\phi)\Delta\phi\leq\varepsilon

The constraint rejects steps whose predicted belief change exceeds the small budget ε. Unless the gradient is already zero, the best step uses the available budget.

\Delta\phi\propto-I(\phi)^{-1}\nabla_\phi F

The solution has this direction: the inverse Fisher transforms the ordinary gradient before its sign is reversed to point downhill. The next calculation shows why.

The proportionality above states the direction of the answer, but it does not yet show how the constraint produces that direction or determine the step's length. The obstacle is that the same candidate $d$ appears in two requirements at once: we want to make $g^{\mathsf T}d$ as negative as possible, while keeping $\tfrac12d^{\mathsf T}I(\phi)d$ within the budget $\varepsilon$ . A standard way to solve an optimization problem with such a constraint is to introduce one temporary number called a Lagrange multiplier, written $\lambda$ . It combines the objective and constraint into a single expression that can be differentiated. $\lambda$ is not a parameter of the probability model and is not chosen independently; solving the equations determines the value that places the optimal step on the budget boundary.

\mathcal L(d,\lambda)=g^{\mathsf T}d+\lambda\!\left(\frac12d^{\mathsf T}I(\phi)d-\varepsilon\right)

The first term rewards a predicted decrease in free energy. The expression in parentheses is zero when the KL budget is used exactly.

\nabla_d\mathcal L=g+\lambda I(\phi)d=0

At the optimum, no infinitesimal change in d can improve the combined objective. Differentiate with respect to the candidate step d and set the result to zero.

\lambda I(\phi)d=-g\qquad\Longrightarrow\qquad d=-\frac1\lambda I(\phi)^{-1}g

Move g to the other side, then undo multiplication by the Fisher matrix. This reveals the natural-gradient direction before its length has been determined.

\varepsilon=\frac12d^{\mathsf T}I(\phi)d=\frac1{2\lambda^2}g^{\mathsf T}I(\phi)^{-1}g

Substitute the proposed step back into the constraint. The remaining unknown λ controls how far to travel along the natural-gradient direction.

\frac1\lambda=\sqrt{\frac{2\varepsilon}{g^{\mathsf T}I(\phi)^{-1}g}}

Solve the budget equation for the scale 1/λ. A larger permitted KL change produces a larger step.

d_{\mathrm{natural}}^*=-\sqrt{\frac{2\varepsilon}{g^{\mathsf T}I(\phi)^{-1}g}}\,I(\phi)^{-1}g

Insert that scale into the natural-gradient direction. The star marks the candidate step that solves this constrained local problem.

The calculation has now produced three related but distinct objects. The ordinary gradient is $g=\nabla_\phi F$ . Applying the inverse Fisher matrix gives the natural gradient $\widetilde{\nabla}_\phi F=I(\phi)^{-1}g$ ; like an ordinary gradient, it points toward increasing free energy, but its direction and scale are defined using the Fisher geometry. Reversing it gives the downhill natural-gradient direction $-I(\phi)^{-1}g$ . Finally, $d_{\mathrm{natural}}^*$ is the particular downhill update obtained by scaling that direction until it uses exactly the chosen local KL budget $\varepsilon$ . Thus the natural gradient supplies the direction, while $d_{\mathrm{natural}}^*$ is the budget-sized step actually selected by this constrained problem.

In the spending analogy, $\varepsilon$ is the available budget and $\lambda$ is the price of using it. A large $\lambda$ makes $1/\lambda$ small and suppresses the step; a small $\lambda$ permits a larger step. If the update is instead written $d=-\rho I(\phi)^{-1}g$ , then the learning rate determines the locally predicted expenditure:

\varepsilon\approx\frac12\rho^2g^{\mathsf T}I(\phi)^{-1}g

The goal is therefore not to find the step that changes the beliefs as little as possible. A zero step would win that trivial contest while learning nothing. The goal is to obtain the greatest predicted reduction in free energy for a fixed small budget of belief change.

Rather than solving for the precise scale that exhausts a prescribed $\varepsilon$ at every iteration, an algorithm may use a positive learning rate $\rho$ to scale the downhill natural-gradient direction. The resulting update can be written in two steps:

NATURAL-GRADIENT STEP

\widetilde{\nabla}_\phi F=I(\phi)^{-1}\nabla_\phi F,\qquad \Delta\phi=-\rho\widetilde{\nabla}_\phi F=-\rho I(\phi)^{-1}\nabla_\phi F

With one parameter, $I(\phi)^{-1}$ is a single positive number. It can only rescale the ordinary gradient: large Fisher information means that a small coordinate movement would cause a large belief change, so its small inverse suppresses the step; small Fisher information has a larger inverse and permits a larger coordinate movement. In one dimension it cannot steer left or right independently of the gradient.

With several parameters, $I(\phi)^{-1}$ is a matrix. It can reweight and mix the gradient components, so the natural gradient can point in a different coordinate direction from the ordinary gradient. Imagine ordinary gradient descent as a skater dropping straight down a ramp according to the ramp's coordinate grid. Fisher geometry recognizes that equal movements in different directions on that grid may produce very unequal changes in the represented beliefs. The natural gradient redirects the local descent toward the greatest drop in free energy allowed by the belief-change budget. Because the gradient and Fisher matrix are recalculated along the way, successive steps may trace a curve across the ramp rather than the ordinary straight-down path. The path can approach a minimum from a different side; the geometry changes the route, not the objective being minimized.

What changes if we use exact KL?

The natural-gradient direction above came from replacing the exact KL boundary with its local Fisher approximation. We can instead retain the exact KL divergence in the step-selection problem:

EXACT-KL TRUST REGION

d_{\mathrm{exact}}^*=\underset{d}{\operatorname{argmin}}\;g^{\mathsf T}d \quad\text{subject to}\quad D_{\mathrm{KL}}\!\left(q_\phi\middle\|q_{\phi+d}\right)\leq\varepsilon

Here “best” means the best nearby endpoint $q_{\phi+d}$ according to the local prediction $g^{\mathsf T}d$ —not the globally optimal posterior found in one step. Applying the same Lagrange multiplier method now gives:

\mathcal L_{\mathrm{exact}}(d,\lambda)=g^{\mathsf T}d+\lambda\!\left[D_{\mathrm{KL}}(q_\phi\|q_{\phi+d})-\varepsilon\right]

The exact, curved KL boundary now replaces the Fisher ellipse.

g+\lambda\nabla_dD_{\mathrm{KL}}(q_\phi\|q_{\phi+d})=0

The derivative of exact KL depends nonlinearly on d through the updated distribution qφ+d, so this equation usually requires an iterative numerical solution.

\nabla_dD_{\mathrm{KL}}(q_\phi\|q_{\phi+d})=I(\phi)d+O(\lVert d\rVert^2)

For a sufficiently small step, the derivative of exact KL is approximated by I(φ)d. The O(‖d‖²) term collects the omitted contributions whose magnitude shrinks at least quadratically as d approaches zero.

g+\lambda I(\phi)d\approx0\qquad\Longrightarrow\qquad d\propto-I(\phi)^{-1}g

With that local replacement, the nonlinear exact-KL condition becomes the solvable natural-gradient condition.

In the three-state visualization, exact KL is easy to evaluate for any proposed step; that is how the solid gold boundary and the reported “actual KL” are calculated. Evaluation and optimization are different tasks, however. Evaluating one candidate tells us its KL cost, while finding the best candidate anywhere on the curved boundary requires solving the constrained problem. The Fisher approximation replaces that boundary locally with an ellipse, which is why the natural-gradient solution has a closed form.

Two approximations at two different levels

Variational inference uses $q_\phi$ to approximate a posterior that may be difficult to compute or represent exactly. Natural-gradient variational inference keeps that same approximate posterior family and the same free-energy objective. Its additional approximation concerns only the next optimization step: it uses the Fisher quadratic to approximate the KL change from $q_\phi$ to a nearby $q_{\phi+d}$ . It is therefore a method for navigating the variational family, not a second target for inference.

An implementation normally does not construct a Taylor polynomial at every iteration or explicitly calculate a matrix inverse. It computes the ordinary gradient $g$ and Fisher matrix $I$ , solves the linear system $Iv=g$ for $v$ , and then updates in the direction $-v$ . A learning rate may set the step directly, or the budget formula above may scale it to a chosen $\varepsilon$ :

IMPLEMENTABLE NATURAL-GRADIENT STEP

I(\phi)v=g,\qquad d=-\rho v\qquad\text{or}\qquad d=-\sqrt{\frac{2\varepsilon}{g^{\mathsf T}v}}\,v

The three-state categorical example

Consider one observation $o$ and one hidden state with three possible values. Write the approximate posterior as $q=(q_1,q_2,q_3)$ . The likelihood times the prior is simply the joint probability $p(o,s)=p(o\mid s)p(s)$ , so its variational free energy is:

THE SAME FREE ENERGY USED BY THE TOY MODEL

F[q]=\sum_s q_s\!\left(\ln q_s-\ln p(o\mid s)-\ln p(s)\right)=\sum_s q_s\!\left(\ln q_s-\ln p(o,s)\right)

Differentiating directly with respect to each probability gives $g_s=\partial F/\partial q_s=\ln q_s+1-\ln p(o,s)$ . An ordinary additive step proposes $v=q-\rho g$ . Because this calculation treats the probabilities as unrestricted coordinates, its proposal need not remain a valid probability distribution: an entry may become negative, or the entries may no longer sum to one.

Aside

Three ways to keep a probability update on the simplex

The probability simplex is the set of valid probability vectors: $\Delta=\{q:q_i\geq0\text{ for every }i,\ \sum_iq_i=1\}$ . Begin with $q=(0.05,0.35,0.60)$ , gradient $g=(0.7,-0.4,-0.3)$ , and learning rate $\rho=0.1$ . The ordinary additive proposal is:

THE UNCONSTRAINED PROPOSAL

v=q-\rho g=(0.05,0.35,0.60)-0.1(0.7,-0.4,-0.3)=(-0.02,0.39,0.63)

The first proposed probability is $-0.02$ , so $v$ lies outside the simplex. The following methods answer the same practical question in different ways: how should this proposed movement be converted into a valid probability distribution?

1. Clip to zero, then renormalize

First replace each negative entry by zero. The surviving entries then sum to $1.02$ , so divide every entry by $1.02$ to restore a total of one:

v=(-0.02,0.39,0.63)\longrightarrow v_{\mathrm{clipped}}=(0,0.39,0.63)

The impossible negative probability is replaced by zero.

q_{\mathrm{new}}=\frac{v_{\mathrm{clipped}}}{\sum_i(v_{\mathrm{clipped}})_i}=\frac{(0,0.39,0.63)}{1.02}\approx(0,0.382,0.618)

Dividing by the new total rescales the surviving entries so they sum to one.

This is the repair used by the toy implementation. It is simple, but the clipping changes the proposed direction, and an exact zero is awkward when later calculations contain $\ln q_i$ . Adding a tiny positive constant avoids evaluating $\ln0$ , but does not change the fact that this is an after-the-step repair.

2. Project to the nearest valid probability vector

Euclidean projection searches the entire simplex and chooses the valid vector closest to the invalid proposal $v$ . “Closest” means the smallest sum of squared coordinate differences:

EUCLIDEAN PROJECTION ONTO THE SIMPLEX

q^*=\underset{q\in\Delta}{\operatorname{argmin}}\;\frac12\lVert q-v\rVert^2=\underset{q\in\Delta}{\operatorname{argmin}}\;\frac12\sum_i(q_i-v_i)^2

Here $q$ ranges over candidate valid distributions. The expression after $\operatorname{argmin}$ assigns each candidate a squared distance from $v$ . Unlike $\min$ , which returns the smallest distance, $\operatorname{argmin}$ returns the candidate $q$ at which that smallest distance occurs. The star labels the selected solution: $q^*$ is the whole projected vector, while $q_i^*$ is its $i$ th entry. Read the equation as: “Among every valid $q$ in $\Delta$ , return the $q$ that makes $\tfrac12\lVert q-v\rVert^2$ smallest.”

The solution subtracts one shared threshold $\tau$ from the entries that can remain positive and sets the rest to zero. The threshold is chosen so the result sums to one:

THE PROJECTION RULE

q_i^*=\max(v_i-\tau,0),\qquad\sum_iq_i^*=1

For $v=(-0.02,0.39,0.63)$ , the first entry remains zero. The other two exceed a total of one by $0.02$ , so the closest solution subtracts $\tau=0.01$ from each of them:

q^*=\bigl(0,0.39-0.01,0.63-0.01\bigr)=(0,0.38,0.62)

Projection therefore changes the coordinates additively by the shared threshold. Clip-and-renormalize instead rescales the surviving coordinates multiplicatively, which is why it produced the slightly different vector $(0,0.382,0.618)$ .

3. Use an exponentiated-gradient update

Exponentiated gradient avoids leaving the simplex rather than repairing an invalid additive proposal. It multiplies each old probability by a positive factor $e^{-\rho g_i}$ and then normalizes:

EXPONENTIATED-GRADIENT UPDATE

q_i^{\mathrm{new}}=\frac{q_i e^{-\rho g_i}}{\sum_jq_j e^{-\rho g_j}}=\operatorname{softmax}(\ln q-\rho g)_i

w_1=0.05e^{-0.07}\approx0.0466

The positive first gradient decreases the first category’s weight.

w_2=0.35e^{0.04}\approx0.3643

The negative second gradient increases the second category’s weight.

w_3=0.60e^{0.03}\approx0.6183

The negative third gradient increases the third category’s weight.

q_{\mathrm{new}}=\frac{(w_1,w_2,w_3)}{w_1+w_2+w_3}\approx(0.045,0.354,0.601)

Every weight is positive; normalizing them returns a valid probability vector.

In log-probability coordinates, multiplication becomes addition: $\ln q_i^{\mathrm{new}}=\ln q_i-\rho g_i-\ln Z$ . This is why the same update can be written as a softmax. It is also called a multiplicative-weights or entropic mirror-descent update: it limits relative change using KL geometry rather than finding the nearest repair by ordinary Euclidean distance.

Simplex projection: J. Duchi et al. (2008), “Efficient Projections onto the ℓ₁-Ball for Learning in High Dimensions”; exponentiated gradient: J. Kivinen and M. Warmuth (1997), “Exponentiated Gradient versus Gradient Descent”.

The toy model's clip-and-renormalize method is therefore one possible way to keep direct probability updates valid. In the tiny one-factor model it serves as a didactic optimization heuristic, but the repair is not the same as solving the constrained gradient problem and need not converge to the exact posterior. That limitation motivates the projection, exponentiated-gradient, and natural-gradient alternatives.

Natural-gradient variational inference takes another route: represent the same categorical beliefs by unconstrained logits $\eta$ and use the distribution's Fisher geometry to scale their update. Softmax converts those logits into probabilities:

CATEGORICAL BELIEFS IN LOGIT COORDINATES

q_s=\operatorname{softmax}(\eta)_s=\frac{e^{\eta_s}}{\sum_j e^{\eta_j}}

The Fisher matrix for these categorical beliefs follows from the score definition already developed earlier. If the sampled category is $S$ , let $e_S$ denote its one-hot vector: for example, $S=2$ gives $e_S=(0,1,0)$ . Differentiating the categorical log-probability with respect to every logit gives the score vector $e_S-q$ . Averaging its outer product produces:

\ln q_\eta(S)=\eta_S-\ln\sum_j e^{\eta_j}

Softmax assigns the observed category its logit and subtracts the shared log-normalizer.

\frac{\partial}{\partial\eta_j}\ln q_\eta(S)=\mathbf 1\{S=j\}-q_j

The indicator is one for the observed category and zero for every other category; differentiating the normalizer contributes -q_j.

\nabla_\eta\ln q_\eta(S)=e_S-q

Collecting those derivatives gives the categorical score vector.

I(\eta)=\mathbb E[(e_S-q)(e_S-q)^{\mathsf T}]

Fisher information is the expected outer product of the score vector.

I(\eta)=\mathbb E[e_Se_S^{\mathsf T}]-\mathbb E[e_S]q^{\mathsf T}-q\mathbb E[e_S]^{\mathsf T}+qq^{\mathsf T}

Expand the outer product before taking its expectation.

\mathbb E[e_S]=q,\qquad \mathbb E[e_Se_S^{\mathsf T}]=\operatorname{diag}(q)

The mean one-hot vector is the probability vector, while its outer product places one only on the sampled category’s diagonal.

I(\eta)=\operatorname{diag}(q)-qq^{\mathsf T}-qq^{\mathsf T}+qq^{\mathsf T}=\operatorname{diag}(q)-qq^{\mathsf T}

Substitute those two expectations and combine the three outer-product terms.

This matrix has diagonal entries $q_i(1-q_i)$ and off-diagonal entries $-q_iq_j$ . It depends on the current distribution $q$ : as the logits move and the probabilities change, the local Fisher geometry changes with them.

The exact Bayesian posterior is the softmax of the prior log-probability plus the log-likelihood. Define its target logits by $\eta_s^*=\ln p(s)+\ln p(o\mid s)=\ln p(o,s)$ . The formula for the ordinary logit gradient can now be derived rather than introduced as a new identity. Exponentiating and normalizing those target logits gives Bayes' rule directly:

TARGET LOGITS PRODUCE THE EXACT POSTERIOR

\operatorname{softmax}(\eta^*)_s=\frac{e^{\ln p(o,s)}}{\sum_j e^{\ln p(o,j)}}=\frac{p(o,s)}{p(o)}=p(s\mid o)

The derivative of free energy with respect to the probabilities is $\nabla_qF=\ln q+\mathbf1-\eta^*$ . Next, because $q$ depends on $\eta$ , the chain rule passes this derivative through the softmax Jacobian. That Jacobian is the same matrix $\operatorname{diag}(q)-qq^{\mathsf T}=I(\eta)$ derived above:

\frac{\partial q_i}{\partial\eta_j}=q_i\left(\mathbf1\{i=j\}-q_j\right)

Differentiating softmax shows how moving logit j changes probability i.

\frac{\partial q}{\partial\eta}=\operatorname{diag}(q)-qq^{\mathsf T}=I(\eta)

Collecting those derivatives produces the softmax Jacobian, which is exactly the categorical Fisher matrix derived above.

\nabla_\eta F=\left(\frac{\partial q}{\partial\eta}\right)^{\mathsf T}\nabla_qF

The chain rule converts a derivative with respect to probabilities into one with respect to the logits that produce them.

\nabla_\eta F=I(\eta)\left(\ln q+\mathbf1-\eta^*\right)

Substitute the Fisher matrix for the symmetric softmax Jacobian and the previously derived probability gradient.

\ln q=\eta-c\mathbf1,\qquad c=\ln\sum_j e^{\eta_j}

Every log-probability is its logit minus the same log-normalizer c.

I(\eta)\mathbf1=0\qquad\Longrightarrow\qquad\nabla_\eta F=I(\eta)(\eta-\eta^*)

The Fisher matrix removes terms shared equally by every logit, leaving only the meaningful relative differences.

The last cancellation follows directly from normalization: $\operatorname{diag}(q)\mathbf1=q$ and $qq^{\mathsf T}\mathbf1=q\sum_iq_i=q$ , so their difference is zero. The ordinary logit gradient is therefore:

ORDINARY LOGIT GRADIENT

\nabla_\eta F=I(\eta)(\eta-\eta^*)

Multiplying by the inverse Fisher removes this coordinate-induced scaling. Strictly, the categorical Fisher matrix has one redundant direction because adding the same constant to every logit leaves the softmax unchanged; fixing one reference logit or using a pseudoinverse gives the same probability update. Up to that irrelevant shared constant:

NATURAL LOGIT UPDATE

\widetilde{\nabla}_\eta F=I(\eta)^{+}\nabla_\eta F=\eta-\eta^*, \qquad \eta_{\mathrm{new}}=(1-\rho)\eta+\rho\eta^*

Ordinary descent in

q

Move, clip, renormalize

The numerical step depends on probability coordinates and may leave the simplex before being repaired.

Natural descent in

\eta

Interpolate toward Bayesian evidence

The softmax always produces a valid distribution, while Fisher scaling measures the step by its effect on beliefs.

Why this example is exceptionally clean

When

\rho=1

, the new logits equal

\eta^*

, so their softmax is exactly

q^*(s)\propto p(o\mid s)p(s)

. A full natural-gradient step therefore reproduces the exact Bayesian update for this single categorical factor. This one-step result comes from the simple categorical exponential-family structure; natural-gradient variational inference does not generally solve every model in one step.

What changes under a mean-field approximation?

For several hidden-state factors, mean field restricts the approximate posterior to a product family:

MEAN-FIELD FAMILY

q_\phi(s_1,\ldots,s_m)=\prod_{i=1}^{m}q_{\phi_i}(s_i)

Because the score of factor $i$ depends only on that factor, score covariances between different factors vanish under the product distribution. The Fisher matrix is therefore block diagonal, with one block for each factor:

PRODUCT-FAMILY FISHER GEOMETRY

I_{\mathrm{MF}}(\phi)=\operatorname{blockdiag}\!\left(I_1(\phi_1),\ldots,I_m(\phi_m)\right)

In a conditionally conjugate exponential-family model, holding every other factor fixed produces an analytic target for factor $i$ . In categorical logit coordinates it has the form:

ONE MEAN-FIELD COORDINATE UPDATE

\eta_i^*(s_i)=\mathbb E_{q_{-i}}\!\left[\ln p(o,s_1,\ldots,s_m)\right]+\text{constant}, \qquad \eta_i\leftarrow(1-\rho)\eta_i+\rho\eta_i^*

With $\rho=1$ , this is the familiar closed-form coordinate-ascent mean-field update. It can also be read as a full natural-gradient step for that factor. Thus conjugate mean-field coordinate ascent already uses the factor's Fisher geometry implicitly; damping with $0<\rho<1$ takes only part of that natural step.

What the product geometry can—and cannot—do

The block-diagonal Fisher matrix makes each retained factor update statistically well scaled. It cannot restore correlations excluded by

q(s_1,\ldots,s_m)=\prod_iq_i(s_i)

. Capturing those dependencies requires a full joint posterior or a structured mean-field family. In nonconjugate, stochastic, or amortized inference, the analytic coordinate target may also disappear; then the natural gradient must be computed or approximated explicitly.

Natural gradients and variational inference: S. Amari (1998), “Natural Gradient Works Efficiently in Learning”; M. Sato (2001), “Online Model Selection Based on the MDL Principle”; M. Hoffman et al. (2013), “Stochastic Variational Inference”.

Outline: the path through this tutorial

Recommended related tutorials

Sources and further reading

The bridge: distributions become points

Why define information like this?

A concrete coin model

Encode coin P using a code designed for coin Q

The price of the wrong codebook

The coin model: 50/50 reality versus a 75/25 model

The local picture: the bottom of a KL bowl

A clean example: known Gaussian noise, unknown mean

From a noisy observation to evidence about the mean

One noisy observation gives a clue; independent observations sharpen it

Why can we rely on the score being centered at zero?

From balanced scores to Fisher information

From one parameter to a parameter space

Insert the nearby model into the cross-entropy definition

Taylor-expand the nearby model's log-density

The linear term vanishes because every distribution normalizes to one

The remaining Hessian is the Fisher information matrix

From local KL to a geometric length

What “zooming in” and “ruler” mean

From points to paths

Two routes between the same distributions

Where KL asymmetry enters

Fisher information: zoom into KL

How distinguishable are two nearby coins?

Fisher information versus the Fisher information metric

Bonus: natural-gradient variational inference

From a free-energy gradient to a belief update

The distributions determine the gradient, Fisher ellipse, and selected update

What changes if we use exact KL?

The three-state categorical example

1. Clip to zero, then renormalize

2. Project to the nearest valid probability vector

3. Use an exponentiated-gradient update

What changes under a mean-field approximation?