Statistical-Computational Trade-off in Variable Selection

• High-dimensional variable selection is a framework for recovering sparse signals when the number of predictors far exceeds the sample size.
• The trade-off reveals that optimal but NP-hard methods like best subset selection contrast with efficient approaches that require stronger design conditions.
• Recent advances, such as graphlet screening and empirical Bayes methods, balance statistical accuracy with computational feasibility to achieve near-optimal recovery.

The statistical-computational trade-off in high-dimensional variable selection concerns the interplay between the minimum statistical sample size, estimation precision, and the algorithmic feasibility of variable selection procedures as the number of predictors p grows much larger than the sample size n. Recent advances have systematically illuminated this interplay, showing that the design structure, signal sparsity, and the choice of selection algorithm fundamentally determine not just statistical optimality but also computational tractability.

1. Formulation of the High-Dimensional Variable Selection Problem

High-dimensional linear regression is generally modeled as $Y = X\beta + z$, where $Y \in \mathbb{R}^n$, $X \in \mathbb{R}^{n \times p}$ with $p \gg n$, and $z \sim N(0, I_n)$. The coefficient vector $\beta$ is assumed sparse, and variable selection aims to recover the support $S = \{j: \beta_j \ne 0\}$.

The minimax Hamming risk is used to quantify selection accuracy: $$\text{Hamm}_p^*(\theta, \kappa, r, a, \Omega) = \inf_{\hat{\beta}} \sup_{\mu \in \Theta_p^*(\tau_p, a)} \mathbb{E}\left[\sum_{j=1}^p 1\{\operatorname{sign}(\hat{\beta}_j) \ne \operatorname{sign}(\beta_j)\}\right],$$ where $\tau_p = \sqrt{2r\log(p)}$ denotes the minimal signal strength, and $\theta, r$ parameterize sparsity and signal regimes (Jin et al., 2012).

2. Statistically Optimal but Computationally Infeasible Procedures

Combinatorial variable selection procedures such as best subset selection (BSS) achieve the optimal information-theoretic sample complexity for support recovery: $$n \gtrsim \max\left\{ \frac{ \log(d-s) + \log(1/\delta)}{\theta^2/\sigma^2 },~ \log\binom{d-s}{s} + \log(1/\delta)\right\}$$ where $d$ is the dimension, $s$ the sparsity, $\theta$ the minimal signal, and $\sigma^2$ the noise variance (Gao et al., 5 Oct 2025, Roy et al., 2022). Subset selection minimizes the projection residual over all $s$-subsets: $$S^{\mathrm{BSS}} = \arg\min_{S \in \mathcal{S}_{d,s}} \| \Pi_S^\perp Y \|^2.$$ This approach is statistically minimax optimal in Hamming risk and support recovery, even under signal heterogeneity.

However, BSS is computationally intractable in general ($\mathsf{NP}$-hard), and no polynomial-time algorithm achieves the minimax sample complexity in the worst-case design, under standard complexity assumptions (Gao et al., 5 Oct 2025, Wang et al., 2014).

3. Efficient Algorithms and Provable Computational Barriers

Efficient algorithms, typically based on convex relaxation (e.g., Lasso or semidefinite programming for sparse PCA), require stronger conditions—usually on the design matrix—than is necessary statistically. For example, their success depends critically on restricted eigenvalues or incoherence. In the worst case, any polynomial-time estimator incurs an additional sample complexity cost inversely proportional to the square of the restricted eigenvalue $\gamma(X)$: $$n \gtrsim \frac{\log d}{\Delta_u } \cdot \frac{1}{\gamma^2}$$ compared to the information-theoretic optimum $n \gtrsim (\log d)/\Delta_l$ achieved by BSS (Gao et al., 5 Oct 2025). The restricted eigenvalue is given by: $$\gamma(X) = \min_{S \in \mathcal{S}_{d,s}} \min_{\|\theta_{S^c}\|_1 \le 3 \|\theta_S\|_1} \frac{ \|X\theta\|^2/n }{ \|\theta\|^2 }.$$ A small $\gamma$ (high covariance or strong dependence) makes the statistical-computational gap large.

For structured matrix estimation problems, such as sparse principal component estimation, convex relaxations are statistically suboptimal in moderate sample regimes. For example, the semidefinite relaxation for sparse PCA achieves error $\sqrt{(k^2 \log p)/(n\theta^2)}$, losing a $\sqrt{k}$ factor relative to the minimax error rate $\sqrt{(k \log p)/(n\theta^2)}$ ((Wang et al., 2014), under the planted clique hypothesis). For linear regression, Lasso can also be rate-suboptimal in correlated design due to "signal cancellation" (Jin et al., 2012).

4. Structure-Exploiting Procedures: Balancing Statistical and Computational Performance

Efficient procedures that exploit additional structure—such as sparsity in the Gram matrix or group dependence—can achieve near-optimal rates with feasible computation:

• Graphlet Screening (GS) utilizes thresholded Gram matrix to form a sparse Graph of Strong Dependence (GOSD), restricting multivariate screening to connected subgraphs. GS's screening cost scales as $O(np \cdot (\log p)^{(m_0+1)\alpha})$ and its cleaning step operates on small disconnected clusters ("graphlets"), which enables it to match the minimax Hamming risk while remaining nearly linear in $p$ (Jin et al., 2012):

$$\text{Hamm}_p^*(\theta, \kappa, r, a, \Omega) \asymp L_p p^{1 - \frac{(\theta + r)^2}{4r}}$$

GS is provably optimal in a wide class of sparse signal/sparse-graph regimes, whereas the global Lasso or subset selection is not.

• Adaptive Subspace Methods and Local Screening: Methods like AdaSub break the high-dimensional space into adaptively chosen low-dimensional problems, "zooming in" on important variables over many iterations. These strategies avoid the combinatorial explosion of BSS, benefiting from favorable selection criteria (e.g., EBIC), and empirically achieve strong accuracy with feasible computation (Staerk et al., 2019).
• Estimate-Then-Screen (ETS): For weak and heterogeneous signals, two-stage frameworks decouple computation and support detection: (1) estimate $\hat{\beta}$ with polylog cost, (2) coordinate-wise screening, often on an independent sample (Roy et al., 2022). ETS achieves model consistency at the information-theoretic optimal signal threshold $r>1$, and, critically, avoids the exponential complexity of BSS.
• Test-Based Selection: Sequential test-based variable selection (e.g., using maximal partial correlation and asymptotic null theory) achieves low error and fast computation, replacing costly cross-validation with statistical testing at each inclusion step (Gong et al., 2017).

5. Empirical Bayes, Variational, and MCMC Approaches

• Empirical Bayes with variational marginalization significantly reduces computation by working directly on the model space, preserving selection consistency and requiring tuning over only $p$ variational parameters (Bernoulli families over indicators), rather than scaling with model size (Tang et al., 14 Feb 2025). Variational solutions match MCMC posterior inclusion probabilities at a fraction of the runtime.
• For fully Bayesian approaches using MCMC, posterior concentration (statistical guarantee) does not imply rapid MCMC mixing (computational guarantee). Truncated sparsity priors (model-size constraints) and carefully designed proposal moves (e.g., single/double flips in Metropolis-Hastings) are essential for both statistical accuracy and polynomial mixing time (Yang et al., 2015). In their absence, the Markov chain may mix exponentially slowly despite correct model-specific concentration.

6. Modern Directions: Synthetic Data, Diffusion Models, and FDR Control

• Diffusion-Driven Variable Selection utilizes pretrained/fine-tuned diffusion models to generate synthetic datasets, applies standard selectors (e.g., Lasso) to each, and aggregates inclusion indicators. This procedure stabilizes variable selection, especially under high correlation, and supports valid inference; it achieves sign consistency under mild conditions while being computationally efficient via parallelism and transfer learning (Wang et al., 19 Aug 2025).
• Synthetic Null Parallelism (SyNPar) controls FDR by parallel model fitting on real and synthetic null data (generated under the null hypothesis), selecting features based on coefficient magnitude comparisons. SyNPar guarantees FDR control and asymptotically full power, outperforming both knockoff- and data-splitting approaches in both statistical accuracy and runtime (Wang et al., 9 Jan 2025).
• Penalized Criterion Calibration: Careful parameter tuning (e.g., penalty constant $K$ in penalized least-squares) can simultaneously control predictive risk and FDR. Non-asymptotic bounds and data-driven selection of $K$ enable practitioners to balance selection conservativeness and predictive performance (Lacroix et al., 2023).

7. Summary Table: Statistical-Computational Trade-off Landscape

Method	Statistical Rate	Computational Complexity	Achievability Condition	Reference
Best subset selection	Minimax optimal	Exponential in $p$, NP-hard	Any design	(Gao et al., 5 Oct 2025)
Lasso / SDP relax.	Suboptimal (worse in $p$ exponent)	Polynomial time	Strong RE/incoherence, large sample	(Gao et al., 5 Oct 2025, Wang et al., 2014, Jin et al., 2012)
Graphlet Screening	Minimax optimal	$O(np ~\mathrm{polylog}(p))$	Sparse Gram, signals decompose on GOSD	(Jin et al., 2012)
AdaSub / Local search	Empirically strong	Polynomial (adaptive)	Depends on selection criterion and OIP	(Staerk et al., 2019)
ETS framework	Minimax optimal	Polylog time per est., linear	Two-stage, accurate initial estimator	(Roy et al., 2022)
SyNPar	Controls FDR, high power	2 model fits (linear in $p$)	Mild (data generative fit, parametric)	(Wang et al., 9 Jan 2025)
Diffusion-agg. methods	Sign consistency	Parallelizable, GPU-efficient	Accurate generative model for $X$ and $Y$	(Wang et al., 19 Aug 2025)

8. Fundamental Insights and Limitations

• The statistical-computational gap manifests when "easy" instances (e.g., i.i.d. designs, strong signals) permit efficient minimax selection, but complex dependencies (e.g., sparse strong dependence, low RE constant) increase the necessary sample size for efficient methods (a $1/\gamma^2$ penalty in $n$ is typical). Only combinatorial algorithms surmount this, at intractable cost.
• Structural assumptions (bandedness, sparsity, block structure) are key for scalable and statistically efficient methods.
• Empirical Bayes and variational approaches can inherit selection consistency if hyperparameters and approximation choices are properly calibrated.
• For model selection with explicit FDR control, recent approaches leveraging synthetic nulls or diffusion-augmented resampling provide both rigorous error control and practical acceleration.

9. Outlook

Despite substantial progress, sharp delineation of the statistical-computational boundary remains an active topic, especially:

• Under alternative signal regimes (rare/weak signals, heterogeneity)
• When dealing with more complex dependence (non-Gaussian, latent factors)
• In integration of generative AI (diffusion models) for robust variable selection and model diagnostics
• For adaptive procedures that can certify statistical guarantees jointly with computational resource bounds.

The statistical-computational trade-off therefore not only structures the attainable rates and feasibility of support recovery but actively informs methodological design and theoretical understanding in contemporary high-dimensional inference.