Worked Multinomial Likelihood and Asymptotic Inference

From Category Counts to a Canonical Likelihood

A multinomial model is a basic setting in which likelihood calculations, Fisher information, and asymptotic inference can all be made explicit. The same calculation also shows why parameterization matters: category probabilities have a sum-to-one constraint, while canonical parameters can vary freely over Euclidean space when all probabilities are positive.

Suppose $n$ independent observations fall into one of $K$ categories, with probabilities $p_1,\ldots,p_K$, where $p_k>0$ and $\sum_{k=1}^K p_k=1$. Let $X_k$ be the number in category $k$, so $\sum_{k=1}^K X_k=n$. Use category $K$ as the reference and define $\theta_j=\log(p_j/p_K)$ for $j=1,\ldots,K-1$.

Definition: The canonical parameter for category $j<K$ is its log-odds relative to category $K$, $\theta_j=\log(p_j/p_K)$. The log-partition function is $A(\theta)=\log(1+\sum_{j=1}^{K-1}e^{\theta_j})$. The probabilities are $p_j=e^{\theta_j-A(\theta)}$ for $j<K$ and $p_K=e^{-A(\theta)}$.

For counts $x_1,\ldots,x_K$ summing to $n$, the multinomial probability mass function can therefore be written

P(X_1=x_1,\ldots,X_K=x_K) = \frac{n!}{\prod_{k=1}^K x_k!} \exp\left\{\sum_{j=1}^{K-1}\theta_jx_j-nA(\theta)\right\}.

The factorial term does not depend on $\theta$. Thus the log-likelihood, apart from that term, is $\ell(\theta)=\sum_{j=1}^{K-1}x_j\theta_j-nA(\theta)$. Differentiating $A$ gives $\partial A/\partial\theta_j=p_j$, so the score has components $x_j-np_j$.

Worked Example: Canonical Parameters and the Multinomial MLE

Suppose $n=100$ observations produce counts $(50,30,15,5)$ in four categories. Take category 4 as the reference. The empirical category probabilities are $(0.50,0.30,0.15,0.05)$, so the canonical estimates are

\hat\theta_1=\log(0.50/0.05)=\log 10,\qquad \hat\theta_2=\log(0.30/0.05)=\log 6,\qquad \hat\theta_3=\log(0.15/0.05)=\log 3.

The fitted normalizing denominator is $1+10+6+3=20$, giving fitted probabilities $(10/20,6/20,3/20,1/20)=(0.50,0.30,0.15,0.05)$, as required. The parameter-dependent part of the maximized log-likelihood is

50\log 10+30\log 6+15\log 3-100\log 20 \approx 185.361-299.573=-114.212.

The omitted log-multinomial coefficient is $\log(100!/(50!30!15!5!))$. It is constant in $\theta$, so it does not affect the MLE, though it is part of the full log-likelihood.

Maximum Likelihood and Boundary Counts

The likelihood equations set $x_j-np_j=0$ for each $j<K$. Since the probabilities sum to one, these equations give $p_j=x_j/n$ for all categories, including the reference category. When every count is positive, these probabilities are in the interior of the simplex and correspond to the finite canonical estimates $\hat\theta_j=\log(x_j/x_K)$.

Theorem: For a multinomial sample with $n>0$, the maximum-likelihood estimate over the closed probability simplex is $\hat p_k=X_k/n$. If every $X_k>0$, the corresponding finite canonical MLE is $\hat\theta_j=\log(X_j/X_K)$. If at least one count is zero, the probability MLE remains the empirical distribution, but it is on the boundary and has no finite canonical representation.

To verify the maximization, write the likelihood, up to its constant coefficient, as $\prod_k p_k^{x_k}$. For positive counts, its log is $\sum_k x_k\log p_k$, a strictly concave function on the interior subject to $\sum_kp_k=1$. A Lagrange multiplier gives $x_k/p_k=\lambda$ for every $k$; summing the resulting probabilities yields $\lambda=n$ and $p_k=x_k/n$. Strict concavity makes this the unique maximum. If some counts are zero, the same maximum over the closed simplex is achieved by assigning zero probability to those categories and $x_k/n$ to the others. A zero probability requires a canonical log-odds parameter to diverge, so no finite $\theta$ attains that boundary point.

Worked Example: A Zero Count and a Boundary Estimate

Suppose $n=20$ and the counts in three categories are $(12,8,0)$. The empirical probabilities $(0.6,0.4,0)$ maximize the likelihood over the closed simplex. If category 3 is the reference, the canonical estimates would require $\hat\theta_1=\log(0.6/0)$ and $\hat\theta_2=\log(0.4/0)$, neither of which is finite. The likelihood can approach its maximum as the reference probability tends to zero, but no finite canonical parameter attains it. This is a boundary MLE, not evidence that the empirical probability estimates are undefined.

Fisher Information and Asymptotic Covariance

For one multinomial observation, record the first $K-1$ category indicators as a vector. Its covariance matrix is $V=\operatorname{diag}(p)-pp^\mathsf{T}$, where $p=(p_1,\ldots,p_{K-1})^\mathsf{T}$. Equivalently, this is the Hessian of the log-partition function. The sample information is $nV$, because the negative Hessian of the sample log-likelihood is $nV$.

For interior probabilities, $V$ is positive definite. Indeed, for any nonzero vector $a\in\mathbb R^{K-1}$, $a^\mathsf{T}Va$ is the variance of a random variable that takes values $a_1,\ldots,a_{K-1},0$ with respective probabilities $p_1,\ldots,p_{K-1},p_K$. Since all probabilities are positive, this variance can be zero only if all those values are equal, which forces $a=0$.

Result: With $D=\operatorname{diag}(p_1,\ldots,p_{K-1})$, the per-observation Fisher information for $\theta$ is $V=D-pp^\mathsf{T}$, and $$ V^{-1}=D^{-1}+\frac{1}{p_K}\mathbf{1}\mathbf{1}^\mathsf{T}, $$ where $\mathbf{1}$ is the $(K-1)$-vector of ones. Thus the sample information is $nV$, and its inverse is $V^{-1}/n$.

To check the inverse, note $D^{-1}p=\mathbf{1}$ and $p^\mathsf{T}D^{-1}p=1-p_K$. Multiplication gives

(D-pp^\mathsf{T})\left(D^{-1}+\frac{\mathbf{1}\mathbf{1}^\mathsf{T}}{p_K}\right) = I-p\mathbf{1}^\mathsf{T} +\frac{p\mathbf{1}^\mathsf{T}-(1-p_K)p\mathbf{1}^\mathsf{T}}{p_K} =I.

This proves the stated inverse. $\square$

Worked Example: Information for a Three-Category Model

Let $K=3$, $n=200$, and $p=(0.2,0.5,0.3)$, with category 3 as the reference. Then

V= \begin{pmatrix} 0.2(0.8)&-0.2(0.5)\\ -0.2(0.5)&0.5(0.5) \end{pmatrix} = \begin{pmatrix} 0.16&-0.10\\ -0.10&0.25 \end{pmatrix}.

The sample information is $200V$. Its inverse, the asymptotic covariance matrix for the canonical MLE, is

\frac{V^{-1}}{200} = \frac{1}{200} \begin{pmatrix} 8.3333&3.3333\\ 3.3333&5.3333 \end{pmatrix} \approx \begin{pmatrix} 0.04167&0.01667\\ 0.01667&0.02667 \end{pmatrix}.

For example, the first diagonal entry also follows from $1/p_1+1/p_3=5+3.3333=8.3333$, divided by $n=200$. The positive off-diagonal covariance reflects the shared reference category in the two log-odds estimates.

Under independent sampling with fixed $K$ and true probabilities strictly inside the simplex, the vector of empirical probabilities for the first $K-1$ categories obeys the multivariate central limit theorem. Applying the smooth log-odds map gives

\sqrt n(\hat\theta-\theta)\xrightarrow{d}N_{K-1}(0,V^{-1}).

For completeness, the derivative of $\theta_j=\log(p_j/p_K)$ with respect to $p_l$, $l<K$, is $\delta_{jl}/p_j+1/p_K$. The resulting derivative matrix is $D^{-1}+\mathbf{1}\mathbf{1}^\mathsf{T}/p_K=V^{-1}$. The empirical probability vector has limiting covariance $V$; hence the delta-method covariance for the log-odds is $V^{-1}VV^{-1}=V^{-1}$. The assumptions ensure the true probabilities are positive and the log-odds map is differentiable in a neighborhood of the truth.

A Probability Contrast by the Delta Method

Often the target is a function of probabilities rather than a canonical parameter. Consider $g(p)=p_1-p_2$, with $K\geq3$, and represent the model by its first $K-1$ probabilities. Its gradient is $a=(1,-1,0,\ldots,0)^\mathsf{T}$. The delta method gives limiting variance $a^\mathsf{T}Va$. Direct calculation yields

a^\mathsf{T}Va = p_1+p_2-(p_1-p_2)^2.

Worked Example: Inference for the Difference of Two Probabilities

Return to counts $(50,30,15,5)$ from $n=100$. The estimated contrast between categories 1 and 2 is $\hat p_1-\hat p_2=0.50-0.30=0.20$. At the fitted probabilities, the estimated asymptotic variance per observation is

$$ 0.50+0.30-(0.50-0.30)^2=0.80-0.04=0.76. $$

Therefore $\sqrt{100}\{(\hat p_1-\hat p_2)-(p_1-p_2)\}$ is approximated by a mean-zero normal variable with variance $0.76$. Equivalently, the plug-in standard error for the estimated contrast is $\sqrt{0.76/100}=\sqrt{0.0076}\approx0.0872$. This quantifies sampling uncertainty in the difference, not in either probability considered separately.

Interpretation and a Common Pitfall

The multinomial model connects three useful representations of uncertainty. Empirical proportions have covariance $V/n$, canonical log-odds have asymptotic covariance $V^{-1}/n$, and a smooth probability contrast has variance obtained by applying its gradient to $V/n$. These are different coordinate descriptions of the same sampling problem.

Pitfall: treating all $K$ probabilities as free parameters. Their sum is fixed at one, so the full $K$-dimensional covariance matrix is singular. Use $K-1$ free probabilities or $K-1$ reference-category log-odds.
Pitfall: using finite log-odds when a category count is zero. The empirical probability MLE still exists, but a zero fitted probability is on the boundary and cannot be represented by finite canonical parameters.
Pitfall: confusing probability and canonical covariance. $V/n$ describes the first $K-1$ empirical probabilities; $V^{-1}/n$ describes the canonical MLE. They are not interchangeable.

Check Your Understanding

Work through these questions using the reference-category parameterization.

For probabilities $(0.4,0.35,0.25)$ with category 3 as reference, what are the two canonical parameters?
Why does a zero count prevent a finite canonical MLE even though the empirical probability estimate is defined?
For $K=3$, write the per-observation information matrix in terms of $p_1,p_2$, and $p_3$.
What is the asymptotic variance per observation of $\hat p_1-\hat p_2$ when the true probabilities are $(0.2,0.5,0.3)$?
Explain why the covariance matrix of all $K$ empirical proportions is singular.