Relaxed Sparsest-Permutation Formulation for Causal Discovery at Scale

Abstract: Despite the growing availability of large datasets, causal structure learning remains computationally prohibitive at scale. We revisit sparsest-permutation learning for linear structural equation models and show that exact Cholesky factorization is unnecessary for structure recovery. This observation motivates a support-level relaxation that searches for sparse triangular factors over a precision-support screening graph. The relaxed formulation can be efficiently evaluated via masked zero-fill incomplete Cholesky factorization, enabling scalable comparison of candidate orderings. At the population level, we establish soundness for Markov equivalence class (MEC) recovery under no-cancellation and sparsest Markov representation assumptions, as well as robustness to ordering misspecification. Motivated by these guarantees, we introduce SCOPE, a sparse-Cholesky pipeline that provides a scalable implementation of the relaxed formulation. Experiments on synthetic and real datasets demonstrate that SCOPE matches the MEC recovery accuracy of substantially slower baselines, while achieving significantly reduced runtime and scaling to 10k variables.

• The paper introduces the ReSP formulation as a relaxation of the sparsest-permutation constraint, using masked feasibility to operate on the support of the precision matrix.
• It leverages a sparse-Cholesky (IC(0)) based pipeline (SCOPE) to efficiently construct candidate factors and select the sparsest model, enhancing scalability.
• Empirical results demonstrate that SCOPE achieves competitive recovery accuracy and reduced runtime in both synthetic settings and high-dimensional breast cancer mRNA data.

Relaxed Sparsest-Permutation Formulation for Scalable Causal Discovery

Motivation and Context

Causal structure learning in high-dimensional settings faces a fundamental tradeoff between statistical precision and computational scalability. Exact methods, such as constraint-based (PC), score-based (GES), and continuous optimization approaches (NOTEARS, GOLEM), rapidly become infeasible as the number of variables ($p$) grows to tens of thousands. The sparsest-permutation (SP) formulation for structure learning in linear SEMs offers theoretical guarantees by searching over orderings to produce the sparsest Cholesky factor, but its computational demands are prohibitive. This paper introduces a principled relaxation—Relaxed Sparsest-Permutation (ReSP)—which circumvents the need for exact factorization and instead operates at the support level of the population precision matrix, enabling scalable causal discovery while preserving key recovery guarantees.

Relaxed Sparsest-Permutation (ReSP) Formulation

Traditional SP learning demands that for each variable ordering $\pi$, the permuted precision matrix $\Omega_\pi$ be represented exactly as $L_\pi L_\pi^\top$, where $L_\pi$ is lower triangular. ReSP replaces the exact constraint with masked feasibility: only entries at positions determined by a mask $M_\pi$ (the support of $\Omega_\pi$) need match in $L_\pi L_\pi^\top$ and $\Omega_\pi$. Zero-fill incomplete Cholesky factorization (IC(0)) is used for efficient construction of sparse triangular factors under the mask, leveraging well-established sparse linear algebra routines. The formulation ensures that no spurious fill-in occurs outside the screened support, and comparisons are made based on the sparsity of $L_\pi$ within $\pi$0. Theoretical analysis confirms that under no-cancellation and sparsest Markov representation (SMR) conditions, ReSP recovers the correct Markov equivalence class (MEC) and is robust to ordering misspecification.

SCOPE: Scalable Sparse-Cholesky Pipeline

SCOPE operationalizes the ReSP relaxation in finite sample settings. The pipeline consists of four stages:

1. Precision Estimation: A sparse estimator for the precision matrix (e.g., DRO-calibrated graphical lasso) provides the support mask.
2. Masked Factor Construction: IC(0) masked by the support constructs candidate $\pi$1 for each ordering in the candidate set.
3. Refit-and-Test Pruning: Local regressions prune small coefficients, eliminating false positives especially from moralization-induced edges.
4. Sparsity-Based Selection: The sparsest factor across candidates is selected as the output.

For maximal scalability, SCOPE defaults to a one-pass strategy with an approximate minimum-degree ordering (AMD), resulting in cubic worst-case complexity dominated by the initial precision estimation.

Theoretical Guarantees

The ReSP formulation yields the following:

• Soundness: With an admissible candidate ordering that aligns with the true structure and under SMR, ReSP minimizes a sparsity objective and recovers a DAG in the correct MEC.
• Robustness: Under arbitrary candidate orderings, the support mask guarantees that all induced edges remain confined to the moral graph, suppressing false discoveries outside the true skeleton.

Empirical Evaluation

Runtime-Accuracy Tradeoff

Figure 1: Runtime-accuracy tradeoffs for MEC recovery across baseline causal discovery methods (means ± error bars) for synthetic Gaussian linear SEMs with bounded in-degree $\pi$2, $\pi$3, and $\pi$4.

Empirical results demonstrate that SCOPE achieves strong recovery accuracy comparable to exact methods, but with drastically reduced runtime. Baselines such as PC, GES, GOLEM, and NOTEARS require extensive computational resources or fail to scale, while SCOPE readily handles $\pi$5 variables.

CPDAG Recovery Performance

Figure 2: CPDAG $\pi$6 versus runtime across problem sizes (bubble size $\pi$7); $\pi$8 denotes timed-out runs. Two one-pass SCOPE results: oracle and AMD heuristic orderings.

SCOPE maintains high $\pi$9 scores for CPDAG recovery across increasing $\Omega_\pi$0. When the true ordering is provided, the recovery approaches optimality ($\Omega_\pi$1), while the heuristic AMD ordering sustains moderate accuracy ($\Omega_\pi$2) even at $\Omega_\pi$3.

Skeleton and Structural Recovery

Figure 3: CPDAG-skeleton $\Omega_\pi$4 scores for bounded in-degree linear SEMs (means ± 1 standard error, 30 replicates).

Figure 4: Normalized structural Hamming distance (nSHD) for CPDAG recovery, indicating consistent structure estimation across replicates.

Across variable in-degree and noise types, SCOPE remains robust, with accuracy only moderately affected by graph densification. Competing methods such as FROSTY and RFD deteriorate rapidly in high-$\Omega_\pi$5 regimes, exhibiting large nSHD and near-zero $\Omega_\pi$6 beyond moderate $\Omega_\pi$7.

Real Data Analysis: Breast Cancer mRNA Expression

SCOPE was applied to TCGA mRNA expression data with $\Omega_\pi$8 and $\Omega_\pi$9, inferring transcription factor regulatory networks. Unlike PC and GES, which required partition-and-merge strategies due to scalability limits, SCOPE operated on the full dataset and achieved substantially higher rates of biologically validated forward edges and two-step paths. The forward-to-reverse ratio for TF → target edge identification was 1.77 for SCOPE, as compared to 1.13–1.32 for PC and GES, with a 2–5x increase in externally validated associations.

Practical and Theoretical Implications

The relaxation of the SP objective to the support level enables scalable MEC recovery, avoiding the combinatorial bottleneck of ordering search and the computational weaknesses of CI testing or dense likelihood optimization. The sparse-Cholesky pipeline is amenable to parallelization and distributed evaluation, further expanding its applicability to large-scale genomics, neuroscience, and socio-economic systems. The mask-based feasibility property offers robustness to misspecification and is a flexible foundation for extensions to interventions and latent variable models.

Future developments may focus on local ordering search in the ReSP context, parallelized evaluation pipelines, and adaptation to settings with measurement noise or incomplete data. The modularity of SCOPE facilitates integration with advanced sparse precision estimation and screening techniques, potentially enhancing statistical power for causal discovery at scale.

Conclusion

The Relaxed Sparsest-Permutation formulation and the SCOPE pipeline significantly improve scalability for causal structure learning in high-dimensional linear SEMs, retaining strong theoretical recovery guarantees and empirical accuracy. The shift from exact factorization to support-level screening and sparse triangular factors offers a practical, robust, and theoretically sound approach to MEC recovery on large graphs, broadening the scope of structure discovery in computational statistics and artificial intelligence (2605.05568).