Variational approximation for heteroscedastic linear models and matching pursuit algorithms

Modern statistical applications involving large data sets have focused attention on statistical methodologies which are both efficient computationally and able to deal with the screening of large numbers of different candidate models. Here we consider computationally efficient variational Bayes approaches to inference in high-dimensional heteroscedastic linear regression, where both the mean and variance are described in terms of linear functions of the predictors and where the number of predictors can be larger than the sample size. We derive a closed form variational lower bound on the log marginal likelihood useful for model selection, and propose a novel fast greedy search algorithm on the model space which makes use of one step optimization updates to the variational lower bound in the current model for screening large numbers of candidate predictor variables for inclusion/exclusion in a computationally thrifty way. We show that the model search strategy we suggest is related to widely used orthogonal matching pursuit algorithms for model search but yields a framework for potentially extending these algorithms to more complex models. The methodology is applied in simulations and in two real examples involving prediction for food constituents using NIR technology and prediction of disease progression in diabetes.

💡 Research Summary

The paper addresses the growing demand for statistical methods that are both computationally efficient and capable of handling massive model spaces, a situation common in modern high‑dimensional data analysis. The authors focus on heteroscedastic linear regression, where both the conditional mean and the conditional variance are expressed as linear functions of the covariates, and where the number of potential predictors can far exceed the sample size.

A Bayesian formulation is adopted: the mean is μ_i = x_i^Tβ and the variance is σ_i^2 = exp(z_i^Tγ). Conjugate priors (multivariate normal for β and inverse‑Gamma for γ) are placed to retain analytical tractability. The core methodological contribution is a variational Bayes (VB) approximation that assumes a factorised q‑distribution of the form Normal‑Inverse‑Gamma. By carefully manipulating the evidence lower bound (ELBO), the authors derive a closed‑form expression that depends only on the first‑ and second‑order moments of β and γ. This enables simple coordinate‑ascent updates that involve only matrix multiplications, making the algorithm scalable to thousands of predictors.

Beyond inference, the paper proposes a novel greedy model‑search strategy that leverages a “one‑step ELBO update.” For each candidate variable, the algorithm computes the change in the ELBO that would result from adding or removing that variable, using only a single VB update. The variable yielding the largest ELBO increase is selected, and the process repeats until a stopping criterion is met. Because the ELBO can be evaluated cheaply, the search scales linearly with the number of candidates, allowing rapid screening of very large predictor sets.

The authors point out that this procedure is closely related to Orthogonal Matching Pursuit (OMP). While OMP selects variables based on the inner product between residuals and predictors, the proposed method selects variables based on the full variational lower bound, thereby incorporating both fit and uncertainty. Consequently, the algorithm inherits the computational simplicity of OMP while retaining a fully probabilistic interpretation and the ability to incorporate prior information.

Simulation experiments demonstrate that, in p ≫ n settings, the VB‑OMP hybrid outperforms standard penalised‑likelihood approaches (LASSO, Elastic Net) and existing VB methods in terms of variable‑selection accuracy and mean‑squared prediction error. The method also successfully recovers the heteroscedastic variance structure, which many competing techniques ignore.

Two real‑world applications illustrate practical impact. In a near‑infrared (NIR) spectroscopy study, the algorithm selects a parsimonious set of wavelengths that dramatically improves prediction of food constituent concentrations. In a longitudinal diabetes cohort, it identifies key clinical and biochemical markers that explain both the mean trajectory and the variability of disease progression, offering interpretable insights for clinicians.

In conclusion, the paper delivers a unified framework that combines variational Bayesian inference with a greedy matching‑pursuit style search, enabling efficient estimation, variable selection, and model comparison for high‑dimensional heteroscedastic linear models. The authors suggest extensions to non‑linear settings, multivariate responses, and online streaming data as promising avenues for future research.