Tutorials › Statistical Theory III - Theories of Inference, Multivariate Methods, and Contingency Tables › Sample Covariance With Known and Estimated Means

Multivariate Normal and Wishart Theory · Tutorial 691 of 972

Sample Covariance With Known and Estimated Means

Learn why knowing the mean gives \(n\) degrees of freedom while estimating it leaves \(n-1\), and how that difference affects covariance estimates and likelihood.

Intermediate 12 min read

What You'll Learn

  • Define the unscaled scatter matrices for known and estimated population means.
  • Derive their Wishart distributions using an orthogonal transformation of a Gaussian sample.
  • Identify the degrees of freedom and covariance normalization in each setting.
  • Distinguish unbiased covariance estimates from maximum-likelihood estimates.
  • Determine when a covariance MLE exists and why singular scatter prevents one.

Why the Mean Changes the Scatter Distribution

A covariance matrix describes variation around a center. If the population mean is known, each observation can be centered at that fixed value. If the mean must be estimated from the same sample, the observations are instead centered at their sample mean, which uses one linear degree of freedom. For a Gaussian sample, this difference changes the exact distribution of the scatter matrix: the known-mean version has \(n\) degrees of freedom, while the estimated-mean version has \(n-1\).

Throughout, let \(X_1,\ldots,X_n\) be independent \(p\)-dimensional random vectors with common distribution \(N_p(\mu,\Sigma)\), where \(\Sigma\) is positive definite. Write \(\bar X=n^{-1}\sum_{i=1}^n X_i\). The scatter matrices below are unscaled: they are sums of outer products, not covariance estimates.

Definition: When \(\mu\) is known, the known-mean scatter is $$ A_\mu=\sum_{i=1}^n(X_i-\mu)(X_i-\mu)^T. $$ When \(\mu\) is estimated by \(\bar X\), the centered scatter is $$ C=\sum_{i=1}^n(X_i-\bar X)(X_i-\bar X)^T. $$ For a positive integer \(\nu\), \(W_p(\nu,\Sigma)\) denotes the Wishart distribution of \(\sum_{j=1}^{\nu}Z_jZ_j^T\), where \(Z_1,\ldots,Z_\nu\) are independent \(N_p(0,\Sigma)\) vectors.

The Two Scatter Distributions

The key fact is that an orthogonal transformation of the observation index separates the sample mean direction from the \(n-1\) contrast directions. Because the sample is Gaussian, the transformed rows remain independent Gaussian vectors. The mean direction contributes to known-mean scatter, but it is removed when the data are centered at \(\bar X\).

Theorem: For \(n\ge2\), under the independent \(N_p(\mu,\Sigma)\) model, $$ A_\mu\sim W_p(n,\Sigma),\qquad C\sim W_p(n-1,\Sigma). $$ In particular, \(E[A_\mu]=n\Sigma\) and \(E[C]=(n-1)\Sigma\). Thus \(A_\mu/n\) and \(C/(n-1)\) are unbiased estimators of \(\Sigma\), when their respective means are used for centering. For \(n=1\), \(C=0\).

Proof. Form the \(n\times p\) matrix \(Y\) whose \(i\)th row is \(X_i^T\). Choose an orthogonal \(n\times n\) matrix \(Q\) whose first row is \(n^{-1/2}(1,\ldots,1)\). The remaining rows are orthonormal contrasts: each is orthogonal to the vector of ones. Let \(Z=QY\), and denote its rows by \(Z_1^T,\ldots,Z_n^T\).

The rows of \(Y\) are independent \(N_p(\mu,\Sigma)\). Orthogonal transformation of the observation index gives jointly Gaussian rows \(Z_j\); their cross-covariances are zero because the rows of \(Q\) are orthogonal. Hence they are independent. Their means are \(E[Z_1]=\sqrt n\,\mu\) and \(E[Z_j]=0\) for \(j\ge 2\), and each has covariance \(\Sigma\).

An orthogonal transformation preserves sums of squared row vectors. Therefore the known-mean scatter, expressed in the transformed coordinates, is $$ A_\mu=\sum_{j=1}^n(Z_j-E[Z_j])(Z_j-E[Z_j])^T. $$ Each centered row in this sum is an independent \(N_p(0,\Sigma)\) vector, so \(A_\mu\sim W_p(n,\Sigma)\).

For the centered scatter, subtracting \(\bar X\) removes the component in the direction of the first row of \(Q\). The contrast rows are unchanged, giving $$ C=\sum_{j=2}^n Z_jZ_j^T. $$ These are \(n-1\) independent \(N_p(0,\Sigma)\) rows, so \(C\sim W_p(n-1,\Sigma)\). Finally, for any \(Z_j\sim N_p(0,\Sigma)\), \(E[Z_jZ_j^T]=\Sigma\). Taking expectations of the sums gives the two stated expectations. \(\square\)

Covariance Normalizations

The Wishart degrees of freedom determine the unbiased divisors. With a known mean, the unbiased estimate is \(A_\mu/n\). With an estimated mean, the conventional sample covariance matrix is \(S=C/(n-1)\), which is unbiased. These formulas do not mean that both estimates have the same likelihood divisor.

If the mean is unknown and estimated jointly with \(\Sigma\), the Gaussian maximum-likelihood estimate of \(\Sigma\), when it exists, is \(C/n\), not \(C/(n-1)\). If the mean is known and fixed, the maximum-likelihood estimate is \(A_\mu/n\). Thus, in the known-mean model the unbiased and maximum-likelihood estimates coincide; in the unknown-mean model they differ by the factor \((n-1)/n\).

Worked Example: Matrix Covariance With Known Mean

Suppose \(n=5\), \(p=2\), the population mean is known, and $$ \Sigma=\begin{pmatrix}4&1\\1&2\end{pmatrix}. $$ The known-mean scatter has distribution \(A_\mu\sim W_2(5,\Sigma)\). Since \(E[A_\mu]=5\Sigma\), $$ E[A_\mu]=\begin{pmatrix}20&5\\5&10\end{pmatrix}. $$ Consequently \(\hat\Sigma_\mu=A_\mu/5\) is unbiased: \(E[\hat\Sigma_\mu]=\Sigma\). If \(A_\mu\) is positive definite, the known-mean Gaussian likelihood is maximized at \(A_\mu/5\), so this estimate is also the MLE. Here \(n=5\ge p=2\), and under the model \(A_\mu\) is positive definite with probability one.

Worked Example: The Univariate Degrees of Freedom

Let \(X_i\sim N(\mu,9)\) independently, with \(n=6\). If \(\mu\) is known, then $$ A_\mu=\sum_{i=1}^6(X_i-\mu)^2,\qquad \frac{A_\mu}{9}\sim\chi^2_6. $$ Thus \(E[A_\mu]=6(9)=54\), and \(A_\mu/6\) is unbiased for the variance \(9\).

If instead \(\mu\) is estimated by \(\bar X\), then $$ C=\sum_{i=1}^6(X_i-\bar X)^2,\qquad \frac{C}{9}\sim\chi^2_5. $$ Therefore \(E[C]=5(9)=45\), and \(C/5\) is unbiased for \(9\). The maximum-likelihood estimate in the unknown-mean model is \(C/6\), whose expectation is \(45/6=7.5\), not \(9\). The lost degree of freedom is exactly why the unbiased divisor is \(5\), while the likelihood divisor remains \(6\).

When Does the Covariance MLE Exist?

The formulas \(A_\mu/n\) and \(C/n\) are covariance MLEs only when the corresponding scatter matrix is positive definite. A covariance parameter must itself be positive definite. If a scatter matrix is singular, the likelihood can increase without bound by shrinking covariance along a direction in which the scatter has no variation. In that case there is no finite positive-definite covariance MLE.

For a Gaussian sample, \(A_\mu\) is positive definite with probability one when \(n\ge p\), and \(C\) is positive definite with probability one when \(n-1\ge p\). If these inequalities fail, the corresponding scatter is necessarily singular by its rank bound. To see the almost-sure full-rank assertion when the number of Gaussian rows is at least \(p\), arrange those rows as a matrix. Its entries have a joint density, and the determinant of a selected \(p\times p\) submatrix is a nonzero polynomial in those entries. The zero set of a nonzero polynomial has Lebesgue measure zero, so that submatrix, and hence the scatter, is nonsingular almost surely.

Theorem: If \(A_\mu\) is positive definite, the known-mean covariance MLE is \(A_\mu/n\). If \(C\) is positive definite, the unknown-mean covariance MLE is \(C/n\), with the mean MLE \(\bar X\). If the relevant scatter matrix is singular, the likelihood over positive-definite covariance matrices is unbounded and has no finite maximizer.

Derivation. For a fixed mean and scatter \(A\), the part of the Gaussian log-likelihood that depends on \(\Sigma\) is $$ \ell(\Sigma)=-\frac n2\log|\Sigma|-\frac12\operatorname{tr}(\Sigma^{-1}A)+\text{constant}. $$ When the mean is unknown, maximizing first over the mean gives \(\bar X\), and the same expression results with \(A=C\). Suppose the relevant scatter \(A\) is positive definite and put \(\Omega=\Sigma^{-1}\). Let \(B=A^{1/2}\Omega A^{1/2}\), which ranges over positive-definite matrices as \(\Omega\) does. Up to a constant independent of \(B\), the log-likelihood becomes $$ \frac n2\log|B|-\frac12\operatorname{tr}(B). $$ If \(\lambda_1,\ldots,\lambda_p\) are the positive eigenvalues of \(B\), this is a sum of terms \((n/2)\log\lambda_j-\lambda_j/2\). Each term is uniquely maximized at \(\lambda_j=n\), since its derivative is \(n/(2\lambda_j)-1/2\), positive for \(\lambda_j<n\) and negative for \(\lambda_j>n\). Thus \(B=nI_p\), so \(\Omega=nA^{-1}\) and \(\hat\Sigma=A/n\).

Now suppose \(A\) is singular with rank \(r<p\). Choose orthonormal coordinates in which \(A\) is diagonal with \(r\) positive diagonal entries and \(p-r\) zeros. Consider diagonal precision matrices \(\Omega\) whose entries in the \(r\) positive-scatter directions are fixed and positive, while each of the remaining \(p-r\) entries equals \(t\). The trace term \(\operatorname{tr}(\Omega A)\) is unchanged as \(t\) grows, but \(\log|\Omega|\) increases by \((p-r)\log t\). Since the log-likelihood in precision form contains \(+(n/2)\log|\Omega|\), it tends to \(+\infty\) as \(t\to\infty\). This proves the nonexistence claim. \(\square\)

Worked Example: A Singular Centered Scatter

Take \(p=3\), \(n=3\), and the observed vectors \(X_1=(0,0,0)^T\), \(X_2=(1,0,0)^T\), and \(X_3=(0,1,0)^T\). Their mean is \(\bar X=(1/3,1/3,0)^T\). Subtracting this mean and summing the outer products gives $$ C=\begin{pmatrix} 2/3&-1/3&0\\ -1/3&2/3&0\\ 0&0&0 \end{pmatrix}. $$ The third coordinate has zero scatter, so \(C\) is singular; also the rank of a centered scatter is at most \(n-1=2<p=3\). There is no finite covariance MLE over positive-definite matrices. In the likelihood, the precision in the third coordinate can increase without affecting \(\operatorname{tr}(\Sigma^{-1}C)\), while increasing the log-determinant term without bound. The formal expression \(C/n\) is singular and therefore is not an admissible positive-definite covariance estimate.

Common Pitfall and Practical Interpretation

A common mistake is to use \(n-1\) in every covariance calculation because the familiar sample covariance divides by \(n-1\). That divisor is appropriate for the unbiased covariance estimate when the mean is estimated. It is not the divisor in the known-mean scatter distribution, and it is not the Gaussian MLE divisor when the mean is estimated. Keep three questions separate: which center is used, how many independent Gaussian contrast vectors remain, and whether the goal is unbiasedness or likelihood maximization.

The same distinction appears directly in the distributions: $$ n\left(\frac{A_\mu}{n}\right)\sim W_p(n,\Sigma), \qquad (n-1)\left(\frac{C}{n-1}\right)\sim W_p(n-1,\Sigma). $$ The scatter distribution identifies the correct degrees of freedom; the estimator's divisor then determines its normalization. For likelihood work, also check positive definiteness before asserting an MLE exists.

Check Your Understanding

Use the definitions and derivations above to answer these questions.

  1. For \(n\) independent \(N_p(\mu,\Sigma)\) observations with known \(\mu\), state the distribution of \(A_\mu\) and its unbiased covariance normalization.
  2. Why does centering at \(\bar X\) reduce the Wishart degrees of freedom by one?
  3. In the unknown-mean Gaussian model, distinguish the unbiased covariance estimate from the covariance MLE.
  4. For \(p=4\) and \(n=4\), can the centered scatter have a finite covariance MLE? Explain using its rank bound.
  5. What happens to the likelihood when a scatter matrix is singular and precision is increased along a null direction?
Clinical Biostats