Tutorials › Multivariate Statistics (Graduate) › Gaussian Mean Estimation With Known Covariance

Gaussian Likelihood and Invariance · Tutorial 168 of 975

Gaussian Mean Estimation With Known Covariance

See why the sample mean is the unique Gaussian maximum-likelihood estimate with known covariance, and derive its exact sampling behavior.

Intermediate 10 min read

What You'll Learn

  • State the multivariate Gaussian model with known positive-definite covariance.
  • Derive the likelihood and log-likelihood for the unknown mean.
  • Prove that the sample mean is the unique maximum-likelihood estimate.
  • Calculate the covariance matrix of the sample mean.
  • Derive the sample mean’s exact multivariate normal distribution.
  • Find the sampling distribution of a linear contrast of the estimated mean.

Model and Estimation Target

Suppose \(X_1,\ldots,X_n\) are independent observations of a \(p\)-dimensional vector. The mean vector \(\mu\) is unknown, while the covariance matrix \(\Sigma\) is known. The goal is to estimate \(\mu\) by maximizing the likelihood, then describe how that estimate varies across repeated samples. Here, “known covariance” means that \(\Sigma\) is fixed and specified, not estimated from these observations.

Definition: The known-covariance Gaussian mean model assumes that \(X_1,\ldots,X_n\) are independent and identically distributed as \(N_p(\mu,\Sigma)\), where \(\mu\in\mathbb{R}^p\) is unknown and \(\Sigma\) is a known, symmetric positive-definite \(p\times p\) matrix. The sample mean is \(\bar X=n^{-1}\sum_{i=1}^n X_i\).

Positive definiteness ensures that \(\Sigma^{-1}\) exists and that the multivariate normal distribution has a density on \(\mathbb{R}^p\). For a single observation \(x\), that density is

$$ f(x;\mu,\Sigma) =(2\pi)^{-p/2}|\Sigma|^{-1/2} \exp\!\left\{-\frac12(x-\mu)^\mathsf{T}\Sigma^{-1}(x-\mu)\right\}. $$

Because the observations are independent, their joint likelihood is the product of these densities. As a function of \(\mu\), its log is a constant minus a sum of squared distances measured using the precision matrix \(\Sigma^{-1}\). The covariance is held fixed throughout this optimization.

Deriving and Verifying the Maximum

The likelihood and log-likelihood, with terms displayed to make their dependence clear, are

$$ L(\mu) =(2\pi)^{-np/2}|\Sigma|^{-n/2} \exp\!\left\{-\frac12\sum_{i=1}^n (X_i-\mu)^\mathsf{T}\Sigma^{-1}(X_i-\mu)\right\}, $$
$$ \ell(\mu) =-\frac{np}{2}\log(2\pi)-\frac n2\log|\Sigma| -\frac12\sum_{i=1}^n (X_i-\mu)^\mathsf{T}\Sigma^{-1}(X_i-\mu). $$

Differentiating with respect to the vector \(\mu\), using symmetry of \(\Sigma^{-1}\), gives the score

$$ \nabla_\mu\ell(\mu) =\Sigma^{-1}\sum_{i=1}^n(X_i-\mu) =n\Sigma^{-1}(\bar X-\mu). $$

Setting this score to zero yields \(\mu=\bar X\). A stationary point is not enough to establish a maximum, so we compare the log-likelihood at any candidate \(\mu\) directly with its value at \(\bar X\).

Theorem: Under the known-covariance Gaussian mean model, the unique maximum-likelihood estimator of \(\mu\) is \(\hat\mu=\bar X\).

Proof. For any fixed data vectors \(X_1,\ldots,X_n\), write \(X_i-\mu=(X_i-\bar X)+(\bar X-\mu)\) and expand the quadratic forms. The cross terms sum to zero because \(\sum_{i=1}^n(X_i-\bar X)=0\). Thus

$$ \sum_{i=1}^n(X_i-\mu)^\mathsf{T}\Sigma^{-1}(X_i-\mu) = \sum_{i=1}^n(X_i-\bar X)^\mathsf{T}\Sigma^{-1}(X_i-\bar X) +n(\bar X-\mu)^\mathsf{T}\Sigma^{-1}(\bar X-\mu). $$

The first term on the right does not depend on \(\mu\). The second is nonnegative because \(\Sigma^{-1}\) is positive definite, and it equals zero only when \(\mu=\bar X\). Since the log-likelihood is minus one half of this sum plus a constant, it is strictly smaller at every \(\mu\ne\bar X\). Therefore \(\bar X\) is the unique global maximizer. \(\square\)

The result does not require the coordinates to be uncorrelated or to have equal variances. In this unrestricted mean model, the estimate is the ordinary componentwise sample mean even when \(\Sigma\) is not diagonal. The known covariance still determines how sharply the likelihood falls away from its maximum.

Worked Example: Estimating a Two-Dimensional Mean

Suppose \(n=4\) independent observations have known covariance \[ \Sigma=\begin{pmatrix}4&1\\1&9\end{pmatrix}, \] and the observed vectors are \((2,1)^\mathsf{T}\), \((4,0)^\mathsf{T}\), \((3,3)^\mathsf{T}\), and \((5,2)^\mathsf{T}\). Find the maximum-likelihood estimate of the mean.

The matrix is positive definite: its first leading principal minor is \(4>0\), and its determinant is \(4(9)-1(1)=35>0\). The coordinate sums are \(2+4+3+5=14\) and \(1+0+3+2=6\). Hence

$$ \hat\mu=\bar X =\begin{pmatrix}14/4\\6/4\end{pmatrix} =\begin{pmatrix}3.5\\1.5\end{pmatrix}. $$

The estimate is the same componentwise average that maximizes the likelihood for any known positive-definite covariance. The off-diagonal entry of \(\Sigma\) does not alter this unrestricted estimate; it does affect the estimate’s sampling covariance, calculated below.

Covariance and Exact Sampling Distribution

The estimator is random because it depends on the sampled observations. Independence and the common covariance give

$$ E[\bar X]=\frac1n\sum_{i=1}^n E[X_i]=\mu, \qquad \operatorname{Cov}(\bar X) =\frac1{n^2}\sum_{i=1}^n\operatorname{Cov}(X_i) =\frac{\Sigma}{n}. $$

The covariance sum has no cross-covariance terms: for \(i\ne j\), independence implies \(\operatorname{Cov}(X_i,X_j)=0\). Thus \(\bar X\) is unbiased and its covariance matrix is exactly \(\Sigma/n\), not an estimate of that matrix.

Theorem: If \(X_1,\ldots,X_n\) are independent and identically distributed as \(N_p(\mu,\Sigma)\) for fixed positive-definite \(\Sigma\), then \[ \bar X\sim N_p\!\left(\mu,\frac{\Sigma}{n}\right). \] In particular, this is the exact finite-sample distribution, not an asymptotic approximation.

Proof. For \(t\in\mathbb{R}^p\), the moment-generating function of one observation is \(E[\exp(t^\mathsf{T}X_i)]=\exp(t^\mathsf{T}\mu+\tfrac12t^\mathsf{T}\Sigma t)\). Independence gives the moment-generating function of the average:

$$ E[\exp(t^\mathsf{T}\bar X)] =\prod_{i=1}^n E[\exp(t^\mathsf{T}X_i/n)] =\left[ \exp\!\left\{\frac{t^\mathsf{T}\mu}{n} +\frac{t^\mathsf{T}\Sigma t}{2n^2}\right\} \right]^n =\exp\!\left\{t^\mathsf{T}\mu+\frac12t^\mathsf{T}\frac{\Sigma}{n}t\right\}. $$

This is the moment-generating function of \(N_p(\mu,\Sigma/n)\). Since the moment-generating function exists in a neighborhood of the origin and determines the distribution there, \(\bar X\) has that multivariate normal distribution. \(\square\)

The diagonal entries of \(\Sigma/n\) are the variances of the estimated coordinates; the off-diagonal entries are their covariances. Increasing \(n\) reduces every entry by the factor \(1/n\), while the standard deviation of each coordinate is reduced by \(1/\sqrt n\). The correlation between two coordinates of \(\bar X\), when both variances are positive, is the same as their correlation under \(\Sigma\), since the common factor \(1/n\) cancels.

Worked Example: Reading the Estimator’s Sampling Covariance

Suppose \(n=10\) and the known covariance is \[ \Sigma=\begin{pmatrix}9&3\\3&4\end{pmatrix}. \] The determinant is \(36-9=27>0\), and the first diagonal entry is positive, so this symmetric matrix is positive definite. The exact covariance matrix of the sample mean is

$$ \operatorname{Cov}(\bar X)=\frac{\Sigma}{10} =\begin{pmatrix}0.9&0.3\\0.3&0.4\end{pmatrix}. $$

The first estimated coordinate has variance \(0.9\) and standard deviation \(\sqrt{0.9}\approx0.949\); the second has variance \(0.4\) and standard deviation \(\sqrt{0.4}\approx0.632\). Their correlation is \[ \frac{0.3}{\sqrt{0.9}\sqrt{0.4}} =\frac{0.3}{\sqrt{0.36}}=\frac{0.3}{0.6}=0.5. \] The sampling distribution is therefore \(N_2(\mu,\Sigma/10)\), with the mean vector \(\mu\) still unknown. Dividing the covariance matrix by \(n\) reduces variances and covariances, but it does not divide the standard deviations by \(n\).

Linear Summaries of the Estimated Mean

The full sampling distribution also gives the distribution of any fixed linear combination of the estimated coordinates. For a fixed vector \(a\in\mathbb{R}^p\), a linear transformation of a multivariate normal vector is normal. Its mean and variance follow from the mean and covariance of \(\bar X\):

$$ a^\mathsf{T}\bar X \sim N\!\left(a^\mathsf{T}\mu,\frac{a^\mathsf{T}\Sigma a}{n}\right). $$

This formula uses the same \(a\) for every sample and treats it as fixed, not as a quantity selected from the observed data. It is useful when the target is a contrast or a weighted combination rather than the entire mean vector.

Worked Example: Sampling Distribution of a Contrast

Let \(n=5\), let the known covariance be \[ \Sigma=\begin{pmatrix}4&2\\2&5\end{pmatrix}, \] and suppose the observed vectors are \((1,0)^\mathsf{T}\), \((2,1)^\mathsf{T}\), \((0,0)^\mathsf{T}\), \((1,1)^\mathsf{T}\), and \((2,0)^\mathsf{T}\). Consider the contrast of the two mean coordinates with \(a=(2,-1)^\mathsf{T}\).

The coordinate sums are \(6\) and \(2\), so \(\bar X=(1.2,0.4)^\mathsf{T}\). The maximum-likelihood estimate of \(2\mu_1-\mu_2\) is \[ a^\mathsf{T}\bar X=2(1.2)-0.4=2.0. \] To find its sampling variance, first multiply:

$$ \Sigma a =\begin{pmatrix}4&2\\2&5\end{pmatrix} \begin{pmatrix}2\\-1\end{pmatrix} =\begin{pmatrix}6\\-1\end{pmatrix}, \qquad a^\mathsf{T}\Sigma a =(2,-1)\begin{pmatrix}6\\-1\end{pmatrix}=13. $$

Consequently, \(2\bar X_1-\bar X_2\) has exact distribution \(N(2\mu_1-\mu_2,13/5)\), with variance \(2.6\). Its observed value is \(2.0\), but the sampling distribution is centered at the unknown contrast \(2\mu_1-\mu_2\), not necessarily at that observed value.

Common Pitfall and Interpretation

A common mistake is to reason that a covariance matrix with unequal variances or nonzero correlations must produce a weighted sample mean for the estimate of \(\mu\). That conclusion confuses estimating an unrestricted mean vector with estimating a mean subject to restrictions. In the present model, every coordinate of \(\mu\) can vary freely. The quadratic decomposition proves that the likelihood is maximized at the ordinary sample mean for every known positive-definite \(\Sigma\). The covariance changes the curvature of the likelihood and the estimator’s sampling uncertainty, not the location of its maximum.

Another mistake is to report \(\Sigma/n\) as the covariance of one observation. It is the covariance of \(\bar X\); one observation has covariance \(\Sigma\). Nor should the exact normal sampling distribution be described as an approximation that requires a large sample: normality here follows directly from the Gaussian model for every \(n\ge1\).

Key takeaway: With independent \(N_p(\mu,\Sigma)\) observations and known positive-definite \(\Sigma\), the unique maximum-likelihood estimate is \(\bar X\). It is unbiased, has covariance \(\Sigma/n\), and has exact distribution \(N_p(\mu,\Sigma/n)\).

Check Your Understanding

Use the model and results developed above to answer these questions.

  1. Why is positive definiteness of the known covariance matrix needed in the uniqueness argument for the maximum-likelihood estimate?
  2. If \(n=20\) and \(\Sigma_{jj}=16\), what are the variance and standard deviation of the \(j\)th coordinate of \(\bar X\)?
  3. Does a nonzero off-diagonal entry of \(\Sigma\) change the unrestricted maximum-likelihood estimate of \(\mu\)? Explain.
  4. For fixed \(a\in\mathbb{R}^p\), state the mean and variance of \(a^\mathsf{T}\bar X\).
  5. Why is the normal sampling distribution of \(\bar X\) exact even for a small sample?
Clinical Biostats