Plug-In LDA Classification From Training Samples

From Labeled Training Data to a Classification Rule

In practice, the population means and shared covariance required by linear discriminant analysis (LDA) are usually unknown. A labeled training sample supplies estimates: average the observations within each class to estimate its mean, and combine the within-class variation to estimate the common covariance. The resulting rule is a plug-in classifier: it substitutes these estimates, and a chosen prior for each class, into the LDA score.

The pooled covariance must reflect variation within classes. Combining all observations around one overall mean would also include differences between class means, which is a different quantity and generally does not estimate the shared within-class covariance.

Definition: Suppose there are $K$ classes and $p$-dimensional labeled training observations $x_{ki}$, where $k=1,\ldots,K$ indexes the class and $i=1,\ldots,n_k$ indexes observations within that class. Write $N=\sum_{k=1}^K n_k$. The class sample mean is $$ \hat{\mu}_k=\frac{1}{n_k}\sum_{i=1}^{n_k}x_{ki}. $$ The pooled within-class covariance estimate is $$ \hat{\Sigma}_p =\frac{1}{N-K}\sum_{k=1}^K\sum_{i=1}^{n_k} (x_{ki}-\hat{\mu}_k)(x_{ki}-\hat{\mu}_k)'. $$ When every class has at least two observations, this is equivalently $\hat{\Sigma}_p=\sum_{k=1}^K(n_k-1)S_k/(N-K)$, where $S_k$ is the usual sample covariance within class $k$.

The denominator $N-K$ counts the within-class degrees of freedom: estimating each class mean uses one degree of freedom per class. This is the unbiased pooled estimator under a common covariance model. It is not the only possible convention; a Gaussian maximum-likelihood estimate uses denominator $N$. In the calculations here, the pooled estimate uses $N-K$.

The standard plug-in choice for the class prior is the observed class proportion, $\hat{\pi}_k=n_k/N$. This choice treats the training sample proportions as representative of the class proportions relevant to classification. If target-population priors are known or deliberately specified, those can instead be used; class counts alone do not determine the scientifically appropriate prior in every study.

The Plug-In LDA Score

Assume the class distributions are modeled as multivariate normal with class-specific means and one common positive-definite covariance matrix, and that each class has positive prior probability. For any new $p$-vector $x$, the LDA discriminant score, with terms common to all classes omitted, is

Definition: Given an invertible pooled estimate $\hat{\Sigma}_p$ and positive prior estimates $\hat{\pi}_k$, the plug-in LDA score for class $k$ is $$ \hat{\delta}_k(x) =x'\hat{\Sigma}_p^{-1}\hat{\mu}_k -\frac{1}{2}\hat{\mu}_k'\hat{\Sigma}_p^{-1}\hat{\mu}_k +\log\hat{\pi}_k. $$ The plug-in LDA rule assigns $x$ to a class with the largest score. If scores tie, the rule needs an explicit tie-breaking convention.

The first two terms are the estimated Gaussian log-density terms that vary by class when the covariance is shared; the final term incorporates the prior. In particular, the prior affects the score even though it is not part of the covariance or mean estimation.

Result: Under the common-covariance Gaussian model, substituting estimated means, a shared covariance estimate, and positive priors into the class log-density gives the score above. When the covariance estimate is invertible, maximizing these scores is equivalent to minimizing $$ (x-\hat{\mu}_k)'\hat{\Sigma}_p^{-1}(x-\hat{\mu}_k)-2\log\hat{\pi}_k. $$

Derivation. For a class with mean $\mu_k$, shared positive-definite covariance $\Sigma$, and prior $\pi_k$, its log joint density for $x$, up to no omissions yet, is

\log\{\pi_k f_k(x)\} =\log\pi_k-\frac{p}{2}\log(2\pi)-\frac{1}{2}\log|\Sigma| -\frac{1}{2}(x-\mu_k)'\Sigma^{-1}(x-\mu_k).

The terms $-\frac{p}{2}\log(2\pi)$ and $-\frac{1}{2}\log|\Sigma|$ are the same for every class, so they do not affect which score is largest. Expanding the remaining quadratic gives

-\frac{1}{2}(x-\mu_k)'\Sigma^{-1}(x-\mu_k)+\log\pi_k =-\frac{1}{2}x'\Sigma^{-1}x +x'\Sigma^{-1}\mu_k-\frac{1}{2}\mu_k'\Sigma^{-1}\mu_k+\log\pi_k.

The term $-\frac{1}{2}x'\Sigma^{-1}x$ is also common to all classes. Removing it leaves the stated score. Replacing the unknown quantities by $\hat{\mu}_k$, $\hat{\Sigma}_p$, and $\hat{\pi}_k$ gives the plug-in rule. Finally, expanding the squared distance in the result gives $x'\hat{\Sigma}_p^{-1}x-2x'\hat{\Sigma}_p^{-1}\hat{\mu}_k+\hat{\mu}_k'\hat{\Sigma}_p^{-1}\hat{\mu}_k$. Its first term is class-independent; subtracting one half of the remaining terms and adding $\log\hat{\pi}_k$ yields $\hat{\delta}_k(x)$. $\square$

Worked Calculations

Worked Example: A Two-Class Classification With Unequal Sample Sizes

Consider two-dimensional observations. Class A has training points $(2,1)'$, $(0,1)'$, and $(1,1)'$; class B has points $(3,3)'$ and $(3,1)'$. Class A's mean is $\hat{\mu}_A=(1,1)'$, and class B's mean is $\hat{\mu}_B=(3,2)'$. The class sizes are 3 and 2, so the empirical priors are $3/5$ and $2/5$.

For class A, the deviations from its mean are $(1,0)'$, $(-1,0)'$, and $(0,0)'$. Their outer products sum to $\begin{pmatrix}2&0\\0&0\end{pmatrix}$. For class B, the deviations are $(0,1)'$ and $(0,-1)'$, whose outer products sum to $\begin{pmatrix}0&0\\0&2\end{pmatrix}$. Therefore the total within-class sum of outer products is $2I$. Since $N-K=5-2=3$,

\hat{\Sigma}_p=\frac{2}{3}I,\qquad \hat{\Sigma}_p^{-1}=\frac{3}{2}I.

Classify $x=(2,1)'$. For A, the linear term is $x'\hat{\Sigma}_p^{-1}\hat{\mu}_A=2(1.5)+1(1.5)=4.5$. The quadratic term is $\frac{1}{2}\hat{\mu}_A'\hat{\Sigma}_p^{-1}\hat{\mu}_A=\frac{1}{2}(1.5+1.5)=1.5$. Thus $\hat{\delta}_A(x)=4.5-1.5+\log(3/5)=3+\log(3/5)\approx2.489$.

For B, the linear term is $x'\hat{\Sigma}_p^{-1}\hat{\mu}_B=2(4.5)+1(3)=12$. Its quadratic term is $\frac{1}{2}(3^2+2^2)(1.5)=9.75$. Hence $\hat{\delta}_B(x)=12-9.75+\log(2/5)=2.25+\log(2/5)\approx1.334$. Since A's score is larger, $x$ is assigned to class A. The calculation uses the full matrix product in the linear term; omitting the covariance inverse would give a different, incorrect score.

Worked Example: How Empirical Priors Can Change the Assignment

Use one-dimensional training observations. Class A consists of $0,1,2,3$, and class B consists of $3,5$. Their sample means are $\hat{\mu}_A=1.5$ and $\hat{\mu}_B=4$. For A, the sum of squared deviations is $5$; for B, it is $2$. There are $N=6$ observations and $K=2$ classes, so

\hat{\Sigma}_p=\frac{5+2}{6-2}=\frac{7}{4}, \qquad \hat{\Sigma}_p^{-1}=\frac{4}{7}.

The empirical priors are $\hat{\pi}_A=4/6=2/3$ and $\hat{\pi}_B=2/6=1/3$. At $x=3$, A's linear term is $3(4/7)(1.5)=18/7$, and its quadratic term is $\frac{1}{2}(4/7)(1.5)^2=9/14$. Therefore $\hat{\delta}_A(3)=18/7-9/14+\log(2/3)=27/14+\log(2/3)\approx1.523$.

For B, the linear term is $3(4/7)(4)=48/7$, and the quadratic term is $\frac{1}{2}(4/7)(4^2)=32/7$. Thus $\hat{\delta}_B(3)=48/7-32/7+\log(1/3)=16/7+\log(1/3)\approx1.187$. The empirical-prior rule assigns $x$ to A. If the priors were equal, the log-prior terms would cancel, and B's score would exceed A's by $16/7-27/14=5/14>0$. Thus the unequal empirical priors reverse the assignment in this example.

This comparison is not a claim that empirical priors are always preferable. It shows precisely where the class proportions enter and why a change in the prior choice can change the result.

Worked Example: Scoring Three Classes

In one dimension, let the training data be A: $0,2$; B: $3,4,5$; and C: $7,9$. The class means are $1,4,8$. Each class has within-class sum of squares $2$, so the total is $6$. Here $N=7$, $K=3$, and the pooled variance is $6/(7-3)=3/2$, with inverse $2/3$. The empirical priors are $2/7,3/7,2/7$.

For the new observation $x=4$, the score for A is $\hat{\delta}_A(4)=4(2/3)(1)-\frac{1}{2}(2/3)(1)^2+\log(2/7)=7/3+\log(2/7)\approx1.081$. For B it is $\hat{\delta}_B(4)=4(2/3)(4)-\frac{1}{2}(2/3)(4)^2+\log(3/7)=16/3+\log(3/7)\approx4.486$. For C it is $\hat{\delta}_C(4)=4(2/3)(8)-\frac{1}{2}(2/3)(8)^2+\log(2/7)=\log(2/7)\approx-1.253$. The largest score is B's, so the observation is assigned to class B.

Why Pooling Is Appropriate, and When It Can Fail

The pooled estimate combines within-class sums of squares and cross-products while allowing each class to have its own mean. Under a common covariance model, this makes it more stable than estimating each class covariance separately, because it uses the within-class information across all classes. The following result explains its unbiasedness under the common-covariance assumption.

Result: If, within each class $k$, the observations are independent with mean $\mu_k$ and covariance $\Sigma$, then $E[\hat{\Sigma}_p]=\Sigma$, provided $N>K$.

Proof. For a fixed class, expand the centered sum as $\sum_i(x_{ki}-\hat{\mu}_k)(x_{ki}-\hat{\mu}_k)'=\sum_i x_{ki}x_{ki}'-n_k\hat{\mu}_k\hat{\mu}_k'$. Independence and the common covariance imply $E[\sum_i x_{ki}x_{ki}']=n_k(\Sigma+\mu_k\mu_k')$. Also, $\operatorname{Var}(\hat{\mu}_k)=\Sigma/n_k$, so $E[\hat{\mu}_k\hat{\mu}_k']=\mu_k\mu_k'+\Sigma/n_k$. The expected centered sum for class $k$ is therefore $n_k(\Sigma+\mu_k\mu_k')-n_k(\mu_k\mu_k'+\Sigma/n_k)=(n_k-1)\Sigma$. Summing over classes gives $(N-K)\Sigma$; dividing by $N-K$ proves the result. $\square$

Unbiasedness does not ensure that the estimate can be inverted. The score requires $\hat{\Sigma}_p^{-1}$, so the pooled matrix must be nonsingular. In particular, the within-class centered vectors span a subspace of dimension at most $N-K$; if $p>N-K$, the pooled estimate cannot have full rank. Even when it is invertible, near-singularity can make scores sensitive to small data changes. The remedy depends on the modeling goal and may involve additional assumptions or regularization; it is not achieved by silently treating a singular matrix as invertible.

Common Mistakes

Pooling around the overall mean. This includes between-class separation in the covariance calculation. Center each training observation at its own class mean before pooling.
Using the wrong denominator without saying so. The pooled unbiased estimate divides by $N-K$; the Gaussian maximum-likelihood estimate divides by $N$. State which estimate is being used and use it consistently.
Forgetting the prior term. If empirical priors are intended, add $\log(n_k/N)$ to each class score. Unequal counts can affect the assigned class, as the second example shows.
Using an incomplete matrix product. Compute $x'\hat{\Sigma}_p^{-1}\hat{\mu}_k$ with the inverse covariance included. A numerical error in this term can change the score even if the final class happens to remain the same.
Ignoring singularity. A pooled covariance estimate with no inverse does not support this plug-in score as written.

Key takeaway: Estimate each class mean from its labeled observations, pool deviations from those class means to estimate the shared covariance, choose priors appropriate to the classification setting, and assign a new observation to the class with the largest plug-in LDA score.

Check Your Understanding

For each question, identify the quantity or calculation that determines the answer.

With $K$ classes and $N$ total labeled observations, why does the pooled within-class covariance estimate use $N-K$ in its denominator?
What variation is incorrectly included if observations are centered at the overall mean rather than at their class means?
For a one-dimensional problem with $\hat{\Sigma}_p=2$, $\hat{\mu}_k=3$, and $\hat{\pi}_k=1/4$, write the plug-in score at $x=2$.
In the second worked example, why would equal priors assign $x=3$ differently from empirical priors?
Give a condition that makes $\hat{\Sigma}_p$ necessarily singular, and explain which part of the LDA rule then cannot be computed as written.