Eigengap Sparsity for Covariance Parsimony

Technical Summary: Eigengap Sparsity for Covariance Parsimony

Problem Statement
Covariance estimation is a fundamental challenge in statistics, particularly in high-dimensional settings where the number of samples $n$ is insufficient to accurately estimate the $p(p+1)/2$ parameters of a $p$-dimensional covariance matrix. The sample covariance matrix (SCM), while the maximum likelihood estimate (MLE), suffers from significant sampling errors, notably the "eigenvalue-spreading" phenomenon described by the Marčenko–Pastur law, where sample eigenvalues disperse widely around their true values.

To address over-parameterization, parsimonious models are preferred. However, traditional model selection for covariance matrices involves discrete choices (e.g., selecting specific sparsity patterns or rank constraints), leading to computationally intractable combinatorial problems. Specifically, the authors focus on covariance parsimony induced by eigenvalue multiplicities, where grouping equal eigenvalues reduces the dimensionality of the parameter space by quotienting rotational invariance within eigenspaces. While this approach (Principal Subspace Analysis) offers a theoretical framework, exact optimization requires solving $2^{p-1}$ distinct unpenalized optimization problems, rendering it computationally prohibitive even for moderate dimensions.

Methodology
The paper proposes a continuous relaxation of the discrete parsimony problem, termed "eigengap sparsity," to make the optimization tractable.

1. Stratification and Penalty Derivation:
   The space of symmetric positive-definite matrices, $\text{Sym}_+(p)$, is stratified into submanifolds based on eigenvalue multiplicities. The dimension of a submanifold corresponding to a specific multiplicity pattern $\gamma$ is derived. The authors prove a key theorem (Theorem 1) stating that this dimension can be expressed as an $\ell_0$-norm of spectral gaps:
   $$\dim(\text{Sym}_+(\gamma)) = 1 + \|\delta(\lambda_s, \lambda_{s+1})\|_0 + \|\delta(\lambda_s, \lambda_t)\|_0$$
   where $\delta(a,b)=0$ iff $a=b$. This formulation reveals that parsimony is achieved not only by equalizing adjacent eigenvalues but also by equalizing non-adjacent ones, effectively "attracting" all eigenvalues toward each other.

2. $\ell_1$-Relaxation:
   To avoid the combinatorial complexity of the $\ell_0$-norm, the authors replace it with an $\ell_1$-norm, creating a convex relaxation. The resulting optimization problem minimizes a penalized likelihood:
   $$\min_{\lambda \in K_+(p)} \left( -2 \ln L(V \text{diag}(\lambda) V^\top) + \alpha \sum_{s=1}^{p-1} \left( \delta(\lambda_s, \lambda_{s+1}) + \sum_{t=s+1}^p \delta(\lambda_s, \lambda_t) \right) \right)$$
   where $K_+(p)$ is the monotone cone of ordered positive eigenvalues.

3. Optimization Algorithm:
   For the Gaussian case, the authors propose a Projected Gradient Descent algorithm on the monotone cone.

   • Gradient Step: Updates eigenvalues to maximize the likelihood while reducing the eigengap penalty.
   • Projection Step: Projects the updated eigenvalues back onto the monotone cone $K_+(p)$. This projection is equivalent to isotonic regression, solved efficiently via the Pools Adjacent Violators Algorithm (PAVA).
   • Mechanism: The algorithm iteratively performs a gradient step and then "equalizes" eigenvalues that violate ordering constraints. The authors note that while close eigenvalues are often a numerical curse in shrinkage literature, in this framework, they are a "blessing" as the projection step naturally equalizes them, inducing parsimony.
   • Eigengap Function: The authors utilize a relative eigengap function $\delta(a, b) = (a-b)/a$ (for $a \ge b$), arguing it yields better results than the absolute difference.

Key Contributions

• Theoretical Reformulation: The paper establishes a novel link between covariance parsimony and eigenvalue gaps, proving that the intrinsic dimension of a covariance manifold can be written as an $\ell_0$-norm of gaps between all pairs of eigenvalues.
• Relaxation Strategy: It introduces "eigengap sparsity" as an $\ell_1$-relaxation of this dimension penalty, transforming a discrete model selection problem into a continuous optimization problem.
• Algorithmic Link: The proposed algorithm draws a direct connection between covariance parsimony and shrinkage methods. It demonstrates that isotonic regression (historically used as a post-processing trick to handle numerical instability) naturally arises as the projection step in a penalized-likelihood framework, automatically grouping eigenvalues to reduce model complexity.
• Computational Efficiency: The method avoids the exponential complexity of exact Principal Subspace Analysis (PSA) by solving a single continuous optimization problem.

Experimental Results
The authors evaluated the method (ESCP) against the Sample Covariance Matrix (SCM), the Ledoit-Wolf (LW) shrinkage estimator, and exact Principal Subspace Analysis (PSA) on synthetic Gaussian datasets with varying dimensions ($p$) and sample sizes ($n$).

• Accuracy vs. Parsimony: ESCP achieved significantly lower penalized-likelihood scores ($L_p$) than SCM and LW in high-dimensional settings, indicating a better trade-off between fit and complexity.
• Estimation Error: Surprisingly, ESCP outperformed the LW estimator in terms of Frobenius norm estimation error ($L_F$), despite not being explicitly designed to minimize this metric.
• Eigenvalue Recovery: In experiments with isotropic and block-diagonal covariance structures, ESCP successfully recovered the true eigenvalue multiplicities (e.g., perfectly recovering $I_{100}$ and the block structure of $\Xi_{80,80,40}$), whereas SCM showed high dispersion and LW showed partial shrinkage without full grouping.
• Speed: ESCP was orders of magnitude faster than exact PSA, which became computationally prohibitive for $p \gtrsim 20$.

Significance and Claims
The paper claims that eigengap sparsity offers a promising perspective on covariance estimation by reframing it through the lens of parsimony induced by eigenvalue equalization. The authors highlight three main bridges to active research areas:

1. Sparse Optimization: The method induces sparsity via isotonic regression rather than classical thresholding (like LASSO).
2. Shrinkage Theory: It addresses the eigenvalue-spreading phenomenon by grouping eigenvalues, extending beyond low-rank models (like PPCA) that only group small eigenvalues as noise.
3. Geometric Statistics: The framework integrates both eigenvalues and the flag of eigenspaces, offering a more complete geometric view than methods focusing solely on eigenvalues.

The authors conclude that while the method shows promise in synthetic experiments and speed, it remains in an early stage of research, with future work needed to evaluate performance on diverse eigenvalue profiles, non-Gaussian distributions, and to explore theoretical guarantees and Riemannian metric extensions.