Multiclass Classification and the Softmax
Logistic regression answers one question: is this a yes or a no? Most real classification problems have more than two answers — a handwritten digit is one of ten, an email is one of several folders, a photograph holds one of a thousand objects. This note generalises the whole of logistic regression to $N$ classes at once. The result is the softmax, and almost nothing else has to change: the derivation is the same argument about log-odds, the loss is the same maximum-likelihood construction, and the gradient comes out in the same prediction minus truth form.
In one line. Give every class its own score $z_j = \vec{w}_j\cdot\vec{x} + b_j$, turn the scores into probabilities with $a_j = e^{z_j} / \sum_k e^{z_k}$, and train on the negative log-probability of the correct class; with $N = 2$ this is logistic regression.
- From two classes to many
- The softmax function
- Where the softmax comes from
- The loss: multiclass cross-entropy
- The gradients
- Convexity
- Softmax as a neural network output layer
- The numerical implementation
- Multi-label is a different problem
- Adam: learning rates that adapt
- Binary and multiclass, side by side
- Where this goes next
From two classes to many
What logistic regression already gives us
Logistic regression produces a single number and derives the second one from it. With $z = \vec{w}\cdot\vec{x} + b$:
Two probabilities that sum to $1$, from one score. The full story — the sigmoid, the log-odds reading, the log loss and why squared error fails — is in logistic-regression.md; this note picks up exactly where that one leaves off.
What multiclass has to produce
Now let the target take $N$ values, $y \in {1, 2, \ldots, N}$ — ten digits, four folders, a thousand objects. The model can no longer emit one number, because $N-1$ of the probabilities are no longer determined by the first. It has to produce the whole distribution at once:
A vector that is positive everywhere and sums to one is a much stronger requirement than "a number between $0$ and $1$", and it is what rules out the obvious ideas. Running $N$ separate sigmoids gives $N$ numbers in $(0,1)$ that sum to whatever they like. The softmax is the construction that satisfies the constraint by design.
The softmax function
One score per class
Give every class its own weight vector and its own bias, and compute one score per class from the same input:
These scores are usually called the logits. Nothing constrains them: a logit can be $-4$ or $+30$. All the work of turning them into probabilities happens in the next step.
Exponentiate and normalise
Both requirements from Equation 2 are met, and it is worth seeing how rather than taking it on trust. The exponential $e^{z}$ is positive for every real $z$, however negative — so every numerator is positive, and so is the denominator, being a sum of them. That gives $a_j > 0$. And the denominator is the sum of all the numerators, so adding up the $a_j$ adds up to the denominator over itself:
That shared denominator is the important structural fact about the softmax, and everything distinctive about it follows from there — including the awkward parts in §7.2 and §8.
A worked example
Take three classes with scores $\vec{z} = \begin{bmatrix} 3.0 & 1.0 & 0.2 \end{bmatrix}$:
| Class $j$ | score $z_j$ | $e^{z_j}$ | $a_j = e^{z_j} / 24.025$ |
|---|---|---|---|
| $1$ | $3.0$ | $20.086$ | $0.836$ |
| $2$ | $1.0$ | $2.718$ | $0.113$ |
| $3$ | $0.2$ | $1.221$ | $0.051$ |
| total | $1.000$ | ||
Read as a prediction: class 1 with about $84\%$ confidence, class 2 at $11\%$, class 3 at $5\%$. Notice how unforgiving the exponential is — class 1 leads class 2 by only $2$ points of score but takes seven times the probability. Softmax turns small score differences into large probability differences, which is what makes it decisive, and also what makes a badly calibrated model badly overconfident.
The classes compete
Because every output shares the one denominator, the probabilities are not computed independently: raising one class's score lowers everyone else's probability, even though their own scores never moved.
It contains the sigmoid
Set $N = 2$ and the softmax should reduce to logistic regression, or the generalisation is not one. It does, exactly:
What the collapse teaches. Only the difference $z_1 - z_2$ matters, not the two scores individually. That is why binary logistic regression needs one score and not two: the second is redundant, and fixing $z_2 = 0$ loses nothing. The same redundancy exists for every $N$ — the softmax has one more parameter set than it strictly needs, which is the subject of §3.4.
Where the softmax comes from
The formula in Equation 4 can look arbitrary — why the exponential, and not squaring, or absolute values, which are also positive? It is not arbitrary. Repeat the argument that produced logistic regression, with $N$ classes instead of $2$, and the softmax is what comes out.
Why the derivation starts at the probability
Both derivations begin at $p$, the quantity to be predicted, and work outwards to a linear model — never the other way round. The reason is that $p$ is trapped in $(0,1)$ while $\vec{w}\cdot\vec{x} + b$ ranges over all of $\mathbb{R}$, so the two cannot simply be set equal. Something has to carry $(0,1)$ onto $\mathbb{R}$ first, and the requirements on that link function are strict enough to pin it down almost uniquely:
| Requirement | Why it is needed |
|---|---|
| Maps $(0,1) \to \mathbb{R}$ | The linear score is unbounded; the probability is not. The link has to absorb the mismatch. |
| Monotonic | A higher score must mean a higher probability, or the model is uninterpretable. |
| Invertible | Training works in score space, prediction is needed in probability space. You must be able to go back. |
| Differentiable | Gradient descent needs a derivative at every point. |
| Statistically grounded | Starting from the Bernoulli (or categorical) likelihood is what makes cross-entropy the loss rather than an arbitrary choice. |
The function that satisfies them is the log-odds, or logit, $\log\left(p / (1-p)\right)$, and modelling it linearly is the assumption from which everything else follows. Starting instead from an unconstrained linear model and hoping the output behaves like a probability gives up all five properties at once.
Relative log-odds against a reference class
With $N$ classes there is no single "odds" to take, because there is no single alternative — so pick one class as the reference and model the log-odds of every other class against it. Taking class $N$ as the reference:
Exponentiating gives every probability in terms of the reference one:
Solving for the probabilities
One condition is still unused: the probabilities have to sum to $1$. Impose it and the reference probability is forced.
Summing Equation 8 over all classes, and writing $z_N = 0$ for the reference (its log-odds against itself is $\log 1 = 0$): $$1 = \sum_{j=1}^{N} P(y = j \mid \vec{x}) = P(y = N \mid \vec{x})\sum_{j=1}^{N} e^{z_j},$$ so $P(y = N \mid \vec{x}) = 1 / \sum_j e^{z_j}$. Substituting that back into Equation 8 gives every other class.
That is Equation 4. The exponential is not a design choice — it is the inverse of the logarithm in the log-odds, and the denominator is not a normalising hack but the unique constant that makes the distribution valid.
The reference class does not matter
Choosing class $N$ as the reference looks like it should bias the model, and it does not, because of a property worth stating on its own: adding the same constant to every score leaves the softmax unchanged.
Two consequences follow, one theoretical and one intensely practical.
- The model has one redundant degree of freedom. Only score differences carry information, so a different reference class is the same model with shifted scores. This is why the two-class case in §2.5 collapses to a single sigmoid.
- You may subtract anything you like before exponentiating. In particular you may subtract $\max_k z_k$, which is the standard defence against overflow — see §8.2.
The loss: multiclass cross-entropy
From maximum likelihood
The loss is not invented either. It is the same five-step maximum-likelihood construction used for the binary log loss, with the Bernoulli distribution replaced by the categorical one.
| Step | Mathematics | In words |
|---|---|---|
| Model | $p_k^{(i)} = \text{softmax}\left(\vec{z}^{(i)}\right)_k$ | The predicted probability of class $k$ for example $i$. |
| Likelihood of one example | $L^{(i)} = \prod_{k=1}^{N} \left(p_k^{(i)}\right)^{y_k^{(i)}}$ | Every exponent is $0$ except the true class's, so only that one factor contributes. |
| Likelihood of the dataset | $L = \prod_{i=1}^{m} L^{(i)}$ | Examples are assumed independent, so probabilities multiply. |
| Log-likelihood | $\log L = \sum_{i=1}^{m}\sum_{k=1}^{N} y_k^{(i)} \log p_k^{(i)}$ | The logarithm turns the product into a sum. |
| Loss to minimise | $J = -\frac{1}{m}\log L$ | Negate to turn "maximise" into "minimise", average to make it independent of $m$. |
One-hot targets, and what survives the sum
The target $y^{(i)}$ is a class label, but the loss wants a vector. One-hot encoding supplies it: a length-$N$ vector that is $1$ in the position of the true class and $0$ everywhere else. For $N = 4$ and a true class of $3$, $\;y = \begin{bmatrix} 0 & 0 & 1 & 0\end{bmatrix}$.
That makes the inner sum in Equation 11 far less intimidating than it looks. Every term is multiplied by a $y_k^{(i)}$ that is zero except one, so the entire sum over $k$ collapses to a single term:
What the loss does to one example
$-\log p$ is worth looking at directly, because its shape is the whole training signal.
The asymmetry is the point. A loss that just counted mistakes would give the same score to a model that was wrong by a hair and one that was catastrophically, confidently wrong. Cross-entropy punishes confident errors without limit, which is precisely the behaviour you want from a probability model.
The gradients
The derivative with respect to the scores
Everything rests on one result, and it is remarkably clean:
Where it comes from. Write $L = -\log p_c$ with $c$ the true class and $p_c = e^{z_c} / S$, where $S = \sum_j e^{z_j}$. Then $L = -z_c + \log S$, which is much easier to differentiate than it first appears. Since $\partial S / \partial z_k = e^{z_k}$, the second term contributes $e^{z_k}/S = p_k$ for every $k$. The first term contributes $-1$, but only when $k = c$ — and "$1$ when $k$ is the true class, $0$ otherwise" is exactly $y_k$. Adding the two gives $p_k - y_k$, with no case analysis left over.
Compare with logistic regression, where the same derivative is $a - y$. It is the identical statement: how wrong were you, in probability. The exponentials and the shared denominator all cancel, which is not luck — it is the pairing of the softmax with its matching log-likelihood loss, and it is the reason those two are never used apart.
The derivative with respect to the parameters
From the scores to the parameters is one application of the chain rule. Since $z_\ell = \vec{w}_\ell\cdot\vec{x} + b_\ell$, the score $z_\ell$ depends on $\vec{w}_k$ only when $\ell = k$, and then the derivative is just $\vec{x}$:
Read what training actually does with these. For the true class, $y_k = 1$ so $p_k - y_k$ is negative, and the update pushes $\vec{w}_k$ toward $\vec{x}$ — raising that class's score on inputs like this one. For every other class, $y_k = 0$ so $p_k - y_k$ is positive, and the update pushes their weights away. Each example simultaneously promotes one class and demotes the rest, in proportion to how much probability each was wrongly holding.
Convexity
For a linear softmax model — no hidden layers — the cross-entropy loss is globally convex in the parameters, so gradient descent cannot get stuck in a bad local minimum. Three facts compose:
- The Hessian in score space is positive semi-definite. For one example with probability vector $\vec{p}$, the second-derivative matrix of the loss with respect to $\vec{z}$ is $H_z = \text{diag}(\vec{p}) - \vec{p}\,\vec{p}^{\mathsf{T}}$ — which is the covariance matrix of a categorical random variable, and every covariance matrix is positive semi-definite.
- An affine map preserves convexity. The scores are affine in the parameters, $\vec{z} = W\vec{x} + \vec{b}$, and a convex function composed with an affine map is still convex.
- Averaging preserves convexity. The cost is a mean over examples, and a mean of convex functions is convex.
The binary case is the same argument with a scalar in place of a matrix: there the second derivative is $\sigma(z)\left(1 - \sigma(z)\right) > 0$, positive at every $z$, and steps 2 and 3 are unchanged.
The caveat that matters later. This guarantee covers the linear model only. Put a softmax on top of a neural network with hidden layers and the composition is no longer affine — the cost surface becomes non-convex, with all the local minima that implies. The loss is still the right loss and the gradient is still $\vec{p} - \vec{y}$; only the promise about global optimality is lost.
Softmax as a neural network output layer
The architecture
In practice the softmax is rarely used on raw inputs. It is the output layer of a network: hidden layers with ReLU do the representation work, and the final layer has one unit per class with softmax as its activation. For handwritten digits, a standard small network is $25$ ReLU units, then $15$, then $10$ softmax units — one per digit.
The last layer works exactly like Equation 3, with the previous layer's activations playing the part of $\vec{x}$: each output unit computes its own score from the same incoming vector, and then the softmax is applied across all of them.
The one activation that is not element-wise
Every other activation in a network is a scalar function applied separately to each unit: $a = \text{ReLU}(z)$ or $a = \sigma(z)$, one unit at a time, with no reference to its neighbours. The softmax is the exception. Because of the shared denominator, computing $a_1$ requires $z_1$ and every other score in the layer:
This is worth knowing for two reasons. It is why the softmax layer's backward pass is a genuine Jacobian rather than an element-wise multiply — the $H_z$ of §6 is that Jacobian in disguise. And it is why frameworks treat the softmax as belonging to the loss rather than to the layer, which is the subject of the next section.
The numerical implementation
Roundoff error is real
Before the softmax specifically, a fact about computer arithmetic that this layer makes painful. Computers store numbers in binary, so a decimal like $0.0002$ has no exact representation — it is stored as the nearest available binary fraction, in the same way $1/3$ has no exact decimal. Two calculations that are equal on paper can therefore differ in the answer:
x1 = 2.0 / 10000
print(f"{x1:.18f}") # 0.000200000000000000
x2 = (1 + 1/10000) - (1 - 1/10000)
print(f"{x2:.18f}") # 0.000199999999999978
Both should be $0.0002$. The second is wrong in the fourteenth decimal place, and it is wrong for a specific reason: it subtracted two nearly equal large-ish numbers, and the digits they had in common cancelled, leaving only the noisy tail. The error is tiny in absolute terms ($2 \times 10^{-17}$) and it does not matter here — but a loss function that repeats such a step millions of times, and then differentiates through it, is a different matter.
Two ways softmax breaks
Implemented literally — compute $\vec{a}$ with Equation 4, then take $\log a_c$ for the loss — the softmax has two failure modes, and both are routine rather than exotic:
- Overflow. A score of $z = 1000$ is unremarkable for an untrained network, and $e^{1000}$ is not representable in double precision at all. The computation returns infinity, and $\infty / \infty$ is not a number. The cure is Equation 10: subtract $\max_k z_k$ from every score first. Nothing about the result changes, and the largest exponent is now $e^0 = 1$.
- Underflow into the logarithm. A correct class given probability $10^{-320}$ rounds to exactly $0.0$, and $\log 0$ is $-\infty$. The loss becomes infinite and the gradients become meaningless — from a prediction that was merely very bad, not impossible.
The fix: keep the logits
Both problems come from computing the probability and then its logarithm, as two separate steps. The fix is to never form the probability at all. Expand the loss algebraically:
This is why frameworks ask you to hand them the raw scores. In practice the output layer is declared with a linear activation and the softmax is folded into the loss:
model = Sequential([
Dense(25, activation='relu'),
Dense(15, activation='relu'),
Dense(10, activation='linear'), # the scores, NOT the probabilities
])
model.compile(loss=SparseCategoricalCrossentropy(from_logits=True))
from_logits=True is the instruction "these are scores; do the softmax and the logarithm
together, the stable way." The cost is that the model's output is no longer a probability, so
predictions have to be pushed through a softmax explicitly afterwards. That is the trade the
whole ecosystem makes: a little inconvenience at prediction time in exchange for a loss that
cannot silently produce infinities during training.
Multi-label is a different problem
There is a problem that looks like multiclass and is not, and mixing them up is a common mistake. Consider a photograph from a car's camera, with three questions: is there a car? a bus? a pedestrian? The answers are not exclusive — a street scene can easily contain all three, or none.
| Multiclass | Multi-label | |
|---|---|---|
| Question | Which one of $N$? | Which of $N$, independently? |
| Example | Which digit is this, $0$ to $9$? | Does the image contain a car, a bus, a pedestrian? |
| Target | one-hot, e.g. $\begin{bmatrix}0&1&0\end{bmatrix}$ | any pattern, e.g. $\begin{bmatrix}1&0&1\end{bmatrix}$ |
| Outputs sum to | $1$, always | anything from $0$ to $N$ |
| Output activation | one softmax over all units | a separate sigmoid on each unit |
| Loss | categorical cross-entropy | binary cross-entropy, summed over units |
You could train three completely separate networks, one per question. Usually you do not: one network with three sigmoid output units shares all the hidden layers, so the features that detect "vehicle-shaped thing in the road" are learned once and reused by every output. The architecture is identical to Figure 3 — only the output activation and the loss change.
The one-line test. Ask whether the answers can be simultaneously true. If they cannot, the classes compete for a fixed budget of probability and you want a softmax. If they can, each unit gets its own independent sigmoid, and the shared denominator that makes softmax work would actively get in the way.
Adam: learning rates that adapt
One last piece of practical machinery, because it is what these networks are actually trained with. Plain gradient descent uses one learning rate for every parameter, which the contour picture in linear-regression.md showed to be the wrong tool for a stretched cost surface: too large in the steep directions, too small in the shallow ones.
Adam — Adaptive Moment estimation — keeps a separate learning rate for every single parameter and adjusts each one during training, using a simple observation about what the recent gradients have been doing:
| What the recent steps look like | What it means | What Adam does |
|---|---|---|
| Consistently in the same direction | The parameter is far from where it needs to be and the path is not in doubt | Increase its learning rate — go faster |
| Oscillating back and forth | The steps are jumping across a narrow valley | Decrease its learning rate — go slower |
In code it is a one-line change, and it is close to a default in modern practice — the initial learning rate still has to be set, but the algorithm is dramatically less sensitive to getting it exactly right:
model.compile(
optimizer=tf.keras.optimizers.Adam(learning_rate=1e-3),
loss=SparseCategoricalCrossentropy(from_logits=True),
)
Binary and multiclass, side by side
Everything above, as one comparison. The right-hand column is the left-hand column with "the other class" replaced by "all the other classes".
| Topic | Binary, $N = 2$ | Multiclass, $N \geq 3$ |
|---|---|---|
| Scores | $z = \vec{w}\cdot\vec{x} + b$, one of them | $z_k = \vec{w}_k\cdot\vec{x} + b_k$, one per class |
| Modelled quantity | the log-odds, $\log\frac{p_1}{p_0} = z$ | pairwise log-odds, $\log\frac{p_k}{p_j} = z_k - z_j$ |
| Scores to probabilities | sigmoid, $p_1 = \frac{1}{1 + e^{-z}}$ | softmax, $p_k = \frac{e^{z_k}}{\sum_j e^{z_j}}$ |
| Likelihood | Bernoulli, $P(y) = p_1^{y}\left(1-p_1\right)^{1-y}$ | categorical, $P(y) = \prod_k p_k^{y_k}$ |
| Loss | binary cross-entropy, $-\left[y\log p_1 + (1-y)\log p_0\right]$ | categorical cross-entropy, $-\sum_k y_k \log p_k$ |
| Gradient at the scores | $\partial L / \partial z = p - y$ | $\partial L / \partial \vec{z} = \vec{p} - \vec{y}$ |
| Convexity, linear model | convex in $\vec{w}, b$ | convex in $W, b$ |
| Special case | — | $N = 2$ collapses to the sigmoid and the binary loss |
Where this goes next
The classifier is now built and, with Adam, trainable. What is missing is everything about whether it is any good — and that turns out to be a bigger subject than the model itself.
- How to tell a good model from a bad one. Splitting data into training, cross-validation and test sets; diagnosing high bias against high variance; reading learning curves; deciding what to try next when the error is too high; error analysis; when adding data helps and when it does not.
- When accuracy lies. On a skewed dataset — a rare disease, a fraud detector — a model that always predicts "no" scores $99.5\%$. Precision, recall and $F_1$ are the metrics that survive that, and the threshold between them is a genuine design decision.
- A different kind of model entirely. Decision trees split on features rather than weighting them, and are trained by maximising information gain rather than by gradient descent. Then ensembles of them, random forests, and XGBoost.
Backwards from here: logistic-regression.md is the two-class model this note generalised, linear-regression.md has the gradient descent and feature-scaling machinery all of it runs on, and neural-networks-reference.md covers the hidden layers that sit underneath the softmax output.