Distributions and MLE

A model in machine learning is often a probability distribution with knobs. Training means turning the knobs until the data you observed looks as probable as possible. That idea, maximum likelihood estimation (MLE), is behind linear regression, logistic regression, and the cross-entropy loss used to train language models.

The normal distribution

$\mu$ moves the bell left and right, and $\sigma$ makes it wider and flatter. The area under the curve is always 1, so a wider curve must also be lower. Try it:

mean μ

0

std dev σ

1

Likelihood

Suppose we observe $n$ independent data points $x_1, \dots, x_n$ and believe they came from $\mathcal{N}(\mu, \sigma^2)$ for some unknown $\mu$ and $\sigma$. The likelihood of a choice of parameters is the probability density of the whole sample under those parameters:

$$L(\mu, \sigma) = \prod_{i=1}^{n} f(x_i;\, \mu, \sigma).$$

Products of many small numbers underflow, and sums are easier to differentiate, so we work with the log-likelihood instead:

$$\ell(\mu, \sigma) = \log L = -\frac{n}{2}\log(2\pi\sigma^2) - \frac{1}{2\sigma^2}\sum_{i=1}^{n}(x_i - \mu)^2.$$

Since $\log$ is increasing, the parameters that maximize $\ell$ also maximize $L$.

Solving it analytically

Set the derivatives to zero. For $\mu$:

$$\frac{\partial \ell}{\partial \mu} = \frac{1}{\sigma^2}\sum_{i=1}^{n}(x_i - \mu) = 0 \quad\Longrightarrow\quad \hat{\mu} = \frac{1}{n}\sum_{i=1}^{n} x_i.$$

Doing the same for $\sigma$ gives the second half of the result:

Here is the log-likelihood of that sample as a function of $\mu$, with $\sigma$ held at $\hat{\sigma}$. It is a downward parabola with its peak exactly at the sample mean:

Solving it numerically

Most models don’t have a closed-form MLE, so in practice we minimize the negative log-likelihood with an optimizer. Doing it here, where we know the answer, is a good check that the machinery works:

import numpy as np
from scipy.optimize import minimize

x = np.array([1.2, 0.4, 2.1, 1.7, 0.9])

def neg_log_lik(params):
    mu, log_sigma = params          # optimize log(sigma) so sigma stays positive
    sigma = np.exp(log_sigma)
    return 0.5 * len(x) * np.log(2 * np.pi * sigma**2) + np.sum((x - mu) ** 2) / (2 * sigma**2)

result = minimize(neg_log_lik, x0=[0.0, 0.0])
mu_hat, sigma_hat = result.x[0], np.exp(result.x[1])

print(mu_hat, sigma_hat)            # ≈ 1.26, 0.595
print(x.mean(), x.std(ddof=0))      # 1.26, 0.5953  (the analytic answer)

Optimizing $\log\sigma$ instead of $\sigma$ is a common trick: it turns a constrained problem ($\sigma > 0$) into an unconstrained one.

Exercise

Draw 1,000 points from $\mathcal{N}(3, 2^2)$ with a fixed seed. Fit $\mu$ and $\sigma$ analytically and with scipy.optimize.minimize, and show that the two agree to at least four decimal places. Then repeat with only 10 points and compare $\hat{\sigma}$ using $n$ and $n - 1$. Which one is closer to 2 on average over 1,000 repetitions?

Key ideas

• A distribution with parameters is a model; the likelihood scores how well parameters explain the data.
• Maximize the log-likelihood, not the likelihood: same answer, better numerics.
• For a Gaussian, the MLE is the sample mean and the (divide-by-$n$) sample variance.
• When there’s no closed form, minimize the negative log-likelihood with an optimizer.