Logistic regression and softmax classification

Why not regress on the label?

If labels are $y\in\{0,1\}$ you could fit a line to them, but a line predicts values below 0 and above 1, and far-away points drag it around even when they are already classified correctly. What you want is a model of the probability $p(y=1\mid\mathbf{x})$ that stays in $(0,1)$ and whose loss stops caring about a point once it is confidently right. Logistic regression is the minimal such model.

Sigmoid and log-odds

Keep the linear score $z = \mathbf{w}^\top\mathbf{x} + b$ and squash it with the sigmoid (logistic function):

Inverting it explains the model: $\log\frac{p}{1-p} = z = \mathbf{w}^\top\mathbf{x}+b$. The log-odds (logit) is linear in the features, so each unit increase in $x_j$ multiplies the odds $p/(1-p)$ by $e^{w_j}$ — that is how you interpret a coefficient. The decision boundary $p = 0.5$ is $z = 0$, a hyperplane: logistic regression is a linear classifier; the nonlinearity only shapes the probabilities around that plane. In the demo, set the hidden layer to zero and you are looking at exactly this boundary; add a hidden layer and you will see what one more layer buys (4.1 From perceptron to MLP; universal approximation).

The loss: cross-entropy from a Bernoulli likelihood

Model each label as a coin flip with probability $p_i = \sigma(z_i)$: $p(y_i\mid\mathbf{x}_i) = p_i^{y_i}(1-p_i)^{1-y_i}$. The likelihood of the dataset is the product over $i$, and its negative log is the binary cross-entropy (1.9 Information theory: entropy, cross-entropy and KL divergence):

Each term is the surprise of the true label under the model: predicting $p=0.99$ for a positive costs $-\log 0.99 \approx 0.01$; predicting $p=0.01$ for it costs $-\log 0.01 \approx 4.6$. The loss is unbounded above, which is why wrong-and-confident is punished so hard, and why numerically you compute it from the logit as $\log(1+e^{-z})$ for $y=1$ (the "with logits" losses in 2.4 PyTorch I: tensors, autograd and the training loop) rather than from $p$, which can round to exactly 0 or 1.

The gradient (and why it looks familiar)

Differentiate one term with respect to $p$: $\frac{\partial}{\partial p}\big[-y\log p-(1-y)\log(1-p)\big] = -\frac{y}{p} + \frac{1-y}{1-p}$. Multiply by $\partial p/\partial z = p(1-p)$:

The messy pieces cancel to the cleanest possible result: the error signal at the logit is prediction minus truth. Chain through $z_i = \mathbf{w}^\top\mathbf{x}_i + b$ and stack over examples:

Compare 3.2 Linear regression: linear regression's gradient was $\frac{2}{n}\mathbf{X}^\top(\hat{\mathbf{y}}-\mathbf{y})$. Same shape — "features times errors" — with the prediction passed through a sigmoid. This is not a coincidence: both losses are negative log-likelihoods of exponential-family models, and for all of them the gradient is $\mathbf{X}^\top(\text{prediction}-\text{target})$.

There is no closed form because $\hat{\mathbf{p}}$ depends on $\mathbf{w}$ nonlinearly, so we iterate. But the problem is convex: the Hessian is $\frac{1}{n}\mathbf{X}^\top\mathbf{S}\mathbf{X}$ with $\mathbf{S} = \operatorname{diag}\big(p_i(1-p_i)\big)$, which is positive semidefinite (1.3 Eigen-things, SVD and low rank), so gradient descent, or Newton's method (called iteratively reweighted least squares here, and what scikit-learn's default solver approximates) reaches the global optimum. One caveat: if the classes are perfectly separable, the likelihood keeps improving as $\|\mathbf{w}\|\to\infty$ and the weights diverge; this is why every practical implementation adds a small L2 penalty by default (3.6 Regularization and generalization; scikit-learn's C is $1/\lambda$).

Worked example

Let $\mathbf{w} = (1, -2)$, $b = 0.5$, and one example $\mathbf{x} = (2, 1)$ with label $y = 1$.

• Logit: $z = 1\cdot 2 + (-2)\cdot 1 + 0.5 = 0.5$.
• Probability: $p = \sigma(0.5) = 1/(1+e^{-0.5}) \approx 0.622$.
• Loss: $-\log 0.622 \approx 0.474$.
• Error signal: $p - y = -0.378$. Gradient w.r.t. $\mathbf{w}$: $(p-y)\,\mathbf{x} = (-0.755, -0.378)$; w.r.t. $b$: $-0.378$.
• One gradient step with $\eta = 0.5$: $\mathbf{w} \leftarrow (1.378, -1.811)$, $b \leftarrow 0.689$; new logit $\approx 1.634$, new $p \approx 0.837$. The step increased the probability of the correct label, as it must.

Many classes: softmax regression

For $K$ classes compute $K$ scores $\mathbf{z} = \mathbf{W}\mathbf{x} + \mathbf{b}$ ($\mathbf{W}\in\mathbb{R}^{K\times d}$) and normalize with the softmax:

Softmax is shift-invariant ($\mathbf{z} + c\mathbf{1}$ gives the same probabilities — subtract the max for stability, 2.2 NumPy and vectorized thinking), reduces to the sigmoid when $K=2$, and its gradient has the same beautiful form: $\partial\ell/\partial\mathbf{z} = \mathbf{p} - \mathbf{y}_{\text{onehot}}$ (1.5 Multivariable and matrix calculus derives it). Every classifier head and the output layer of every language model is softmax regression on learned features; the vocabulary of a GPT is a $K \approx 10^5$-way softmax (7.6 Decoder-only: GPT and autoregressive generation).

Thresholds, imbalance and calibration

The output is a probability; turning it into a decision requires a threshold, and $0.5$ is just the default that minimizes 0–1 error when the classes and costs are balanced. If missing a fraud case costs 50× a false alarm, lower the threshold; if positives are 1% of the data, a model can be useful at thresholds far from 0.5. Sweep thresholds and read precision/recall or the ROC curve (3.5 Evaluating models: metrics, validation and the bias–variance trade-off). For heavy imbalance, class weights (reweighting terms in $J$) or resampling shift the fit; thresholds shift the decision. A separate question is calibration: among examples with $p \approx 0.8$, are about 80% positive? Logistic regression is usually well calibrated because it optimizes the likelihood directly; trees and modern neural networks often are not (3.9 Ensembles: bagging, random forests and gradient boosting, 4.3 Loss functions catalog), and the standard fix — Platt scaling — is a logistic regression fitted on the model's scores. Reliability diagrams and the Brier score / log loss measure it.

import numpy as np
def sigmoid(z): return 1 / (1 + np.exp(-z))
def train_logreg(X, y, lr=0.1, steps=5000, l2=1e-3):
    w, b = np.zeros(X.shape[1]), 0.0
    for _ in range(steps):
        p = sigmoid(X @ w + b)
        gw = X.T @ (p - y) / len(y) + l2 * w
        gb = (p - y).mean()
        w -= lr * gw; b -= lr * gb
    return w, b

In PyTorch this is nn.Linear(d, 1) with BCEWithLogitsLoss, or nn.Linear(d, K) with CrossEntropyLoss (2.4 PyTorch I: tensors, autograd and the training loop) — a one-layer network. Insert nn.Linear(d, h) and a ReLU before it and you have a multilayer perceptron (4.1 From perceptron to MLP; universal approximation).

Where this shows up later

• The sigmoid and softmax are the output activations of 4.2 Activation functions; cross-entropy is the loss of 4.3 Loss functions catalog; the gradient $\mathbf{p}-\mathbf{y}$ starts every backward pass in 4.4 Backpropagation, fully derived.
• Softmax over a vocabulary is how a language model predicts the next token — 6.2 Classical NLP and n-gram language models, 7.6 Decoder-only: GPT and autoregressive generation; temperature scaling of the same softmax is decoding — 8.5 Decoding and sampling.
• Logistic regression on frozen features is the linear probe of 5.8 Self-supervised and contrastive learning and 8.11 Interpretability: what is happening inside; the Bradley–Terry reward model of RLHF is logistic regression on preference pairs — 8.3 Alignment: RLHF, PPO, DPO, GRPO and Constitutional AI.
• Thresholds, imbalance and calibration are evaluated formally in 3.5 Evaluating models: metrics, validation and the bias–variance trade-off.