Concave Penalized Estimation of Sparse Gaussian Bayesian Networks

TL;DR

This paper introduces a fast, score-based algorithm for estimating sparse Gaussian Bayesian networks using concave regularization, capable of handling high-dimensional data without search space restrictions, with theoretical guarantees and superior empirical performance.

Contribution

The authors propose a novel concave penalized likelihood framework and a fast algorithm for high-dimensional Gaussian Bayesian network structure learning, without search space restrictions and with theoretical guarantees.

Findings

Algorithm is significantly faster than competitors.

Achieves higher sensitivity with similar false discovery rates.

Can handle networks with thousands of variables in about one hour.

Abstract

We develop a penalized likelihood estimation framework to estimate the structure of Gaussian Bayesian networks from observational data. In contrast to recent methods which accelerate the learning problem by restricting the search space, our main contribution is a fast algorithm for score-based structure learning which does not restrict the search space in any way and works on high-dimensional datasets with thousands of variables. Our use of concave regularization, as opposed to the more popular $\ell_0$ (e.g. BIC) penalty, is new. Moreover, we provide theoretical guarantees which generalize existing asymptotic results when the underlying distribution is Gaussian. Most notably, our framework does not require the existence of a so-called faithful DAG representation, and as a result the theory must handle the inherent nonidentifiability of the estimation problem in a novel way. Finally, as a matter of independent interest, we provide a comprehensive comparison of our approach to several standard structure learning methods using open-source packages developed for the R language. Based on these experiments, we show that our algorithm is significantly faster than other competing methods while obtaining higher sensitivity with comparable false discovery rates for high-dimensional data. In particular, the total runtime for our method to generate a solution path of 20 estimates for DAGs with 8000 nodes is around one hour.