Why High-Dimensional Precision Estimation Needs Regularization

Why Inverting the Sample Covariance Can Fail

A precision matrix is the inverse covariance matrix. In a Gaussian model, its entries describe conditional relationships among variables, which makes it useful for studying dependence and graphical structure. But the usual unregularized estimate requires inverting the sample covariance matrix. When the dimension is large relative to the sample size, that matrix may be singular; even when it is invertible, it may have very small eigenvalues that make inversion highly sensitive to data variation.

The first issue is algebraic: the observations constrain the sample covariance to a subspace whose dimension is limited by the number of observations. The second is numerical and statistical: a small eigenvalue becomes a large eigenvalue of the inverse. These are related problems, but they are not the same. Rank gives an exact condition for singularity; parameter counting and small eigenvalues help explain why instability can arise more broadly.

Definition: For observations $X_1,\ldots,X_n$ in $p$ dimensions, let $\bar X=n^{-1}\sum_{i=1}^n X_i$. The centered sample covariance is $$ S=\frac{1}{n}\sum_{i=1}^n (X_i-\bar X)(X_i-\bar X)^\mathsf{T}. $$ The precision matrix is $\Omega=\Sigma^{-1}$, where the population covariance $\Sigma$ is positive definite. The unregularized Gaussian precision estimate, when it exists, is $S^{-1}$.

The divisor $n$ is convenient for the Gaussian likelihood; using $n-1$ instead rescales the sample covariance but does not change its rank or whether it is invertible. When the mean is unknown, centering the observations imposes a linear constraint: the centered vectors sum to zero. Consequently, there can be at most $n-1$ linearly independent centered vectors.

Theorem: The centered sample covariance satisfies $$ \operatorname{rank}(S)\leq \min(p,n-1). $$ In particular, if $p\geq n$, then $S$ is singular.

Proof. Form the $p$-by-$n$ matrix $Y$ whose $i$th column is $X_i-\bar X$. Then $S=n^{-1}YY^\mathsf{T}$, so $\operatorname{rank}(S)\leq\operatorname{rank}(Y)\leq p$. Also, $Y\mathbf{1}=0$, where $\mathbf{1}$ is the $n$-vector of ones. Its columns therefore satisfy a nontrivial linear relation, which gives $\operatorname{rank}(Y)\leq n-1$. Combining these inequalities proves the result. $\square$

What Singularity Does to the Likelihood

For Gaussian observations with unknown mean, profiling the mean at $\bar X$ gives a log-likelihood for a positive-definite precision matrix $\Omega$ of the form $$ \ell(\Omega)=C+\frac{n}{2}\{\log\det(\Omega)-\operatorname{tr}(S\Omega)\}, $$ where $C$ does not depend on $\Omega$. Thus, maximizing the likelihood requires balancing the log determinant against the trace term. If $S$ is positive definite, the balance has a finite maximizer. If $S$ is singular, precision can grow without limit in a direction that the observed scatter does not constrain.

Theorem: If $S$ is singular and positive semidefinite, the Gaussian precision log-likelihood above is unbounded above over positive-definite $\Omega$. Hence there is no finite unregularized maximum-likelihood estimate in that parameter space.

Proof. Since $S$ is singular, there is a unit vector $v$ such that $Sv=0$. For $t\geq1$, define $\Omega_t=I+(t-1)vv^\mathsf{T}$. Its eigenvalue in direction $v$ is $t$, and its eigenvalues in directions orthogonal to $v$ are $1$, so $\Omega_t$ is positive definite and $\log\det(\Omega_t)=\log t$. Moreover, $$ \operatorname{tr}(S\Omega_t) =\operatorname{tr}(S)+(t-1)\operatorname{tr}(Svv^\mathsf{T}) =\operatorname{tr}(S)+(t-1)v^\mathsf{T}Sv =\operatorname{tr}(S). $$ Therefore $\ell(\Omega_t)=C+\frac{n}{2}\{\log t-\operatorname{tr}(S)\}$, which tends to infinity as $t$ tends to infinity. No finite positive-definite matrix can maximize this likelihood. $\square$

For comparison, if $S$ is positive definite, the maximizer is $S^{-1}$. To see this, set $A=S^{1/2}\Omega S^{1/2}$, which is positive definite. Then $$ \log\det(\Omega)-\operatorname{tr}(S\Omega) =-\log\det(S)+\log\det(A)-\operatorname{tr}(A). $$ If the eigenvalues of $A$ are $\lambda_1,\ldots,\lambda_p>0$, the part depending on $A$ is $\sum_{j=1}^p(\log\lambda_j-\lambda_j)$. Each term is maximized uniquely at $\lambda_j=1$, since $\log\lambda-\lambda\leq-1$ for every $\lambda>0$. The maximum occurs at $A=I$, or $\Omega=S^{-1}$. This also shows exactly why the ordinary estimate cannot be used when $S$ is singular.

Worked Example: A Singular Covariance from Three Observations

Consider $n=3$ observations in $p=4$ dimensions: $X_1=(1,0,0,0)^\mathsf{T}$, $X_2=(0,1,0,0)^\mathsf{T}$, and $X_3=(-1,-1,0,0)^\mathsf{T}$. Their sample mean is zero. The centered covariance is $$ S=\frac{1}{3}\sum_{i=1}^3X_iX_i^\mathsf{T} =\begin{pmatrix} 2/3&1/3&0&0\\ 1/3&2/3&0&0\\ 0&0&0&0\\ 0&0&0&0 \end{pmatrix}. $$ Its rank is $2$, so it is singular. In particular, $v=(0,0,1,0)^\mathsf{T}$ is a unit vector with $Sv=0$.

There are $p(p+1)/2=4(5)/2=10$ distinct entries in a symmetric $4$-by-$4$ precision matrix. For $\Omega_t=I+(t-1)vv^\mathsf{T}$, the determinant is $t$, while $\operatorname{tr}(S\Omega_t)=\operatorname{tr}(S)=4/3$. The log-likelihood is therefore $$ \ell(\Omega_t)=C+\frac{3}{2}\left(\log t-\frac{4}{3}\right). $$ It increases without bound as $t$ increases. The problem is not merely that matrix inversion is inconvenient: under the unregularized Gaussian likelihood, no finite estimate exists.

Rank Is Exact; Parameter Count Is a Warning

A symmetric $p$-by-$p$ matrix has $p(p+1)/2$ distinct entries, so an unrestricted precision matrix has that many free parameters. This count is a useful warning: the number of covariance relationships to estimate grows quadratically with dimension. But comparing this count directly with $n$ is not an exact test for invertibility. Invertibility of $S$ is governed by its rank, not by whether the parameter count is larger than the number of observations.

For centered observations, $p\geq n$ forces singularity. If $p<n$, full rank is possible, but it is not guaranteed: the centered data must actually span all $p$ dimensions. Conversely, a parameter count greater than $n$ does not itself force singularity. The following example makes that distinction concrete.

Worked Example: More Covariance Entries Than Observations, Yet Full Rank

Let $p=10$ and $n=20$. An unrestricted symmetric precision matrix has $10(11)/2=55$ distinct entries, more than the 20 observations. Nevertheless, take the observations to be $X_{2j-1}=\sqrt{10}\,e_j$ and $X_{2j}=-\sqrt{10}\,e_j$, for $j=1,\ldots,10$, where $e_j$ is the $j$th coordinate vector in ten dimensions.

The paired observations sum to zero, so $\bar X=0$. For each $j$, the two outer products contribute $10e_je_j^\mathsf{T}+10e_je_j^\mathsf{T}=20e_je_j^\mathsf{T}$. Thus $$ S=\frac{1}{20}\sum_{j=1}^{10}20e_je_j^\mathsf{T}=I_{10}. $$ This covariance is positive definite, and the unregularized estimate is $S^{-1}=I_{10}$. The example does not say that estimating 55 parameters from 20 observations is generally reliable. It shows that parameter count alone is not the algebraic criterion for whether the sample covariance can be inverted.

Why an Invertible Estimate Can Still Be Unstable

Full rank only rules out exact singularity. If the smallest eigenvalue of $S$ is close to zero, the corresponding eigenvalue of $S^{-1}$ is large. Small changes in the data can then cause large changes in the estimated precision. A useful way to make this precise is to compare $S^{-1}$ with the inverse of a perturbed covariance matrix.

Theorem: Let $S$ be positive definite, and let $E$ be symmetric with $\|E\|_2<\lambda_{\min}(S)$, where $\|\cdot\|_2$ is the spectral norm and $\lambda_{\min}(S)$ is the smallest eigenvalue of $S$. Then $S+E$ is positive definite and $$ \|(S+E)^{-1}-S^{-1}\|_2 \leq \frac{\|E\|_2}{\lambda_{\min}(S)\{\lambda_{\min}(S)-\|E\|_2\}}. $$

Proof. For every unit vector $u$, $u^\mathsf{T}(S+E)u\geq\lambda_{\min}(S)-\|E\|_2>0$, so $S+E$ is positive definite and $\|(S+E)^{-1}\|_2\leq\{\lambda_{\min}(S)-\|E\|_2\}^{-1}$. Also, the inverse identity $$ (S+E)^{-1}-S^{-1}=-(S+E)^{-1}ES^{-1} $$ follows by multiplying out the right-hand side. Taking spectral norms and using $\|S^{-1}\|_2=1/\lambda_{\min}(S)$ gives the stated bound. $\square$

Worked Example: A Small Eigenvalue Amplifies a Small Change

Suppose a positive-definite sample covariance has eigenvalues $1$, $0.25$, and $0.0004$. In the same eigenvector coordinates, its inverse has eigenvalues $1$, $4$, and $2500$. Now change only the smallest covariance eigenvalue from $0.0004$ to $0.0001$. The covariance perturbation has spectral norm $0.0003$, but the corresponding inverse eigenvalue changes from $2500$ to $10000$, a difference of $7500$.

The covariance remains invertible in both cases. Yet a change of $0.0003$ along one direction produces a very large change in the estimated precision along that direction. This is why checking only whether a numerical inverse exists is not enough; the eigenvalues, or equivalently the conditioning, matter as well.

What Regularization Addresses

Regularization changes the estimation problem so that directions unsupported by the data do not receive arbitrarily large precision, or so that poorly determined directions are stabilized. Depending on the model and scientific goal, this can mean imposing structure on the precision matrix or modifying the likelihood objective. The specific choice has consequences: it can introduce bias, encode assumptions, and affect which conditional relationships are estimated. The key motivation here is more basic: when the unrestricted likelihood has no finite maximizer, or inversion is extremely sensitive, an unmodified inverse sample covariance cannot provide a dependable estimate.

Common Mistakes

Assuming $n>p$ guarantees invertibility. It makes full rank possible for centered data, but the observed vectors still must span all $p$ dimensions.
Treating parameter count as the rank condition. The count $p(p+1)/2$ describes the number of unrestricted matrix entries. The exact singularity condition comes from the rank of the centered data.
Calling a singular-likelihood problem a numerical failure. When $S$ is singular, the unregularized Gaussian likelihood is genuinely unbounded over positive-definite precision matrices; a more accurate matrix inversion cannot fix that.
Assuming every invertible covariance gives a stable precision estimate. A very small positive eigenvalue is still inverted into a very large one, amplifying perturbations.

Key takeaway: Centering limits the rank of the sample covariance to at most $\min(p,n-1)$. If it is singular, the unregularized Gaussian precision likelihood has no finite maximizer; if it is nearly singular, its inverse may be unstable. Parameter count signals complexity, but rank and eigenvalue size diagnose the immediate problems.

Check Your Understanding

Use the rank and likelihood arguments to answer these questions.

For $n=12$ centered observations in $p=12$ dimensions, what is the largest possible rank of the sample covariance? What does that imply for its invertibility?
Why does a null vector of $S$ allow the Gaussian precision log-likelihood to increase without bound?
In what sense is comparing $p(p+1)/2$ with $n$ informative, and why is it not an exact rank test?
If the smallest eigenvalue of $S$ is positive but very small, what happens to the corresponding eigenvalue of $S^{-1}$?
For centered data with $p<n$, is full rank guaranteed? State the additional condition needed.