Korean ver.

Interactive Probability

Date: |Estimated Reading Time: 52 min|Author: Seungheon Doh

Reference: Goodfellow, Bengio & Courville, Deep Learning — Chapter 3: Probability and Information Theory. Information theory (entropy, KL, cross-entropy) is covered in a separate article.

In the earlier linear algebra article, data was represented by vectors and matrices, and a model was a linear map transforming those representations. But that language alone cannot explain why feeding the same input to the same model twice can produce a different answer. A language model generates a different sentence each time from the same prompt; a diffusion model produces a different image each time under the same condition. The reason is that we are handling distributions, not deterministic functions.

There are three broad reasons machine learning needs probability. First, uncertainty always remains in a prediction — because the world is inherently stochastic (inherent randomness), because our observations are incomplete (partial observability), or because the model is a simplification (model misspecification). Second, only by quantifying that uncertainty can we judge how much to trust a prediction. Third, if we learn the distribution that generates the data itself, we can produce new samples — which is exactly what a generative model does.

Probability even determines the shape of the loss function. Cross-entropy for classification, MSE for regression, MAE for robust regression, and L2/L1 regularization each follow as the negative log-likelihood of a different probabilistic model. Instead of memorizing loss functions, ask "what distribution is this model assuming?" and the formula follows.

This article starts from the definition of a random variable and moves through PMF/PDF and sampling, joint, marginal and conditional probability, the chain rule and independence, expectation, variance and covariance, and on to maximum likelihood estimation, Bayes' rule and MAP. Rather than listing distributions like a catalogue, it treats in depth only the three that touch machine learning losses directly — the Categorical distribution of classification, the Gaussian of regression, and its alternative the Laplace. Each concept is presented as definition, intuition, algebraic form, and use in machine learning.

Random Variables

A random variable is the device that translates "what could happen" into numbers. A die roll, a classifier label, a person's height — objects of entirely different kinds are all handled through this single format.

Definition. A random variable \(X\) is a function \(X:\Omega\to\mathbb{R}\) assigning a real number to each outcome \(\omega\) in the sample space \(\Omega\). Despite the name, a random variable is not itself random — it is a function; what is random is which \(\omega\) gets drawn. By convention the variable is written with a capital \(X\) and the value it actually takes with a lowercase \(x\).

Why we need this device. "The die shows an even number," "this photo is a cat," and "this person is about 170cm tall" are statements of different kinds. A random variable turns all of them into statements about numbers. Once something is a number we can add it, average it, and differentiate it. Machine learning can work with probability because of this translation.

Discrete and continuous. If the values \(X\) can take are finite or countable, it is a discrete random variable; if they fill a real interval, it is continuous. This distinction is what separates the PMF from the PDF in the next section. In the diagram, the die and the class label are drawn as bars while height is drawn as a smooth curve — that is the difference.

In machine learning. Almost everything in training is a random variable: the input features \(X\), the label \(Y\), unobserved latent variables \(Z\) (a VAE's code, a mixture component index), and the action \(A\) and reward \(R\) of reinforcement learning. Building a model means specifying the joint or conditional distribution among them.

\[\begin{aligned} X &: \Omega \to \mathbb{R} \\ P(X = x) &= P(\{\omega \in \Omega : X(\omega) = x\}) \end{aligned}\]

Press the three buttons. You can compare how an outcome of the sample space (a die face, an image) is carried over to a single number, and why for a continuous variable it is area rather than bar height that becomes probability.

Probability Distributions

If a random variable decides what gets observed as a number, a probability distribution decides how often each number appears. A distribution describes both the set of possible values and how probability is spread across them.

Definition. The probability distribution of a random variable \(X\) is the rule assigning \(P(X\in A)\) to each set \(A\). Put simply, it collects the answers to every question of the form "how likely is \(X\) to land in this value or range?" For a discrete variable we can list the probability of each value; for a continuous one we describe the distribution by the probability of falling in an interval.

A random variable and a distribution are different things. The random variable \(X\) is the rule translating outcomes of the sample space into numbers; the distribution is the behavioural pattern describing how probability sits on the results of that translation. By analogy with a program, the random variable is the variable or function determining what gets returned, and the distribution is the statistical pattern its return values make when the code is run many times.

Different random variables can share a distribution. The variable \(X\) mapping heads to 1 and tails to 0, and the variable \(Y\) mapping exam pass to 1 and fail to 0, are defined on different experiments. But if the probability of heads and of passing are both 0.7, the two follow the very same \(\mathrm{Bernoulli}(0.7)\) distribution. A distribution discards the specific story behind the outcome and keeps only the structure of values and probabilities.

Variables of the same kind can have different distributions. Two dice may both take values 1 through 6, but if one is fair and the other favours 6, their distributions differ. Knowing the list of possible values is not enough; we must know the probability attached to each value before we can predict the variable's behaviour.

Notation. \(X\sim p\) or \(X\sim\mathcal D\) reads "\(X\) follows the distribution \(p\) (or \(\mathcal D\))." For \(X\sim\mathcal N(\mu,\sigma^2)\) the family is Gaussian, and \(\mu\) and \(\sigma^2\) are the parameters fixing its location and spread. A collection of distributions sharing one name and differing only in parameters is called a distribution family.

How a distribution is written down. A distribution is one concept, and the PMF, PDF and CDF are different tools for expressing it. A discrete variable writes the probability of each value directly with a PMF; a continuous variable expresses interval probabilities as areas under a PDF. The CDF can be used for both. The next section looks concretely at the differences among these three representations.

\[\begin{aligned} X\sim\mathcal D &\quad\Longleftrightarrow\quad P(X\in A)=\mathcal D(A) \\ X\sim\mathrm{Bernoulli}(0.7) &\quad\Rightarrow\quad P(X=1)=0.7 \\ X\sim\mathcal N(\mu,\sigma^2) &\quad\Rightarrow\quad \mu:\text{ center},\;\sigma^2:\text{ spread} \end{aligned}\]

PMF, PDF and Probability

A discrete random variable assigns probability directly to each value (PMF), but for a continuous variable the probability of any single point is zero, so we assign a density and integrate over an interval to get a probability (PDF). The two worlds are not disconnected: narrowing the bin width joins them naturally.

Definition (PMF). The probability mass function \(P(x)\) of a discrete random variable assigns probability directly to each value. There are only two conditions: \(0\le P(x)\le 1\) for every \(x\), and the total \(\sum_x P(x)=1\).

Definition (PDF). For a continuous random variable the probability of a single point is zero — the probability that a height is exactly 170.000…cm is zero. So instead of probability we assign a density \(p(x)\ge 0\) and obtain probability by integrating over an interval. It suffices that the total integral is \(\int p(x)\,dx=1\).

Density is not probability. This is the most common misconception. \(P(x)\) can never exceed 1, but \(p(x)\) may exceed 1 freely. The uniform distribution on \([0,\,0.1]\), for instance, has density 10. Density is probability per unit length, so it only becomes a probability once multiplied by a width \(\Delta\): \(P(x\le X\le x+\Delta)\approx p(x)\,\Delta\).

The two worlds are connected. Slice a continuous distribution into bins of width \(\Delta\) and collect the probability of each cell, and you get a discrete distribution. In the diagram, shrinking \(\Delta\) makes the bars ever finer until they cling to the curve — this is the limit taking a PMF to a PDF. Conversely, discretizing a continuous random variable into bins is always possible, and in practice image pixel values and audio samples are handled exactly that way.

The CDF: one expression binding both. The cumulative distribution function \(F(x)=P(X\le x)\) is defined identically whether the variable is discrete or continuous. It increases monotonically, going to 0 at \(-\infty\) and 1 at \(+\infty\). An interval probability comes out of a single subtraction: \(P(a\le X\le b)=F(b)-F(a)\). The sampling in the next section uses precisely this function.

\[\begin{aligned} \textstyle\sum_x P(x) &= 1 && \text{(PMF)} \\ \int_{-\infty}^{\infty} p(x)\,dx &= 1 && \text{(PDF)} \\ P(a\le X\le b) &= \int_a^b p(x)\,dx = F(b)-F(a) \end{aligned}\]

Shrink the Δ slider. The bars (a discretized probability) converge onto the curve (a density). Watching the "bar height × Δ = probability" relation in the readout makes the difference between density and probability plain. Use "Show CDF" to see the cumulative form as well.

Sampling

Knowing a distribution and being able to draw samples from it are different abilities. The way a histogram converges to the underlying density as samples accumulate is the fact on which Monte Carlo estimation, SGD minibatches, and generative sampling all rest.

Definition. To draw a sample from a distribution \(p\) is to generate values so that, in the long run, the relative frequency of each value matches \(p\). We write this \(x^{(i)}\sim p(x)\). Knowing the formula for a distribution and being able to sample from it are separate abilities, and without the latter a generative model cannot exist.

Samples reveal the distribution. A histogram of \(N\) samples approximates the underlying density. For small \(N\) it is lumpy, but it smooths out as \(N\) grows, and the gap between the sample mean and the true expectation shrinks at roughly the rate \(1/\sqrt{N}\). The slowness of this convergence matters: raising precision tenfold requires a hundred times as many samples.

How to draw: inverse transform sampling. Most pseudorandom generators supply uniform integers or \(\mathrm{Uniform}(0,1)\) numbers as raw material, and numbers from other distributions are made by transforming those. Draw a uniform \(u\) and compute \(x=F^{-1}(u)\), and that \(x\) follows \(p\) — because the CDF spreads the probability mass of the \(x\)-axis evenly onto \([0,1]\) on the \(y\)-axis. For distributions where \(F^{-1}\) is not easy to compute, other methods such as rejection sampling are used.

In machine learning. Replacing a hard-to-compute expectation with a sample mean is Monte Carlo estimation. The canonical example is SGD, which approximates the full loss by the average over a minibatch drawn uniformly from the dataset. Under such standard uniform-sampling conditions the minibatch gradient is an unbiased estimator of the full-data gradient. The process by which a diffusion model builds an image from noise, and a language model picks tokens one at a time, are both sampling too.

\[\begin{aligned} x^{(i)} &\sim p(x), \quad i=1,\ldots,N \\ \mathbb{E}_{p}[f(X)] &\approx \frac{1}{N}\sum_{i=1}^{N} f(x^{(i)}) \\ u &\sim \mathrm{Uniform}(0,1) \;\Rightarrow\; F^{-1}(u) \sim p \end{aligned}\]

Raise the N slider from 1 to 10,000. The histogram converges onto the orange density curve and the error in the readout shrinks. Press "Resample" and the result differs even at the same N — a sample always carries chance with it. The lower figure shows how a single uniform number becomes a sample from the desired distribution.

Inverse transform sampling: push a uniform u back through the CDF to get x

Joint Probability

Handle two or more random variables at once and probability sits on a grid rather than in a table. The joint distribution describes exhaustively what values the variables take together, and marginalization and conditioning are both just different ways of reading this grid.

Definition. The joint distribution \(P(X=x,\,Y=y)\) of two random variables \(X, Y\) is the probability that both take those values simultaneously. In the discrete case it becomes a table on a grid whose cells sum to 1: \(\sum_x\sum_y P(x,y)=1\). In the continuous case the double integral of the joint density \(p(x,y)\) equals 1.

The joint distribution holds everything. The answer to every question one can ask about \(X\) and \(Y\) is inside this table. The distribution of \(X\) alone, of \(Y\) alone, and of one given the other are all obtained by reading the grid differently. The next three sections — marginalization, conditioning, and the chain rule — are exactly those ways of reading.

The cost is size. If each variable takes \(K\) values, the grid has \(K^2\) cells, and with \(n\) variables it needs \(K^n\). For a language model with a vocabulary of 50,000, the joint distribution of just two words already runs to 2.5 billion cells. This exponential blow-up is the fundamental problem of probabilistic modelling, and independence assumptions, conditional independence, factorization and neural parameterization are all attempts to dodge it.

In machine learning. What a generative model is ultimately trying to learn is the joint distribution of the data \(p(x_1,\ldots,x_n)\). A latent-variable model inserts an unobserved \(Z\) to work with \(p(x,z)\), then sums \(z\) away when only the observed part is needed. That operation is the marginalization of the next section.

\[\begin{aligned} P(x,y) &\ge 0 \\ \textstyle\sum_x\sum_y P(x,y) &= 1 \\ \iint p(x,y)\,dx\,dy &= 1 \end{aligned}\]

Drag cells up and down to change their values (arrow keys work too). Raising one cell automatically renormalizes the rest so the total stays 1. Use the preset buttons to compare "correlated / independent / uniform" — this grid keeps reappearing over the next four sections.

Marginalization

Summing an unwanted variable away in a joint distribution leaves the distribution of what remains. This operation — erasing latent variables to obtain the distribution of the observed ones — is where the computational burden of latent-variable models and Bayesian inference arises.

Definition (sum rule). Summing a joint distribution over all values of one variable leaves the distribution of the other. The \(P(x)\) so obtained is called the marginal distribution, and the operation marginalization. In the continuous case the sum becomes an integral; the logic is unchanged.

Why "marginal"? The name comes from the practice of drawing a grid on paper and writing the sum of each row and column in the margin. In the diagram, summing over \(Y\) leaves each column's total in the bottom margin as \(P(X)\), and summing over \(X\) leaves each row's total in the side margin as \(P(Y)\).

It is an information-destroying operation. Marginalization cannot be undone. Knowing \(P(x)\) and \(P(y)\) does not recover \(P(x,y)\) — there are infinitely many joint distributions with the same marginals (in the "independence" demo of the next section, moving \(t\) leaves the marginals unchanged, which is an example). Information about the variable summed away is gone.

In machine learning: this is the bottleneck. A latent-variable model is defined by \(p(x)=\int p(x,z)\,dz\), and for complex models this integral is often intractable to compute exactly. That is why a VAE optimizes a lower bound called the ELBO instead, and why Bayesian inference depends on MCMC or variational inference because of the evidence \(p(x)=\int p(x\mid\theta)p(\theta)\,d\theta\). Diffusion models likewise introduce many intermediate latent states to convert an intractable data likelihood into stagewise probabilities and a computable training objective.

\[\begin{aligned} P(x) &= \sum_y P(x,y) \\ p(x) &= \int p(x,y)\,dy \\ p(x) &= \int p(x,z)\,dz \quad\text{(marginalizing out the latent)} \end{aligned}\]

Use the buttons to change the axis being summed over. The bars in the readout show the sum along that direction, and the total is always 1. Drag the grid cells to change the joint distribution and watch how the marginals follow.

Conditional Probability

A conditional probability is one slice taken out of the joint distribution and renormalized so that it sums to 1 again. This is precisely the distribution supervised learning learns, and it takes the form of an output distribution given an input.

Definition. When \(P(X=x)>0\), the conditional probability of \(Y\) given \(X=x\) is the joint probability divided by the probability of the condition: \(P(y\mid x)=P(x,y)/P(x)\). What the denominator does is normalize — it makes the detached slice sum to 1 again.

Geometrically: cut, then measure again. From the grid, keep only the column where \(X=x\) and erase the rest. The values in that column are proportions relative to the whole grid, so they do not sum to 1. Dividing by the column total \(P(x)\) turns them into proportions within that column, and now they sum to 1. To condition is to narrow the sample space and then measure probability afresh inside it.

Change the condition and the distribution changes. Moving the \(X\) slider changes the shape of \(P(Y\mid X=x)\). That change is exactly the information \(X\) carries about \(Y\). If the shape were identical for every \(x\) you chose, \(X\) would tell us nothing about \(Y\) — and that is the independence of two sections from here.

In machine learning: this is what supervised learning learns. A classifier learns \(p(y\mid x)\), and so does a regression model (usually as a Gaussian). Image captioning is \(p(\text{text}\mid\text{image})\), and conditional generation is \(p(x\mid c)\). What a neural network does, instead of memorizing the whole joint grid, is learn a function that takes an input \(x\) and directly outputs the normalized distribution of that column — which is why it need not store \(K^n\) cells.

\[\begin{aligned} P(y\mid x) &= \frac{P(x,y)}{P(x)}, \qquad P(x)>0 \\ \textstyle\sum_y P(y\mid x) &= 1 \quad\text{(for fixed } x \text{)} \end{aligned}\]

Change the condition with the X slider. The selected column is highlighted and the rest dimmed, and the readout shows that column normalized into \(P(Y\mid X=x)\) as bars. Check that the total always returns to 1.

The Chain Rule

A joint distribution can always be decomposed into a product of conditional probabilities, and the order of that decomposition is not unique. The fact that any order recovers the same joint distribution is why autoregressive models and next-token prediction in language models are legitimate.

Definition. Multiplying both sides of the definition \(P(y\mid x)=P(x,y)/P(x)\) by \(P(x)\) and turning it around gives \(P(x,y)=P(x)P(y\mid x)\). With \(n\) variables, repeating this decomposes the joint distribution completely into a product of conditional probabilities. This is the chain rule.

We choose the order. \(P(x,y)=P(x)P(y\mid x)\) and \(P(x,y)=P(y)P(x\mid y)\) hold simultaneously. The two decompositions store entirely different numbers yet multiply out to the same joint distribution. With \(n\) variables, all \(n!\) orderings are valid. No order is wrong; it is only that some are easy to compute and others hard.

Why this is decisive. The chain rule splits a joint distribution into \(n\) conditional distributions. However, storing conditional probability tables outright for every possible context still requires an exponential number of values. What lets a language model avoid an enormous table is the chain rule plus a single neural network — shared across every position and context — that computes "the distribution of the next item given what has been seen so far." The chain rule supplies the exact decomposition, and neural parameterization compresses and approximates it.

In machine learning: this is a language model. Decomposing the probability of a sentence \(P(w_1,\ldots,w_T)\) from left to right gives \(\prod_t P(w_t\mid w_{<t})\). During training the model is only ever asked to predict "the distribution of the next token," and thanks to the chain rule that alone amounts to representing the joint distribution of the whole sentence. PixelCNN, WaveNet and autoregressive diffusion models share the same skeleton, differing only in the axis.

\[\begin{aligned} P(x,y) &= P(x)\,P(y\mid x) = P(y)\,P(x\mid y) \\ P(x_1,\ldots,x_n) &= \prod_{i=1}^{n} P(x_i \mid x_1,\ldots,x_{i-1}) \\ \log P(w_{1:T}) &= \sum_{t=1}^{T} \log P(w_t \mid w_{<t}) \end{aligned}\]

Switch the decomposition-order buttons. The readout shows the numbers each decomposition actually stores (one marginal plus one conditional table) — the numbers differ but the products agree. Click the token strip below and the autoregressive decomposition advances one token at a time, accumulating log-probability.

Autoregressive decomposition: multiply token by token, growing the condition

Independence

If two variables are independent, the joint distribution factorizes completely into the product of the individual distributions. Instead of remembering the whole grid we need only the per-axis distributions, so an independence assumption dramatically reduces the number of parameters at the cost of expressiveness.

Definition. \(X\) and \(Y\) are independent when \(P(x,y)=P(x)P(y)\) for all \(x,y\), written \(X\perp Y\). Rewritten conditionally, \(P(y\mid x)=P(y)\) — learning \(x\) changes our belief about \(y\) not at all.

Geometrically: the grid becomes an outer product. An independent joint distribution is the outer product of a column vector and a row vector. Every row is therefore a constant multiple of every other, and likewise for the columns. Sending \(t\) to 1 in the demo shows the grid settling into exactly that shape. The marginals do not change at all in the process — a visual confirmation of the previous section's remark that many joint distributions share the same marginals.

Algebraically: it is compression. A general \(K\times K\) joint distribution has \(K^2-1\) free parameters, but assuming independence reduces that to \(2(K-1)\). At \(K=5\), 24 becomes 8, and with \(n\) variables \(K^n\) becomes \(nK\). Exponential turns into linear. In exchange, no interaction between the variables can be expressed — compression always has a price.

In machine learning. Naive Bayes assumes all features are conditionally independent given the class, writing \(p(x\mid y)=\prod_j p(x_j\mid y)\). It is plainly a false assumption, but the parameter count drops dramatically so it trains on little data, and classification performance is often surprisingly serviceable. Treating a dataset as i.i.d., making each latent dimension independent in a VAE (\(p(z)=\prod_j \mathcal{N}(z_j;0,1)\)), and adding independent noise per pixel in a diffusion model are all assumptions of the same kind.

\[\begin{aligned} X \perp Y \iff P(x,y) &= P(x)\,P(y) \quad \forall x,y \\ \iff P(y\mid x) &= P(y) \\ \text{parameter count}:\quad K^2-1 &\;\longrightarrow\; 2(K-1) \end{aligned}\]

Move the t slider from 0 to 1, or press "Force independence." The upper grid converges to the outer product of the marginals, and the residual grid below, \(P(x,y)-P(x)P(y)\), vanishes to all zeros. A nonzero residual means the joint distribution cannot be recovered from the two marginals alone.

Residual P(x,y) − P(x)P(y): independence when it reaches 0

Expectation

Expectation summarizes a distribution into a single representative value. Because the sample mean converges to the expectation as samples accumulate, expected losses and expected rewards that cannot be integrated can be replaced by minibatch averages.

Definition. The expectation of a function \(f\) is \(f(x)\) averaged with probability as the weight: \(\mathbb{E}[f(X)]=\sum_x P(x)f(x)\) in the discrete case and \(\mathbb{E}[f(X)]=\int p(x)f(x)\,dx\) in the continuous case. The case \(f(x)=x\) is what is usually called the mean \(\mu\).

Geometrically: a centre of mass. Regard the density curve as the mass distribution of a plank with no thickness, and the expectation is the position of the fulcrum where it balances. So an expectation need not be a value the distribution actually takes. The expectation 3.5 of a die roll is not on any face, and the expectation 0.7 of a binary label cannot be a label.

Linearity — the most-used property. \(\mathbb{E}[aX+bY]=a\mathbb{E}[X]+b\mathbb{E}[Y]\) holds always, even when \(X\) and \(Y\) are not independent. This property lets us split the expectation of a sum of losses term by term. Moreover, when the function is smooth enough for differentiation and integration to be exchanged, we may write \(\nabla_\theta\mathbb{E}[f_\theta(X)]=\mathbb{E}[\nabla_\theta f_\theta(X)]\). By contrast \(\mathbb{E}[XY]=\mathbb{E}[X]\mathbb{E}[Y]\) holds only under independence, and that difference is precisely the covariance two sections from here.

In machine learning. What training minimizes is not the loss on one example but the expected loss \(\mathbb{E}_{(x,y)\sim p_{\text{data}}}[\ell(f(x),y)]\). The objective of reinforcement learning is the expected reward \(\mathbb{E}_\pi[\sum_t \gamma^t r_t]\). Neither true distribution is known, so the integral cannot be done directly and is estimated by a sample mean instead — the law of large numbers is what justifies the substitution.

\[\begin{aligned} \mathbb{E}[f(X)] &= \sum_x P(x) f(x) = \int p(x) f(x)\,dx \\ \mathbb{E}[aX+bY] &= a\,\mathbb{E}[X] + b\,\mathbb{E}[Y] \quad \text{(always)} \\ \bar{x}_N = \frac{1}{N}\sum_{i} x^{(i)} &\;\xrightarrow[N\to\infty]{}\; \mathbb{E}[X] \quad \text{(law of large numbers)} \end{aligned}\]

Press "Draw samples." Samples pile up as ticks in the upper figure, while the lower figure traces the running sample mean converging to the true expectation (orange dashed line). The way it swings wildly at first and gradually settles is the same phenomenon as the noise in a minibatch gradient.

Convergence of the sample mean: it approaches E[X] as samples accumulate

Variance

Two distributions with the same mean can look completely different. Variance measures how far values are scattered from the mean, and it expresses the uncertainty of a prediction, the stability of an estimator, and the noise in optimization alike.

Definition. Variance is the expectation of the squared deviation from the mean: \(\mathrm{Var}(X)=\mathbb{E}[(X-\mu)^2]\). Expanding gives the computationally convenient form \(\mathbb{E}[X^2]-\mathbb{E}[X]^2\). Because a square was used the units are those of the original value squared, and undoing that gives the standard deviation \(\sigma=\sqrt{\mathrm{Var}(X)}\).

Why square? The expectation of the plain deviation \(X-\mu\) is always 0 by definition — the two sides cancel. To prevent that cancellation the sign must go, and the reason for squaring rather than taking an absolute value is that squares are differentiable and fall out in closed form. Choosing the absolute value gives the mean absolute deviation, and that choice leads later to MAE and the Laplace distribution.

Geometrically. The two curves in the diagram have the same mean and differ only in variance. Small variance makes a tall narrow peak; large variance a low broad one. The total area is always 1, so width and height must be traded against each other. The \(\pm\sigma\) markers show that width directly.

Why averaging reduces the wobble. The last line of the formula above, \(\mathrm{Var}(\bar{x}_N)=\mathrm{Var}(X)/N\), is the fact reused most often in this article. The reason takes two lines. By the scaling rule \(\mathrm{Var}(aX)=a^2\mathrm{Var}(X)\), the leading coefficient \(1/N\) of \(\bar{x}_N=\frac{1}{N}\sum_i x^{(i)}\) comes out squared as \(1/N^2\), and if the samples are independent the variance of the sum simply adds to \(N\cdot\mathrm{Var}(X)\). Multiplying the two gives \(\frac{1}{N^2}\cdot N\,\mathrm{Var}(X)=\mathrm{Var}(X)/N\). Intuitively, the errors of the individual samples point in different directions and partly cancel one another — if they were all wrong in the same direction (perfectly correlated), no cancellation would occur and this gain would vanish.

Which is where √N comes from. If the variance falls as \(1/N\), the standard deviation falls as its square root \(1/\sqrt{N}\). What we actually feel as error is the standard deviation, not the variance, so the conclusion of the sampling section — that raising precision tenfold requires a hundred times as many samples — follows from here. It is the same \(\sqrt{N}\) that makes quadrupling the minibatch size reduce gradient noise only by half rather than fourfold, and that gives diminishing returns to larger batches.

In machine learning: three faces. First, the predictive uncertainty a model outputs — learning \(\sigma\) alongside a Gaussian output lets the model say how confident it is. Second, the variance of an estimator — in the bias–variance decomposition, overfitting is the high-variance state. Third, the noise of optimization — a minibatch gradient is an unbiased estimator of the true gradient but carries variance, and increasing the batch size by a factor \(B\) reduces that variance by roughly \(1/B\). This is why the learning rate and batch size have to be tuned together.

\[\begin{aligned} \mathrm{Var}(X) &= \mathbb{E}[(X-\mu)^2] = \mathbb{E}[X^2]-\mathbb{E}[X]^2 \\ \sigma &= \sqrt{\mathrm{Var}(X)} \\ \mathrm{Var}(aX+b) &= a^2\,\mathrm{Var}(X) \\ \mathrm{Var}(\bar{x}_N) &= \mathrm{Var}(X)/N \end{aligned}\]

Move the two σ sliders. The mean is pinned at 0, so only the spread changes. Making one very small and the other very large makes clear how different a story "same mean, different distribution" can tell.

Covariance and Correlation

Covariance carries the direction and magnitude with which two variables move together; the correlation coefficient strips the magnitude and keeps only the direction. The axes of the ellipse drawn by the covariance matrix are its eigenvectors, and that is what PCA does.

Definition. Covariance is the average of the product of the two deviations: \(\mathrm{Cov}(X,Y)=\mathbb{E}[(X-\mu_X)(Y-\mu_Y)]=\mathbb{E}[XY]-\mathbb{E}[X]\mathbb{E}[Y]\). As the second form says, covariance measures how badly \(\mathbb{E}[XY]=\mathbb{E}[X]\mathbb{E}[Y]\) fails — exactly the difference deferred in the previous section.

Definition (correlation). Covariance depends on units. Changing height from cm to m changes the value a hundredfold. Dividing by each standard deviation removes that dependence, giving the correlation coefficient \(\rho=\mathrm{Cov}(X,Y)/(\sigma_X\sigma_Y)\), which always satisfies \(-1\le\rho\le 1\).

Geometrically: the ellipse's axes are eigenvectors. The contours of a two-dimensional Gaussian are ellipses, and the shape of that ellipse is determined by the covariance matrix \(\Sigma\). Eigendecomposing \(\Sigma\) makes the eigenvectors the principal axis directions of the ellipse and the eigenvalues the variance along those directions. Raising \(\rho\) tilts the ellipse; returning it to 0 aligns the axes with the coordinate axes. This is exactly what PCA does — finding the direction of greatest spread in the data is the same problem as finding the leading eigenvector of the covariance matrix.

Zero correlation is not independence. Independence implies zero correlation, but not conversely. If \(X\sim\mathcal{N}(0,1)\) and \(Y=X^2\), then \(Y\) is completely determined by \(X\) and yet \(\mathrm{Cov}(X,Y)=\mathbb{E}[X^3]=0\). The reason is that correlation captures only linear relationships. This is why "the features have been decorrelated" does not mean "the features are independent" in representation learning.

In machine learning. The entire shape of a multivariate Gaussian \(\mathcal{N}(\mu,\Sigma)\) is held in \(\Sigma\). PCA changes coordinates to the eigenvectors of \(\Sigma\) to reduce dimension, and whitening turns \(\Sigma\) into the identity. Techniques that regularize the covariance between features to prevent representation collapse in contrastive learning (Barlow Twins, VICReg) also work on this matrix directly.

\[\begin{aligned} \mathrm{Cov}(X,Y) &= \mathbb{E}[XY]-\mathbb{E}[X]\mathbb{E}[Y] \\ \rho &= \frac{\mathrm{Cov}(X,Y)}{\sigma_X \sigma_Y} \in [-1, 1] \\ \Sigma &= \begin{pmatrix} \sigma_X^2 & \rho\sigma_X\sigma_Y \\ \rho\sigma_X\sigma_Y & \sigma_Y^2 \end{pmatrix} = Q\Lambda Q^\top \end{aligned}\]

Move the ρ slider from −0.95 to 0.95. The point cloud and the contour ellipse tilt together, and the red and violet axes (the eigenvectors) follow that direction. Leaving ρ at 0 and changing only σx and σy keeps the axes on the coordinate axes and changes just the ellipse's aspect ratio.

The Categorical Distribution

A distribution over a finite set of choices is determined by a single probability vector, and the space that vector can occupy is the simplex of vectors summing to 1. Softmax is the map carrying whatever logit vector a network emits onto that simplex, and both classification and token prediction happen there.

Definition. The distribution of a random variable taking \(K\) values is determined by a single probability vector \(\boldsymbol{p}=(p_1,\ldots,p_K)\). The conditions are \(p_k\ge 0\) and \(\sum_k p_k=1\). The special case \(K=2\) is the Bernoulli distribution, covering coin flips and binary classification.

Geometrically: the simplex. The set of points satisfying both conditions is a \((K-1)\)-dimensional figure sitting in \(K\)-dimensional space. For \(K=3\) it is the triangle joining three vertices, called the probability simplex. The vertices are one-hot distributions (certainty) and the exact centre is the uniform distribution (complete ignorance). A probability vector a model outputs must always lie on this figure, and its degrees of freedom are \(K-1\), not \(K\).

Softmax: the map onto the simplex. The last layer of a network emits an unconstrained real vector (the logits) \(\boldsymbol{z}\). Softmax \(p_k=e^{z_k}/\sum_j e^{z_j}\) carries it onto the simplex. It does two things — \(\exp\) makes every component positive, then dividing by the total (the L1 norm) forces the sum to 1. That adding a constant to \(\boldsymbol{z}\) leaves the result unchanged (shift invariance) corresponds exactly to the fact that there are \(K-1\) degrees of freedom.

Temperature. In \(p_k \propto \exp(z_k/T)\), the parameter \(T\) controls how peaked the distribution is. As \(T\to 0\) the probability piles onto the single argmax and becomes one-hot (greedy decoding); as \(T\to\infty\) it flattens to the uniform distribution. This is the very knob that tunes generation diversity in a language model.

Connection to the loss function. When the correct answer is class \(c\), the negative log-likelihood of this distribution is \(-\log p_c\). That is the cross-entropy loss itself — because for a one-hot target \(y\), \(-\sum_k y_k \log p_k = -\log p_c\). In other words, training a classifier is solving the MLE of a categorical distribution. Take the whole vocabulary as the classes and it becomes next-token prediction in a language model.

\[\begin{aligned} P(Y=k) &= p_k, \quad p_k \ge 0, \quad \textstyle\sum_{k} p_k = 1 \\ p_k &= \mathrm{softmax}(\boldsymbol{z}/T)_k = \frac{e^{z_k/T}}{\sum_j e^{z_j/T}} \\ -\log P(Y=c) &= -\log p_c \quad \text{(cross-entropy)} \end{aligned}\]

Lower the temperature T to 0.1 — the probability piles onto one class. Raise it to 3 and the five classes become nearly uniform. Note that the logits are unchanged while only the distribution moves, and that the sum is always 1.

The Gaussian Distribution

Determined by just two parameters, a mean and a variance, it arises naturally whenever many factors add together. Using MSE in regression is precisely the same statement as assuming the noise is Gaussian.

Definition. The density of a Gaussian (normal) distribution with mean \(\mu\) and variance \(\sigma^2\) is \(p(x)=\frac{1}{\sqrt{2\pi\sigma^2}}\exp\!\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)\). Two numbers determine the whole distribution, and the leading coefficient is only the normalizing constant making the integral 1 — it plays no part in the shape.

Why this distribution in particular? There are two reasons. First, the central limit theorem — adding together many independent random variables makes the sum approach a Gaussian whatever the individual distributions were. This is why quantities accumulated from countless small factors, such as measurement noise, tend to follow a Gaussian. Second, maximum entropy — when we know only the mean and the variance and nothing else, the distribution embodying the fewest assumptions is the Gaussian. It is the choice that refrains from inventing what we do not know.

Geometrically. \(\mu\) fixes the location of the peak and \(\sigma\) its width. Changing \(\mu\) slides the shape left and right unchanged; changing \(\sigma\) widens or narrows it in place. For any \(\mu,\sigma\), about 68% falls within \(\pm 1\sigma\) and about 95% within \(\pm 2\sigma\) — the two shaded bands in the diagram.

MSE is a Gaussian assumption. Assuming \(y=f(x)+\varepsilon\) with \(\varepsilon\sim\mathcal{N}(0,\sigma^2)\) in regression gives \(p(y\mid x)=\mathcal{N}(y; f(x), \sigma^2)\). Taking the negative log brings the square inside the exponent straight down: \(-\log p(y\mid x)=\frac{(y-f(x))^2}{2\sigma^2}+\text{const}\). With \(\sigma\) fixed this is identical to MSE up to a constant factor. When we use mean squared error, we are saying the noise is Gaussian.

In machine learning. A VAE's latent prior is \(\mathcal{N}(0,I)\) and its encoder outputs \(\mathcal{N}(\mu(x),\sigma(x)^2)\). A diffusion model is itself the process of progressively adding Gaussian noise to data and undoing it, and the fact that each step's transition distribution is Gaussian is what makes training closed-form. Weight initialization, Gaussian noise as an alternative to dropout, and the approximate posterior of a Bayesian neural network all come from here too.

\[\begin{aligned} p(x) &= \frac{1}{\sqrt{2\pi\sigma^2}}\exp\!\left(-\frac{(x-\mu)^2}{2\sigma^2}\right) \\ -\log p(x) &= \frac{(x-\mu)^2}{2\sigma^2} + \tfrac{1}{2}\log(2\pi\sigma^2) \\ \text{MSE} &\;\equiv\; \text{Gaussian NLL} \;(\sigma \text{ fixed}) \end{aligned}\]

Move μ and σ. The shading marks the ±1σ (68%) and ±2σ (95%) intervals, and those proportions do not change as σ grows. Use "Draw samples" to check with tick marks what fraction of actual samples falls in those intervals.

The Laplace Distribution

It centres in the same place as a Gaussian but has heavier tails and a sharper peak. From the loss-function point of view, that single difference appears as the difference between the parabola of MSE and the V of MAE.

Definition. The density of a Laplace distribution with location \(\mu\) and scale \(b\) is \(p(x)=\frac{1}{2b}\exp\!\left(-\frac{|x-\mu|}{b}\right)\). Set beside a Gaussian, the difference is in exactly one place inside the exponent — the square has become an absolute value. Its variance is \(2b^2\).

Geometrically: a sharp peak, heavy tails. The figure overlays the two distributions with their variances matched, so only shape is being compared, not width. The Laplace has a sharp, non-differentiable corner at the centre, and in exchange both tails decay more slowly. Turning on "View tails on a log scale" makes it plain that the Gaussian falls as a parabola (fast, like a square) while the Laplace falls as a straight line (slow, like an absolute value).

What heavy tails mean. Under a Gaussian, an observation \(5\sigma\) away is effectively impossible, so the model drags its mean hard toward it in order to explain it. Under a Laplace the same observation is "rare but plausible," and the fit moves far less. Robustness to outliers comes from the thickness of the tail.

MAE is a Laplace assumption. Simply repeat the calculation of the previous section. Since \(-\log p(x)=\frac{|x-\mu|}{b}+\log 2b\), up to a constant factor this is the mean absolute error (MAE, the L1 loss). Switching a regression loss from MSE to MAE is not tinkering with a hyperparameter; it is replacing the noise model from Gaussian with Laplace.

In machine learning. The first use is MAE and the Huber loss in regression with outliers; the second is sparsity. Placing a Laplace prior on the parameters makes the penalty \(\lambda|w|\), that is L1 regularization (Lasso), and the sharp corner gives rise to solutions that are exactly zero. That story is verified directly in the final section.

\[\begin{aligned} p(x) &= \frac{1}{2b}\exp\!\left(-\frac{|x-\mu|}{b}\right), \quad \mathrm{Var}(X) = 2b^2 \\ -\log p(x) &= \frac{|x-\mu|}{b} + \log 2b \\ \text{MSE} \leftrightarrow \text{Gaussian}, &\qquad \text{MAE} \leftrightarrow \text{Laplace} \end{aligned}\]

Move b and compare the two curves — the variances are always matched, so all that changes is shape. The lower figure is the key one: plot −log p and the Gaussian becomes a parabola (MSE) while the Laplace becomes a V (MAE). The two loss functions turn out to have been two distributions.

As a loss: −log p(x). The Gaussian is a parabola (MSE), the Laplace a V (MAE)

Likelihood

Read the same formula as a function of the data and it is a probability; read it as a function of the parameters and it is a likelihood. Choosing the parameters that make the observed data most plausible is maximum likelihood estimation, and most training amounts to that.

Definition. Take the single formula \(p(x\mid\theta)\). Holding \(\theta\) fixed and reading it as a function of \(x\) gives a probability (density); holding \(x\) fixed and reading it as a function of \(\theta\) gives the likelihood \(L(\theta)=p(x\mid\theta)\). The formula is the same; only which argument is treated as the variable differs.

A likelihood is not a distribution. This is the crucial distinction. \(\int p(x\mid\theta)\,dx=1\) holds, but there is no reason whatever for \(\int L(\theta)\,d\theta\) to be 1. That is why the area under the lower curve is not 1. A likelihood is only a score expressing relative preference among values of \(\theta\); it is not a probability distribution over \(\theta\). To give \(\theta\) genuine probability we need a prior, and that is the Bayes section further on.

With several data points. If the observations are i.i.d., the joint likelihood is the product of the individual ones: \(L(\theta)=\prod_i p(x^{(i)}\mid\theta)\). Choosing the \(\theta\) that maximizes this product is maximum likelihood estimation (MLE). It is the principle "choose the explanation that makes the observed data most plausible," and training in supervised learning, language models and generative models is almost always of this form.

The Gaussian case: the answer is the sample mean. With \(p(x\mid\mu)=\mathcal{N}(x;\mu,\sigma^2)\) we have \(\log L(\mu)=-\frac{1}{2\sigma^2}\sum_i (x^{(i)}-\mu)^2+\text{const}\), so differentiating in \(\mu\) and setting to zero gives \(\hat\mu=\frac{1}{N}\sum_i x^{(i)}\). We do not use the mean because "of course it is the mean," but because the MLE solution under a Gaussian assumption is the mean. Change the noise model to a Laplace and the same calculation answers with the median.

The question MLE answers. MLE asks, "what should the parameters be for the data we just observed to be as plausible as possible?" Because it compares likelihoods only, it uses no knowledge or preference held before seeing the data. It is powerful and simple when there is enough data, but with little data or many parameters it can fit the observed sample too closely.

\[\begin{aligned} L(\theta) &= p(x \mid \theta) \quad \text{(} x \text{ fixed, } \theta \text{ as the variable)} \\ L(\theta) &= \prod_{i=1}^{N} p(x^{(i)} \mid \theta) \quad \text{(i.i.d.)} \\ \hat{\theta}_{\text{MLE}} &= \arg\max_{\theta} L(\theta) \end{aligned}\]

Click the upper figure to add data points. The likelihood curve \(L(\mu)\) below changes in real time, and its maximum (the green line) always sits exactly at the sample mean of the data. Piling points on one side pulls the curve that way.

Likelihood curve L(μ) over the parameter μ — the maximum is the MLE solution

Maximum Likelihood Estimation

If a likelihood is the score a parameter receives, maximum likelihood estimation adopts as the answer whichever parameter scores highest. Depending on which distribution family is chosen, the same principle yields estimators as different as the sample mean, the median, and the relative frequency.

Definition. Given observed data \(\mathcal{D}=\{x^{(1)},\ldots,x^{(N)}\}\) and a model \(p(x\mid\theta)\) indexed by a parameter \(\theta\), the maximum likelihood estimator (MLE) is the parameter maximizing the likelihood: \(\hat\theta_{\text{MLE}}=\arg\max_{\theta\in\Theta}L(\theta)\), where \(L(\theta)=p(\mathcal{D}\mid\theta)\). If the observations are i.i.d. then \(L(\theta)=\prod_i p(x^{(i)}\mid\theta)\), and since \(\log\) is strictly increasing, the same \(\hat\theta\) may be found as the maximizer of the log-likelihood \(\ell(\theta)=\sum_i \log p(x^{(i)}\mid\theta)\).

What is held fixed and what moves. The distinction of the previous section now does real work. The data \(\mathcal{D}\) has already been observed and is a fixed constant; the only thing that moves is \(\theta\). So MLE does not ask "how might data come out" but "what explanation makes the data that already came out least surprising?" It shows the same data to many candidate values of \(\theta\), asks each "what probability does your worldview assign to this data?", and adopts the candidate that names the largest number.

The procedure. In practice one almost always walks the same four steps. ① Choose a distribution family \(p(x\mid\theta)\) (this is the assumption, and it determines every answer that follows). ② Write down the log-likelihood \(\ell(\theta)=\sum_i\log p(x^{(i)}\mid\theta)\). ③ Solve \(\nabla_\theta\ell=0\) (this equation is called the likelihood equation). ④ Confirm the solution is a maximum via a second-order condition or the concavity of the function. When there is no closed form, step ③ is replaced by gradient descent — and that is neural network training.

Three standard examples. The same procedure produces statistics we already know. For a Bernoulli, \(p(x\mid\theta)=\theta^x(1-\theta)^{1-x}\), solving \(\ell'(\theta)=0\) gives \(\hat\theta=\frac{1}{N}\sum_i x^{(i)}\), the observed relative frequency — flip a coin ten times and see seven heads, and the MLE is 0.7. For a Gaussian, \(\hat\mu\) is the sample mean and \(\hat\sigma^2=\frac{1}{N}\sum_i(x^{(i)}-\hat\mu)^2\) the sample variance. For a Laplace, one must minimize \(\sum_i|x^{(i)}-\mu|\), so the answer is the median. The values we took to be "obvious summaries" were in fact MLE solutions under particular distributional assumptions.

Intuition: choosing the family is choosing the estimator. The difference among those three results comes not from the data but from the assumption. A Gaussian penalizes squared error, so a single distant point protests in proportion to \((x-\mu)^2\) and the mean is dragged toward it. A Laplace penalizes only the absolute value, so however far away the point is the strength of its protest is constant, and the median does not budge. Drag an outlier in the demo and you can watch the two estimates separate. The message of the earlier sections — that choosing a loss function is choosing a noise model — becomes, under the single principle of MLE, another way of saying the same sentence.

Good properties. There are reasons MLE is used so widely. Consistency — if the model contains the true distribution, \(\hat\theta_N\) converges to the true value as the sample grows. Asymptotic normality — for large \(N\) the distribution of \(\hat\theta\) approaches a Gaussian around the true value whose width is set by the Fisher information, which yields confidence intervals. Invariance — to estimate \(g(\theta)\) there is no need to solve again; \(g(\hat\theta_{\text{MLE}})\) will do. That the MLE of the standard deviation is the square root of the MLE of the variance is an instance.

Three limitations. First, it can be biased. The MLE of a Gaussian variance divides by \(N\), and because it measures deviations against the \(\hat\mu\) already fitted from the data, it comes out systematically smaller than the true value (\(\mathbb{E}[\hat\sigma^2]=\frac{N-1}{N}\sigma^2\)). The familiar sample variance dividing by \(N-1\) is a different estimator that corrects this bias. Second, overfitting — with little data or many parameters the likelihood conforms too closely to the observed sample. In the extreme, if an event is never observed, MLE assigns it probability exactly 0, and as the information theory article shows, the loss becomes infinite the moment that event actually occurs. Third, it is a point estimate — it returns a single value and says nothing about how confident it is.

Which leads to the following sections. Maximizing the product directly breaks down numerically, so we take a log to turn it into a sum and flip the sign — that is the NLL of the next section, the quantity actually minimized in implementations. And addressing the three limitations above requires putting a prior on \(\theta\), which is Bayes' rule and MAP. In fact, setting the prior to be uniform in MAP returns exactly MLE — MLE is the special case of MAP that declares "all parameter values are equally preferred."

In machine learning. Nearly every standard supervised loss is a conditional MLE. For a conditional model \(p(y\mid x,\theta)\) one solves \(\hat\theta=\arg\max_\theta\sum_i\log p(y^{(i)}\mid x^{(i)},\theta)\), and inserting a categorical distribution gives cross-entropy, a Gaussian gives MSE, and a Laplace gives MAE. Next-token prediction in a language model is the MLE of a categorical distribution over the vocabulary, extended to a whole sentence by the chain rule, and normalizing flows and autoregressive generative models likewise maximize the log-likelihood directly. "Training a model" is usually shorthand for this sentence.

\[\begin{aligned} \hat{\theta}_{\text{MLE}} &= \arg\max_{\theta} \prod_{i=1}^{N} p(x^{(i)}\mid\theta) = \arg\max_{\theta} \sum_{i=1}^{N} \log p(x^{(i)}\mid\theta) \\ \text{Bernoulli}(\theta) &\;\longrightarrow\; \hat\theta = \tfrac{1}{N}\textstyle\sum_i x^{(i)} \quad (\text{relative frequency}) \\ \mathcal{N}(\mu,\sigma^2) &\;\longrightarrow\; \hat\mu = \bar{x} \quad (\text{sample mean}) \\ \mathrm{Laplace}(\mu,b) &\;\longrightarrow\; \hat\mu = \mathrm{median}(x) \quad (\text{median}) \end{aligned}\]

Switch between Gaussian and Laplace with the family buttons. The data is the same, yet the MLE solution moves from the sample mean to the median. Click near the right edge of the upper figure to add an outlier, and watch the Gaussian solution get dragged while the Laplace solution barely moves — a distributional assumption is an attitude toward outliers.

Log-likelihood ℓ(μ) — the maximum (green) is the MLE solution

View Python code
import numpy as np
from scipy.optimize import minimize_scalar

x = np.array([-1.2, 0.3, 0.9, 1.8, 2.4, 9.0])   # the last point is an outlier

# Gaussian MLE: minimize squared error → sample mean
gauss_nll = lambda mu: np.sum((x - mu) ** 2)
# Laplace MLE: minimize absolute error → median
laplace_nll = lambda mu: np.sum(np.abs(x - mu))

print(minimize_scalar(gauss_nll).x,   np.mean(x))      # 2.20 2.20
print(minimize_scalar(laplace_nll).x, np.median(x))    # 1.35 1.35

# the Bernoulli MLE is just the relative frequency
flips = np.array([1, 0, 1, 1, 1, 0, 1, 1, 0, 1])
print(flips.mean())                                     # 0.7

Log-Likelihood and NLL

A likelihood is a product of probabilities, so it collapses numerically to zero as soon as the data grows even slightly. This single transformation — taking a log to turn the product into a sum — is why machine learning loss functions are almost without exception of negative log-likelihood form.

Definition. Taking the log of the likelihood gives the log-likelihood \(\ell(\theta)=\log L(\theta)=\sum_i \log p(x^{(i)}\mid\theta)\), and flipping the sign to make it a minimization problem gives the negative log-likelihood (NLL) \(-\ell(\theta)\). Since by optimization convention a loss is minimized, it is the latter we use.

Reason 1: products collapse numerically. If each \(p(x^{(i)})\) is about 0.2, then at \(N=100\) the product is \(10^{-70}\), and at \(N=500\) it passes the smallest magnitude float64 can represent (\(\approx 10^{-308}\)) and becomes exactly 0. The moment it hits 0 we can no longer compare which \(\theta\) is better, and the gradient dies too. Taking a log turns the product into a sum and the problem disappears wholesale.

Reason 2: the argmax is preserved. Because \(\log\) is strictly increasing, \(L(\theta_1)>L(\theta_2)\) and \(\log L(\theta_1)>\log L(\theta_2)\) are equivalent. So taking the log does not change the optimum in the slightest. It is a nearly free transformation that swaps in an easier objective while leaving the answer alone.

Reason 3: the derivative decomposes. The derivative of a product tangles terms together through the product rule, but the derivative of a sum is the sum of the term derivatives. So \(\nabla_\theta \ell = \sum_i \nabla_\theta \log p(x^{(i)}\mid\theta)\), and approximating that sum with a minibatch is SGD. The earlier logic — linearity of expectation leading to minibatch estimation — is reused verbatim here.

And so every loss function is an NLL. Gathering what we have seen produces the table below — choosing a loss function is really choosing a distribution.

\[\begin{aligned} \text{Categorical} &\;\longrightarrow\; -\log p_c \;=\; \text{cross-entropy} \\ \text{Gaussian} &\;\longrightarrow\; \tfrac{1}{2\sigma^2}(y-\hat{y})^2 \;=\; \text{MSE} \\ \text{Laplace} &\;\longrightarrow\; \tfrac{1}{b}|y-\hat{y}| \;=\; \text{MAE} \\ \mathcal{L}(\theta) &= -\sum_{i=1}^{N} \log p(x^{(i)} \mid \theta) \end{aligned}\]

Raise the N slider. The "product" in the readout shrinks rapidly and at some point collapses to 0 (underflow), while the "log sum" descends linearly as though nothing had happened. The red vertical line marks where float64 gives up on the product.

Bayes' Rule

Multiply a prior by a likelihood and normalize, and you have a posterior. It is the rule describing how belief held before observing data is updated by observation, and it is nothing more than the definition of conditional probability turned around.

Definition. Decomposing \(P(x,\theta)\) in two different orders (the chain rule) gives \(P(\theta\mid x)P(x)=P(x\mid\theta)P(\theta)\), and dividing both sides by \(P(x)\) yields Bayes' rule: \(P(\theta\mid x)=\frac{P(x\mid\theta)P(\theta)}{P(x)}\). It is not a new axiom but merely the definition of conditional probability turned around.

The names of the four pieces. \(P(\theta)\) is the belief before seeing data, the prior; \(P(x\mid\theta)\) is the likelihood of the previous section; \(P(\theta\mid x)\) is the belief after seeing data, the posterior; and the denominator \(P(x)=\int P(x\mid\theta)P(\theta)\,d\theta\) is the evidence. The evidence does not depend on \(\theta\) and is a normalizing constant, so it is usually written only as \(P(\theta\mid x)\propto P(x\mid\theta)P(\theta)\) — but that integral is precisely the marginalization bottleneck seen earlier.

How it differs from MLE. MLE looks only at the likelihood and picks a single \(\theta\). Bayes multiplies by a prior and obtains not one point but a whole distribution. So it can answer not only "the most plausible \(\theta\)" but also "how certain that estimate is." The less data there is, the broader the posterior and the closer to the prior; as data accumulates it narrows and the likelihood takes over.

The example: estimating a coin's bias. Let \(\theta\) be the probability of heads and take the prior to be \(\mathrm{Beta}(2,2)\) ("probably somewhere near half"). Observing \(h\) heads and \(t\) tails makes the posterior exactly \(\mathrm{Beta}(2+h,\,2+t)\) — a pairing in which prior and posterior stay in the same family is called a conjugate prior, and thanks to it the update finishes with parameter addition and no integral at all.

In machine learning. Bayesian neural networks put a posterior on the weights, and uncertainty estimation and active learning are built on top of that. The more common practical use is the MAP of the next section — take only the maximum instead of computing the whole posterior, and the problem becomes a familiar optimization with a regularization term attached.

\[\begin{aligned} P(\theta \mid x) &= \frac{P(x \mid \theta)\,P(\theta)}{P(x)}, \qquad P(x) = \int P(x\mid\theta)P(\theta)\,d\theta \\ \text{posterior} &\;\propto\; \text{likelihood} \times \text{prior} \\ \mathrm{Beta}(a,b) \;\xrightarrow{\;h,\,t\;}\; &\mathrm{Beta}(a+h,\, b+t) \end{aligned}\]

Press "Accumulate observations." The violet posterior starts from the grey dashed prior and narrows around the true value (green line) as observations pile up. Compare how far the posterior mean and the MLE differ when there are five or fewer observations, and how close they become by about fifty.

Priors and MAP

MAP estimation, which takes the maximum of the posterior, is an optimization with one extra term — the prior — added to the likelihood. A Gaussian prior becomes L2 regularization and a Laplace prior becomes L1, and only the latter drives solutions exactly to zero.

Definition. Taking only the single maximum of the posterior rather than handling the whole distribution is maximum a posteriori (MAP) estimation: \(\hat\theta_{\text{MAP}}=\arg\max_\theta P(x\mid\theta)P(\theta)\). The evidence \(P(x)\) does not involve \(\theta\), so it has no effect on the argmax and can be dropped.

The difference between MLE and MAP. Both are point estimates picking a single parameter, but the material for the judgement differs. MLE maximizes \(P(x\mid\theta)\) alone and picks the value fitting the data best; MAP multiplies in the prior \(P(\theta)\) and picks a value that both fits the data and is plausible under prior knowledge. With little data the prior does much to stabilize MAP, and when there is a lot of data and the likelihood becomes very peaked, MLE and MAP usually come close. Note, though, that MAP keeps only a single maximum rather than the entire posterior, so it does not preserve the full uncertainty the way complete Bayesian inference does.

Take a log and regularization appears. Following the habit of the earlier section, taking the negative log splits the product into a sum: \(-\log P(x\mid\theta)-\log P(\theta)\). The first term is the NLL we already know (the data term), and the second is the regularization term. Regularization is not a device bolted on afterwards for the sake of optimization; it is a term that follows automatically the moment we hold a prior.

Gaussian prior = L2 = ridge. Setting \(P(\theta)=\mathcal{N}(0,\tau^2)\) gives \(-\log P(\theta)=\frac{\theta^2}{2\tau^2}+\text{const}\), that is \(\frac{\lambda}{2}\|\theta\|^2\). In the one-dimensional squared-error example on the right the solution shrinks proportionally as \(\hat\theta=\hat\theta_{\text{MLE}}/(1+\lambda)\). This closed form is specific to that example's scaling and is not a formula that applies unchanged to every model. In general, though, L2 tends to shrink parameters smoothly toward zero.

Laplace prior = L1 = lasso. Setting \(P(\theta)\propto e^{-|\theta|/b}\) makes the penalty \(\lambda|\theta|\). In the same one-dimensional squared-error example the solution takes the soft-thresholding form \(\hat\theta=\mathrm{sign}(\hat\theta_{\text{MLE}})\max(0,\,|\hat\theta_{\text{MLE}}|-\lambda)\), reaching exactly zero once \(\lambda\) is large enough. The exact scale of the coefficients depends on how the objective is defined, but the property that L1 induces sparse solutions because of the sharp corner at the origin holds regardless.

MLE is a special case of MAP. Make the prior flat (uniform) and \(-\log P(\theta)\) becomes a constant and disappears, and MAP falls back to MLE. Turning regularization off (\(\lambda=0\)) is the same statement as declaring "I prefer all parameter values equally."

\[\begin{aligned} \hat\theta_{\text{MAP}} &= \arg\min_\theta \big[\underbrace{-\log P(x\mid\theta)}_{\text{NLL}} \;\underbrace{-\log P(\theta)}_{\text{regularization}}\big] \\ \mathcal{N}(0,\tau^2) &\;\longrightarrow\; \tfrac{\lambda}{2}\|\theta\|_2^2 \;\;(\text{L2, Ridge}) \\ \mathrm{Laplace}(0,b) &\;\longrightarrow\; \lambda\|\theta\|_1 \;\;(\text{L1, Lasso}) \end{aligned}\]

Raise λ from 0 and alternate between the two priors. Under the Gaussian the green solution only ever approaches zero, but under the Laplace it reaches exactly zero the moment λ passes 1.6. The kink in the violet objective curve at the origin is the cause.

← Back to all posts