With the softmax function in our toolkit, we can now construct the complete multinomial logistic regression model—a principled extension of binary logistic regression to any number of classes. This model is foundational: it underpins the output layers of virtually every classification neural network and remains a powerful baseline for multi-class problems.

Multinomial logistic regression is known by several names:

• Softmax regression — emphasizing the softmax output transformation
• Maximum entropy classifier — from its information-theoretic derivation
• Multiclass logit / Polytomous logit — from econometrics
• Multinomial logit — the formal statistical name

Despite this nomenclature diversity, all these terms describe the same model. In this page, we develop multinomial logistic regression completely—from model specification through parameterization choices, geometric interpretation, and practical considerations.

The Setting

We have:

• Features: $\mathbf{x} \in \mathbb{R}^d$ — a $d$-dimensional feature vector
• Labels: $y \in {1, 2, \ldots, K}$ — one of $K$ possible classes
• Training data: ${(\mathbf{x}_1, y_1), (\mathbf{x}_2, y_2), \ldots, (\mathbf{x}_n, y_n)}$

Our goal is to model the conditional probability distribution $P(y | \mathbf{x})$ over all $K$ classes given features $\mathbf{x}$.

The Model

Multinomial logistic regression models the conditional probabilities as:

$$P(y = k | \mathbf{x}; \boldsymbol{\Theta}) = \frac{\exp(z_k)}{\sum_{j=1}^{K} \exp(z_j)} = \text{softmax}(\mathbf{z})_k$$

where the logit (or score) for class $k$ is a linear function of features:

$$z_k = \mathbf{w}k^T \mathbf{x} + b_k = \sum{i=1}^{d} w_{ki} x_i + b_k$$

Here:

• $\mathbf{w}_k \in \mathbb{R}^d$ is the weight vector for class $k$
• $b_k \in \mathbb{R}$ is the bias (intercept) for class $k$
• $\boldsymbol{\Theta} = {(\mathbf{w}_1, b_1), \ldots, (\mathbf{w}_K, b_K)}$ is the full parameter set

The Modeling Assumption

Multinomial logistic regression makes a specific assumption about how class probabilities depend on features: the log-odds between any two classes is a linear function of features.

For classes $k$ and $l$:

$$\log \frac{P(y=k|\mathbf{x})}{P(y=l|\mathbf{x})} = \log \frac{e^{z_k}}{e^{z_l}} = z_k - z_l = (\mathbf{w}_k - \mathbf{w}_l)^T \mathbf{x} + (b_k - b_l)$$

This is the log-linear model assumption: log-odds are linear in features. The assumption is strong but often reasonable and leads to a model with many desirable properties:

1. Interpretable parameters: Each weight $w_{ki}$ represents how feature $i$ affects the log-odds of class $k$
2. Convex optimization: Maximum likelihood estimation is a convex problem
3. No distributional assumptions on X: Unlike discriminant analysis, we don't assume feature distributions
4. Calibrated probabilities: Under correct specification, predicted probabilities match true frequencies

Formal Distributional Form

The model specifies that $y | \mathbf{x}$ follows a categorical distribution (generalization of Bernoulli to $K$ outcomes):

$$y | \mathbf{x}; \boldsymbol{\Theta} \sim \text{Categorical}(p_1, p_2, \ldots, p_K)$$

where $p_k = P(y=k|\mathbf{x}) = \text{softmax}(\mathbf{z})_k$.

Equivalently, writing $y$ as a one-hot vector $\mathbf{y} \in {0,1}^K$ with $y_k = 1$ iff $y = k$:

$$P(\mathbf{y}|\mathbf{x}; \boldsymbol{\Theta}) = \prod_{k=1}^{K} p_k^{y_k} = \prod_{k=1}^{K} \left(\frac{e^{z_k}}{\sum_j e^{z_j}}\right)^{y_k}$$

This product form is crucial for likelihood-based learning.

Due to the translation invariance of softmax, multinomial logistic regression is fundamentally overparameterized: adding a constant to all logits doesn't change the probability distribution. This creates both opportunities and complications.

The Identifiability Problem

Consider parameters $\boldsymbol{\Theta} = {(\mathbf{w}k, b_k)}{k=1}^K$ and shifted parameters $\boldsymbol{\Theta}' = {(\mathbf{w}k + \mathbf{c}, b_k + d)}{k=1}^K$ for any $\mathbf{c} \in \mathbb{R}^d$, $d \in \mathbb{R}$:

$$z'_k = z_k + \mathbf{c}^T\mathbf{x} + d$$

Since the shift is identical for all classes: $$\text{softmax}(\mathbf{z}') = \text{softmax}(\mathbf{z})$$

These infinitely many parameter settings produce identical predictions. The model is not identifiable without constraints.

Resolution 1: Reference Class Parameterization

Fix one class (say class $K$) as the reference with $\mathbf{w}_K = \mathbf{0}$ and $b_K = 0$:

$$z_K = 0 \quad \text{(always)}$$ $$z_k = \mathbf{w}_k^T \mathbf{x} + b_k \quad \text{for } k = 1, \ldots, K-1$$

Probabilities under reference class parameterization:

$$P(y=k|\mathbf{x}) = \frac{e^{z_k}}{1 + \sum_{j=1}^{K-1} e^{z_j}} \quad \text{for } k = 1, \ldots, K-1$$

$$P(y=K|\mathbf{x}) = \frac{1}{1 + \sum_{j=1}^{K-1} e^{z_j}}$$

Note the similarity to binary logistic regression:

• Binary: $P(y=1) = \frac{e^z}{1+e^z} = \sigma(z)$ — exactly the $K=2$ case

Advantages:

1. Identifiability: Parameters uniquely determined (no redundancy)
2. Direct interpretation: $\mathbf{w}_k$ represents log-odds of class $k$ vs. reference class $K$
3. Fewer parameters: $(K-1)(d+1)$ instead of $K(d+1)$

Disadvantages:

1. Asymmetric interpretation: Reference class is special, others are relative
2. Coefficient comparison: $\mathbf{w}_k - \mathbf{w}_l$ for non-reference classes requires extra computation
3. Reference choice matters for interpretation (but not for predictions)

Resolution 2: Full Parameterization with Regularization

Alternatively, keep all $K$ parameter sets but add regularization that implicitly breaks symmetry:

$$\mathcal{L}{\text{regularized}} = \mathcal{L}{\text{NLL}} + \lambda \sum_{k=1}^{K} |\mathbf{w}_k|^2$$

The L2 penalty prefers solutions where parameters are small, effectively centering them around zero. While still technically overparameterized, regularization:

1. Selects a unique solution (the minimum-norm one among equivalent solutions)
2. Improves generalization through the usual regularization benefits
3. Maintains symmetric treatment of all classes

Neural Network Convention:

Deep learning frameworks typically use full parameterization:

• Final layer: $\mathbf{z} = W\mathbf{h} + \mathbf{b}$ where $W \in \mathbb{R}^{K \times d_{\text{hidden}}}$
• Softmax: $\mathbf{p} = \text{softmax}(\mathbf{z})$
• Weight decay applies to $W$

This symmetric treatment is convenient for implementation and doesn't cause optimization issues thanks to regularization.

The log-linear structure of multinomial logistic regression enables rich interpretation of parameters. Understanding log-odds is essential for model diagnostics and domain insights.

Log-Odds Definition

The log-odds (or logit) of class $k$ versus class $l$ is:

$$\log \frac{P(y=k|\mathbf{x})}{P(y=l|\mathbf{x})} = z_k - z_l = (\mathbf{w}_k - \mathbf{w}_l)^T \mathbf{x} + (b_k - b_l)$$

This is the natural logarithm of the odds ratio—how many times more likely class $k$ is than class $l$.

Interpreting Coefficients (Reference Class)

With reference class $K$ (where $\mathbf{w}_K = 0$, $b_K = 0$):

$$\log \frac{P(y=k|\mathbf{x})}{P(y=K|\mathbf{x})} = \mathbf{w}_k^T \mathbf{x} + b_k$$

For a unit increase in feature $x_i$ (holding others constant):

$$\log \frac{P(y=k|\mathbf{x} + \mathbf{e}_i)}{P(y=K|\mathbf{x} + \mathbf{e}i)} - \log \frac{P(y=k|\mathbf{x})}{P(y=K|\mathbf{x})} = w{ki}$$

Thus:

• $w_{ki} > 0$: Unit increase in $x_i$ increases log-odds of class $k$ vs. reference by $w_{ki}$. The odds of class $k$ are multiplied by $e^{w_{ki}}$.
• $w_{ki} < 0$: Unit increase in $x_i$ decreases log-odds of class $k$ vs. reference.
• $w_{ki} = 0$: Feature $x_i$ doesn't affect relative probability of class $k$ vs. reference.

Odds Ratio Interpretation

Exponentiating gives the odds ratio:

$$\text{OR}{ki} = e^{w{ki}}$$

For a unit increase in feature $x_i$:

• $\text{OR}{ki} > 1$: Odds of class $k$ (vs. reference) multiply by $\text{OR}{ki}$
• $\text{OR}{ki} < 1$: Odds of class $k$ (vs. reference) divide by $1/\text{OR}{ki}$
• $\text{OR}_{ki} = 1$: Feature has no effect on relative class probability

Example:

Suppose we're classifying customer churn into: Stayed (reference), Churned to Competitor A, Churned to Competitor B.

If the coefficient for 'years as customer' for 'Churned to A' is $w = -0.3$:

$$\text{OR} = e^{-0.3} \approx 0.74$$

Interpretation: For each additional year as a customer, the odds of churning to Competitor A (versus staying) multiply by 0.74—i.e., they decrease by 26%.

If the same coefficient for 'Churned to B' is $w = -0.1$:

$$\text{OR} = e^{-0.1} \approx 0.90$$

Tenure reduces churn to B by only 10% per year—less protective than against A.

Comparing Non-Reference Classes

To compare classes $k$ and $l$ (neither is reference):

$$\log \frac{P(y=k|\mathbf{x})}{P(y=l|\mathbf{x})} = (\mathbf{w}_k - \mathbf{w}_l)^T \mathbf{x} + (b_k - b_l)$$

The effect of feature $i$ on log-odds of $k$ vs. $l$ is $w_{ki} - w_{li}$.

Statistical Inference:

To test whether classes $k$ and $l$ differ in their relationship with feature $i$:

• Null hypothesis: $H_0: w_{ki} = w_{li}$
• Test statistic: Wald test on $(w_{ki} - w_{li})$ requires estimating the covariance between $\hat{w}{ki}$ and $\hat{w}{li}$

This is more complex than testing single coefficients and often requires simultaneous inference adjustments (Bonferroni, etc.) when comparing multiple pairs.

The geometry of multinomial logistic regression reveals a beautiful structure: the feature space is partitioned into regions by linear decision boundaries (hyperplanes).

Pairwise Decision Boundaries

The boundary between classes $k$ and $l$ occurs where their probabilities are equal:

$$P(y=k|\mathbf{x}) = P(y=l|\mathbf{x})$$

This implies $z_k = z_l$, so:

$$(\mathbf{w}_k - \mathbf{w}_l)^T \mathbf{x} + (b_k - b_l) = 0$$

This is a hyperplane in $\mathbb{R}^d$:

• Normal vector: $\mathbf{w}_k - \mathbf{w}_l$
• Offset: $-(b_k - b_l)$

For $K$ classes, there are $\binom{K}{2} = \frac{K(K-1)}{2}$ pairwise boundaries.

The Voronoi-like Partition

Even though there are $\binom{K}{2}$ pairwise boundaries, the actual decision regions are simpler. Each class $k$ 'owns' the region where it has the highest logit:

$$\text{Predict class } k \text{ where } z_k > z_j \text{ for all } j \neq k$$

This creates a partition of feature space into (at most) $K$ convex polytopes—regions bounded by intersecting hyperplanes.

Properties of the Decision Regions:

1. Convexity: Each class region is convex (intersection of half-spaces)
2. Simply connected: No class region is disconnected
3. May be unbounded: Regions can extend to infinity
4. May be empty: With certain parameters, some classes may never be predicted

Soft vs. Hard Boundaries:

The model outputs probability distributions, not hard predictions. Near a decision boundary:

• Probabilities are close to uniform (high uncertainty)
• Moving perpendicular to the boundary shifts probability smoothly
• Temperature scaling controls the sharpness: high T → gradual transitions, low T → sharp boundaries

Binary logistic regression is not merely a special case of multinomial logistic regression—it is multinomial logistic regression with $K=2$. Let's verify this correspondence precisely.

Binary Setup Recap

In binary logistic regression with classes ${0, 1}$:

$$P(y=1|\mathbf{x}) = \sigma(\mathbf{w}^T \mathbf{x} + b) = \frac{1}{1 + e^{-(\mathbf{w}^T \mathbf{x} + b)}}$$ $$P(y=0|\mathbf{x}) = 1 - P(y=1|\mathbf{x}) = \frac{1}{1 + e^{\mathbf{w}^T \mathbf{x} + b}}$$

Parameters: Single weight vector $\mathbf{w} \in \mathbb{R}^d$ and bias $b \in \mathbb{R}$.

Multinomial with $K=2$

Using reference class parameterization with class 0 as reference:

• $z_0 = 0$ (reference)
• $z_1 = \mathbf{w}_1^T \mathbf{x} + b_1$

Then: $$P(y=1|\mathbf{x}) = \frac{e^{z_1}}{e^{z_0} + e^{z_1}} = \frac{e^{z_1}}{1 + e^{z_1}} = \frac{1}{1 + e^{-z_1}}$$

With $\mathbf{w}_1 = \mathbf{w}$ and $b_1 = b$, this is exactly the sigmoid formula.

Similarly: $$P(y=0|\mathbf{x}) = \frac{1}{1 + e^{z_1}} = \frac{1}{1 + e^{\mathbf{w}^T \mathbf{x} + b}}$$

The correspondence is exact.

Gradient Correspondence

The gradient structures also correspond. For binary logistic regression with $y \in {0, 1}$:

$$\frac{\partial \mathcal{L}}{\partial \mathbf{w}} = \frac{1}{n} \sum_{i=1}^{n} (\hat{p}_i - y_i) \mathbf{x}_i$$

For multinomial with one-hot $\mathbf{y}$ and class $k$ for sample $i$:

$$\frac{\partial \mathcal{L}}{\partial \mathbf{w}k} = \frac{1}{n} \sum{i=1}^{n} (\hat{p}{ik} - y{ik}) \mathbf{x}_i$$

The form $$(\text{predicted probability} - \text{true label}) \times \text{features}$$ is universal. This elegant gradient structure arises from the canonical link function property of logistic regression within the generalized linear model framework.

Understanding when multinomial logistic regression is appropriate—and when it fails—requires examining its underlying assumptions.

Assumption 1: Log-linear Relationship

The model assumes log-odds are linear functions of features:

$$\log \frac{P(y=k|\mathbf{x})}{P(y=l|\mathbf{x})} = \text{linear in } \mathbf{x}$$

This is violated when:

• Feature effects are nonlinear (quadratic, threshold, etc.)
• Features interact in complex ways
• The decision boundary is inherently curved

Remedies:

• Feature engineering: Add polynomial features, interactions
• Kernel methods: Implicit feature expansion
• Neural networks: Learned nonlinear feature extraction

Assumption 2: Independence of Irrelevant Alternatives (IIA)

This famous property (from econometrics) states that the ratio of probabilities for any two classes doesn't depend on other classes:

$$\frac{P(y=k|\mathbf{x})}{P(y=l|\mathbf{x})} = \frac{e^{z_k}}{e^{z_l}} = e^{z_k - z_l}$$

This ratio only involves $k$ and $l$—adding or removing class $m$ doesn't affect it.

When IIA is Problematic:

• Classes have 'nested' structure (buses vs. car)
• Classes share unobserved attributes
• Adding a new class should mostly draw from similar existing classes

Remedies for IIA Violation:

• Nested logit: Hierarchical choice structure
• Mixed logit: Random coefficients allowing correlation
• Probit models: Correlated error terms via covariance matrix

Assumption 3: No Complete Separation

Complete separation occurs when a hyperplane perfectly separates all instances of some classes. When this happens:

• MLEs don't exist: coefficients → ±∞
• Predicted probabilities → 0 or 1 exactly
• Standard errors are undefined

Detection:

• Warnings from optimization (convergence failure)
• Very large coefficient estimates
• Predicted probabilities at 0 or 1

Remedies:

• Regularization (L1/L2): Most common, prevents coefficient explosion
• Firth's penalized likelihood: Bias-reduced estimates
• Exact logistic regression: For small samples
• Bayesian methods: Proper priors prevent unbounded posteriors

Assumption 4: Correctly Specified Linear Predictor

The model assumes we've included the 'right' features in the linear predictor. Misspecification leads to:

• Biased coefficient estimates
• Poor calibration (predicted probabilities don't match observed frequencies)
• Reduced predictive power

Remedies:

• Model selection: AIC, BIC, cross-validation
• Regularization: Lasso for feature selection
• Residual diagnostics: Check for systematic patterns

Practical Robustness:

Despite these assumptions, multinomial logistic regression is remarkably robust in practice:

• Often performs well even with mild violations
• Regularization handles many edge cases
• Provides interpretable baselines
• Fast to train on large datasets

It remains the go-to method for multi-class classification before trying more complex models.

Implementing multinomial logistic regression correctly requires attention to several practical details. Let's examine production-quality patterns.

Data Preparation

1. Feature scaling: Unlike tree-based methods, logistic regression benefits from standardized features (mean 0, std 1). This helps optimization converge faster and makes coefficients comparable.

2. One-hot encoding: Categorical features must be one-hot encoded. Drop one level per feature to avoid multicollinearity (or rely on regularization).

3. Label encoding: Convert string labels to integers ${0, 1, \ldots, K-1}$. For loss computation, convert to one-hot format.

4. Train/validation split: Essential for hyperparameter tuning (regularization strength).

We have developed multinomial logistic regression as the complete framework for probabilistic multi-class classification. Let's consolidate the key insights:

What's Next:

With the model structure established, we now turn to the crucial question: How do we train this model? The next page develops the cross-entropy loss function, deriving it from maximum likelihood principles and understanding why it's the natural choice for training softmax-based classifiers.