Robust Latent Variable Graphical Models: Exact Recovery and Fundamental Limits

Chuansen Peng; Xiaojing Shen

stat.MLOct 19, 2026

Robust Latent Variable Graphical Models: Exact Recovery and Fundamental Limits

Authors: Chuansen Peng, Xiaojing Shen

Organizations: School of Mathematics, Sichuan University, Chengdu 610064, China

Abstract

We study robust estimation of a Gaussian precision matrix that decomposes as a sparse matrix minus a symmetric, possibly indefinite, low-rank matrix, a structure that arises when latent-variable marginalization is combined with conditioning on linear measurements. The proposed convex estimator replaces the sample covariance by the centered Gaussian measure-transformed (MT) covariance, which preserves the signed graph of the sparse component under the nominal model, and penalizes the low-rank component by the nuclear norm of a congruence transform, an atomic norm that encodes structural side information. Our main result is a primal--dual witness theorem: under restricted invertibility and irrepresentability conditions on the Fisher information, the raw convex optimizer is unique, has no false edges, recovers the signs of all sufficiently strong edges, and has exactly the true rank and inertia. These conditions reduce to the graphical-lasso conditions when the low-rank part is absent, and we verify them from primitive incoherence and near-isotropy assumptions. Over a class on which they hold uniformly, we prove information-theoretic limits for Frobenius estimation, signed graph recovery, and inertia recovery. In particular, inertia recovery requires $N\gtrsim p/ν_{\min}^2$ samples, where $ν_{\min}$ is the smallest nonzero eigenvalue modulus of the low-rank part; the spectral beta-min condition of the witness theorem is necessary up to a $\sqrt{\log p}$ factor; and the Frobenius risk attains the minimax rate up to a factor $\sqrt{μ(1+\min\{d,r\})\log p}$ , where $μ$ is the incoherence parameter, $d$ the maximal degree, and $r$ the rank. We further quantify robustness under Huber contamination, with a bias bound that is sharp for the MT functional and vanishes for remote outliers, and give a globally convergent linearized ADMM with closed-form updates.

Figures & tables

Task	Converse (every estimator)	Guarantee (Theorem 1 )	Ratio
Frobenius risk	$\sqrt{\frac{pd}{N}}+\sqrt{\frac{pr}{N}}$ (Thm. 3 )	$\sqrt{\mu(1+\min\{d,r\})\,\ell}\,\Big(\sqrt{\frac{pd}{N}}+\sqrt{\frac{pr}{N}}\Big)$ (Cor. 2 )	$\sqrt{\mu(1+\min\{d,r\})\log p}$
Edge magnitude $\beta_{\min}$	$\sqrt{\frac{\log p}{N}}$ (Thm. 4 )	$\sqrt{\frac{\mu r\,\ell}{N}}$ (C4)	$\sqrt{\mu r}$
Eigenvalue modulus $\nu_{\min}$	$\sqrt{\frac{p}{N}}$ (Thm. 5 )	$\sqrt{\frac{p\ell}{N}}$ (C1); $\sqrt{\frac{p\ell}{N}}\sqrt{\frac{p}{\mu r}}$ (C3)	$\sqrt{\log p}$ (C1); $\sqrt{\frac{p\log p}{\mu r}}$ (C3)
Sample size for inertia	$\frac{p}{\nu_{\min}^{2}}$ (Thm. 5 )	$\frac{p^{2}\ell}{\mu r\,\nu_{\min}^{2}}$ ( 24 )	$\frac{p\log p}{\mu r}$

TABLE I: Information-theoretic limits versus the witness guarantees over the common class

\mathcal{F}

(

\mathbf{B}=\mathbf{I}

\ell\asymp\log p

\nu_{\min}<\omega_{-}

; absolute constants and the constants

\bar{\kappa},\alpha

of Lemma 1 suppressed).

Fig. 1: Relative error

\|\hat{\Omega}-\Omega^{*}\|_{F}/\|\Omega^{*}\|_{F}

versus the sample size

N

(

p=100

, ER graph,

r=3

with inertia

(2,1)

\tau=2

). (a) Nominal Gaussian data. (b)

10\%

remote outliers. The dotted line has slope

-1/2

. Shaded bands are

95\%

confidence intervals over

50

trials.

Fig. 2: Relative error versus the contamination ratio

\epsilon

(

p=100

N=800

). (a) Remote outliers. (b) Intermediate-radius outliers at the least favourable radius of Proposition 8 .

	$F_{1}$ , $\epsilon=0.1$ (remote)				$F_{1}$ , $\epsilon=0$	ER, $\epsilon=0.1$
Method	Chain	ER	Scale-free	Lattice	ER, $t_{3}$	MCC	Rank (%)	Inertia (%)
GLasso	0.28 $\pm$ 0.04	0.31 $\pm$ 0.02	0.26 $\pm$ 0.01	0.29 $\pm$ 0.03	0.41 $\pm$ 0.04	0.27	—	—
LVGLasso	0.27 $\pm$ 0.04	0.29 $\pm$ 0.01	0.25 $\pm$ 0.03	0.28 $\pm$ 0.05	0.47 $\pm$ 0.05	0.25	2.0	0.0
tlasso	0.63 $\pm$ 0.04	0.58 $\pm$ 0.04	0.49 $\pm$ 0.04	0.55 $\pm$ 0.01	0.64 $\pm$ 0.03	0.57	—	—
Trimmed GLasso	0.61 $\pm$ 0.05	0.56 $\pm$ 0.07	0.47 $\pm$ 0.05	0.53 $\pm$ 0.04	0.58 $\pm$ 0.02	0.55	—	—
LV-SKEPTIC	0.79 $\pm$ 0.03	0.71 $\pm$ 0.03	0.60 $\pm$ 0.08	0.69 $\pm$ 0.02	0.76 $\pm$ 0.03	0.70	64.0	58.0

TABLE II: Edge-recovery

F_{1}

score (mean

\pm

standard deviation over

50

trials), MCC, and percentages of exact rank and inertia recovery of

\hat{\mathbf{L}}

(

p=100

N=400

). “—” indicates a method without a low-rank component. Best values in bold.

Fig. 3: Precision–recall curves for edge recovery obtained by sweeping

\lambda

(

p=100

N=400

\epsilon=0.1

remote). (a) ER graph. (b) Scale-free graph. Bracketed values in the legend are the AUPRC.

Fig. 4: Dependence on the dimension

p

(

N=400

\epsilon=0.1

remote). (a) Relative estimation error. (b) Edge-recovery

F_{1}

score.

Fig. 5: Empirical probability that

\hat{\mathbf{L}}

has the rank and inertia of

\mathbf{L}^{*}

versus

N/\log p

(

\epsilon=0.05

remote,

200

trials per point). (a) GR-LVGM (solid) and MT-LVGM (dashed) for three dimensions. (b) Comparison of methods at

p=100

Variant	$\Omega$ err.	$\mathbf{L}$ err.	$F_{1}$	Inertia (%)	$\Omega$ err.	$\mathbf{L}$ err.	$F_{1}$	Inertia (%)
	$\epsilon=0$				$\epsilon=0.1$ (remote)
GR-LVGM (full)	0.225	0.29	0.86	96	0.250	0.33	0.83	90
w/o MT (sample covariance, $\tau>0$ )	0.215	0.27	0.87	97	0.910	0.98	0.29	0
w/o structural prior ( $\tau=0$ )	0.255	0.36	0.82	86	0.280	0.40	0.78	77
PSD-constrained $\mathbf{L}\succeq\mathbf{0}$ (trace penalty)	0.272	0.58	0.79	0	0.300	0.61	0.76	0
w/o low-rank part ( $\gamma\to\infty$ )	0.385	1.00	0.63	—	0.405	1.00	0.60	—

TABLE III: Ablation study (

p=100

N=400

, ER graph). Relative errors of

\hat{\Omega}

and

\hat{\mathbf{L}}

, edge-recovery

F_{1}

score, and percentage of trials with exact inertia recovery. For the misspecified operators,

\tau

is fixed at the value tuned for the correct

\mathbf{K}

Fig. 6: Sensitivity of GR-LVGM (

p=100

N=400

). (a) Prior strength

\tau

for the correct, partially rewired, and random structural operators; the dashed line marks

\tau=0

. (b) MT scale under nominal, remote, and intermediate-radius contamination (

\epsilon=0.1

). (c) Perturbation of the congruence operator: relative error (left axis) and

F_{1}

score (right axis), with the

\tau=0

reference levels.

Fig. 7: Point-mass contamination

\Psi=\delta_{\mathbf{r}}

. (a) Population bias of the MT covariance: exact values from ( 34 ) (solid), the bound ( 33 ) (dashed), and Monte Carlo estimates (circles); the dotted line marks

\|\mathbf{r}\|_{2}^{2}=2s^{2}

. (b) Mean estimated rank (

\epsilon=0.05

N=800

) for

\tau=0

and

\tau=2

, with a generic and a smooth spurious direction.

Fig. 8: Computational behaviour. (a) Relative objective suboptimality of Algorithm 1 for several

\rho

(

p=200

), and of the unlinearized four-block splitting. (b) Wall-clock time to reach a relative tolerance of

10^{-6}

versus

p

; the dotted line has slope

3

Fig. 9: Single-trial estimates ( $p=50$ , $N=400$ , $\epsilon=0.1$ remote). Top: off-diagonal entries of $\mathbf{S}^{*}$ and of the estimated sparse components, with $F_{1}$ scores and numbers of false positives (FP) and false negatives (FN). Bottom: $\mathbf{L}^{*}$ and the estimated low-rank components, with their inertia and relative Frobenius errors.

Site	TR (s)	$T$	ADHD	TDC
NYU	2.00	172	98	86
Peking	2.00	236	52	88
KKI	2.50	124	20	58
NeuroIMAGE	1.96	261	27	22
OHSU	2.50	74 ( $\times$ 3 )	18	51
Total	—	—	215	305

TABLE IV: ADHD-200 sites used in this study (after quality control).

T

is the number of acquired frames per run.

Fig. 10: Held-out NLL per region on scrubbed test frames versus the number of training frames

N

(ADHD-200, Peking and NeuroIMAGE). (a) Scrubbed training frames. (b) Unscrubbed training frames. Shaded bands are

95\%

confidence intervals over scans.

Fig. 11: Sensitivity to head motion. (a) Correlation between the edge weights estimated from unscrubbed and scrubbed data versus the quintile of mean FD (median FD in parentheses). (b) Percentage of edges whose weight is significantly associated with mean FD across subjects (

p<0.05

); the dashed line is the null level. Error bars are

95\%

bootstrap intervals.

Fig. 12: Split-half reproducibility on scrubbed data. (a) Dice coefficient of the edge sets (violins: distribution over scans; boxes: interquartile range). (b) Median edgewise intraclass correlation. (c) Percentage of scans for which

\hat{\mathbf{L}}

has the same rank, or the same inertia, in both halves.

Fig. 13: Effect of global signal regression (GSR) on

\hat{\mathbf{L}}

. (a) Percentage of scans with at least one negative eigenvalue. (b) Distribution of the estimated rank of GR-LVGM . (c) Alignment of the eigenvector of the negative eigenvalue with the constant vector under GSR; the dashed line is the expected alignment of a random unit vector.

Fig. 14: ADHD versus TDC classification from the estimated networks. (a) Pooled sites,

10\times 10

-fold stratified cross-validation (mean and standard deviation over repetitions). (b) Leave-one-site-out (mean and standard deviation over held-out sites; dots are individual sites). Horizontal lines mark the full Pearson correlation matrix and chance level.

Fig. 15: Sensitivity of GR-LVGM on unscrubbed data (

N=140

). (a) Held-out NLL and (b) split-half Dice versus

\tau

for three structural operators; the dashed line marks

\tau=0

. (c) Held-out NLL versus the MT scale for scrubbed and unscrubbed training data; dotted lines mark the sample-covariance limit

s\to\infty

, and the vertical line the default scale.

Fig. 16: Agreement with established features of rs-fMRI networks at a common edge density of

5\%

. (a) Fraction of homotopic pairs recovered. (b) Modularity with respect to the canonical networks, subcortex, and cerebellum. (c) Fraction of non-homotopic edges longer than

75

mm. Error bars are

95\%

confidence intervals over scans.

Fig. 17: ADHD versus TDC edgewise comparison. (a) Number of edges with

p<0.005

, with and without mean FD as a covariate; the dashed line is the number expected by chance. (b) Jaccard index between the two edge sets. (c) Share of the edges (FD covariate included) that connect only sensorimotor and visual regions, that involve the default mode, frontoparietal, subcortical, or cerebellar systems, or neither.

Fig. 18: Consensus TDC networks (edges selected in at least

50\%

of the TDC scans). Regions are arranged by hemisphere (left: L, right: R) and lobe (outer ring). Red and blue edges denote positive and negative partial correlations. Below each panel: number of edges, percentage of homotopic pairs recovered, and percentage of long-range edges between non-adjacent lobes.

Fig. 19: System-level summary of the ADHD versus TDC differences (

p<0.005

, FD covariate included). Chord width is proportional to the number of differential edges between two systems; loops denote within-system edges. DMN: default mode; FPN: frontoparietal; DAN: dorsal attention; VAN: ventral attention; LIM: limbic; SUB: subcortical; CB: cerebellum; VIS: visual; SMN: somatomotor.

Robust Estimation

Robust estimation infers reliable parameters and states despite outliers, heavy-tailed noise, data contamination, or model uncertainty.

See all 83 papers

Sep 15, 2026

Near-Optimal Nonconvex Matrix Completion

We study nonconvex methods for matrix completion, the problem of recovering a low-rank matrix from a subset of its entries. Convex methods achieve sample complexity linear in the matrix dimension and the rank, up to logarithmic factors, whereas global guarantees for commonly used nonconvex methods require a higher polynomial dependence on the rank. We close this gap by analyzing Riemannian gradient descent (RGD) and Riemannian Gauss--Newton (RGN) methods. For an

n\times n

matrix of rank

r

with incoherence parameter

μ

and condition number

κ

, the two methods achieve exact recovery with high probability from

O(μnr\log n\log(nκ))

and

O(μnr\log n\log(2μrκ))

observations, respectively. The methods use a multiscale residual initialization, while the analysis simultaneously controls the spectral error and incoherence. The resulting RGD iterates converge linearly, whereas RGN eventually converges Q-quadratically.

Jian-Feng Cai, Xiliang Lu, Juntao You

Department of Mathematics, Hong Kong University of Science and Technology. · School of Mathematics and Statistics, Hubei Center for Applied Mathematics, and Hubei Key Laboratory of Computational Science, Wuhan University, Wuhan 430072, China. · School of Artificial Intelligence, Hubei Center for Applied Mathematics, and Hubei Key Laboratory of Computational Science, Wuhan University, Wuhan 430072, China.

Riemannian Optimization Nonconvex Optimization Matrix Completion

Ask PDF

May 6, 2026

Convexity in Disguise: A Theoretical Framework for Nonconvex Low-Rank Matrix Estimation

Nonconvex methods have emerged as a dominant approach for low-rank matrix estimation, a problem that arises widely in machine learning and AI for learning and representing high-dimensional data. Existing analyses for these methods often require additional regularization to mitigate nonconvexity, even though such regularization is often unnecessary in practice. Moreover, most analyses rely on problem-specific arguments that are difficult to generalize to more complex settings. In this paper, we develop a theoretical framework for studying nonconvex procedures across a broad class of low-rank matrix estimation problems. Rather than focusing on a specific model, we reveal a fundamental mechanism that explains why nonconvex procedures can behave well in low-rank estimation. Our key device is a {\it benign regularizer} that does not alter the original update rule, but yields an equivalent locally strongly convex formulation of the algorithm. This perspective uncovers a disguised convexity inherent in the nonconvex procedure and provides a new route to theoretical guarantees for nonconvex low-rank matrix estimation.

Chengyu Cui, Gongjun Xu

Department of Statistics, University of Michigan

Nonconvex Optimization

Ask PDF

Jul 7, 2026

The Regularization Parameter: Sparse Precision Matrix Estimation

Sparse precision matrix estimation provides an interpretable and computationally efficient framework for modeling conditional dependencies in high-dimensional, low-sample-size data. A recurring challenge is appropriately selecting the regularization parameter that controls estimator sparsity and strikes a balance between underfitting and overfitting. We propose a closed-form, matrix-valued regularization parameter derived from the sampling distribution of the first-order optimality conditions of the

\ell_1

-regularized Gaussian maximum-likelihood estimator. By prescribing the probability that each nonzero entry of the estimator satisfies its optimality condition under resampling, we eliminate the need for cross-validation. The resulting regularization parameter is shown to attain asymptotic scaling properties that, under standard conditions, provide consistency and sparsistency of the estimator. On synthetic Gaussian and non-Gaussian datasets, as well as real-world gene microarray and neuroimaging applications, the proposed approach achieves estimation accuracy comparable to cross-validation, delivers superior support recovery, and reduces runtime by several orders of magnitude.

Aryan Eftekhari, Daniel Sergio Vega, Ernst-Jan Camiel Wit +1

Institute of Computing, Faculty of Informatics, Universit`a della Svizzera italiana (USI), Lugano, Switzerland

Graph Structure Learning Parameter Estimation Sparse Recovery Hyperparameter Optimization

Ask PDF

Abstract

Robust Estimation

Explore similar work

Near-Optimal Nonconvex Matrix Completion

Convexity in Disguise: A Theoretical Framework for Nonconvex Low-Rank Matrix Estimation

The Regularization Parameter: Sparse Precision Matrix Estimation