Comparing Likelihood And Neighborhood Selection

Two Routes to a Sparse Gaussian Graph

For a multivariate Gaussian distribution, a zero in the precision matrix corresponds to conditional independence between two variables given the others. The graphical lasso and neighborhood selection use this same relationship but approach estimation differently. The graphical lasso estimates one symmetric precision matrix by penalized likelihood. Neighborhood selection fits a separate penalized regression for each variable, then combines the resulting supports into a graph.

The distinction matters in practice. Joint likelihood estimation enforces symmetry automatically, whereas separate regressions can disagree about whether an edge is present. Their computational costs and support-recovery guarantees also depend on different features of the problem. Neither method dominates in every setting.

The Objectives Estimate Different Objects Directly

With a common penalty convention, the graphical-lasso estimate minimizes the negative Gaussian log-likelihood, up to constants and a positive factor, plus an \(L_1\) penalty on off-diagonal entries:

Here \(\Theta\succ0\) means that \(\Theta\) is symmetric positive definite. The sum counts both \((j,k)\) and \((k,j)\); an alternative convention sums only over \(j<k\), changing the numerical penalty scale by a factor of two. Because the optimization variable is one symmetric matrix, the selected edge support is symmetric by construction.

Neighborhood selection instead regresses each variable on all the others. For each \(j\), with centered data, it solves

The estimate \(\hat{\beta}^{(j)}_k\) describes the fitted contribution of variable \(k\) in the regression for \(j\). In general, \(\hat{\beta}^{(j)}_k\) and \(\hat{\beta}^{(k)}_j\) need not agree in whether they are zero. An edge rule is therefore needed: the OR rule includes an edge if either regression selects it, while the AND rule includes it only if both do.

Proof. The Gaussian density is proportional to \(\exp(-x^{\mathsf T}\Theta x/2)\). Fix \(x_{-j}\) and collect the terms involving \(x_j\) in the exponent's quadratic form. Those terms are

Thus the conditional distribution of \(X_j\) given \(X_{-j}\) is normal with mean \(-\sum_{k\ne j}(\Theta_{jk}/\Theta_{jj})X_k\) and variance \(1/\Theta_{jj}\). Its regression coefficients are therefore \(-\Theta_{jk}/\Theta_{jj}\). Since \(\Theta_{jj}>0\), a coefficient is zero exactly when the corresponding off-diagonal precision entry is zero. \(\square\)

Symmetry Is Automatic in One Method, a Decision in the Other

The likelihood estimate has one entry for each unordered pair, so \(\hat{\Theta}_{jk}=\hat{\Theta}_{kj}\). Neighborhood selection has one coefficient from each of two different regressions for a candidate edge. Finite samples, separate tuning choices, and numerical optimization can produce different support decisions in those regressions, even though the corresponding population coefficients are both zero or both nonzero.

The OR and AND rules have a simple set relationship. Let \(A_j\) be the set of variables selected in the regression for \(j\). For each unordered pair \(\{j,k\}\), the OR rule selects it when \(k\in A_j\) or \(j\in A_k\); the AND rule selects it when both conditions hold. Therefore every AND edge is an OR edge, but not conversely. The OR rule is more inclusive; the AND rule is more conservative. Neither rule repairs an incorrect selection made by both regressions.

Computational Trade-Offs

Neighborhood selection solves \(p\) lasso problems, each with \(p-1\) candidate predictors. The regressions can be solved in parallel and can use standard lasso algorithms. This can be attractive when parallel computing is available or when the main goal is a list of candidate neighbors for each node. Its cost is not automatically negligible: for large \(p\), there are many regressions, and a dense implementation may repeatedly process much of the same data.

The graphical lasso solves one coupled positive-definite matrix optimization problem. Its iterates and matrix operations involve \(p\times p\) objects, and enforcing positive definiteness is part of the optimization. Joint estimation can be computationally demanding for large \(p\), although sparse structure, algorithm choice, and warm starts can substantially affect runtime.

A count of the number of regressions versus the number of matrix entries is not a reliable runtime comparison. Sample size, sparsity, stopping tolerances, software, available parallelism, and whether intermediate matrices are dense all matter. The practical comparison should measure runtime and memory for the intended dimensions and implementation.

Support-Recovery Guarantees Are Conditional

Support recovery means correctly identifying which off-diagonal entries of \(\Theta\) are zero and which are nonzero. It is distinct from estimating all entries accurately: a small estimation error does not by itself guarantee correct edge decisions when a true nonzero entry is very small. Both methods therefore need a separation condition, commonly expressed as a beta-min condition, as well as assumptions on the design or precision structure.

For the graphical lasso, standard support-consistency results for Gaussian models use conditions that include bounded eigenvalues, sufficiently sparse rows (with maximum degree \(d\)), an appropriate incoherence condition for the likelihood's Hessian, and a minimum nonzero precision entry large enough relative to the estimation error. A common sufficient scaling in such results is \(n\) of order at least \(d^2\log p\), up to constants and problem-specific conditions. The estimate is symmetric, so exact support recovery returns an undirected edge set directly.

For neighborhood selection, support-consistency results are typically applied to each conditional regression. They require suitable design regularity and an irrepresentability condition for the active predictors, along with sparsity and a beta-min condition for the nonzero regression coefficients. Under Gaussian assumptions and appropriate boundedness conditions, a common sufficient sample-size scaling is of order \(d\log p\), with constants depending on the design conditions. Recovery across all nodes also requires controlling errors across the collection of regressions and choosing a rule for combining their supports.

These are representative sufficient regimes, not universal thresholds or a guarantee that neighborhood selection is always more sample-efficient. The assumptions differ, constants can be consequential, and violations can invalidate the guarantees. The regression coefficients are scaled versions of precision entries, \(-\Theta_{jk}/\Theta_{jj}\); translating a minimum-signal condition between the methods requires control of the diagonal entries as well.

Common Mistake

A common mistake is to say that the two methods must select exactly the same finite-sample graph because they encode the same Gaussian conditional independences. The population relationship does not imply equality of penalized estimates: the likelihood couples all precision entries and enforces positive definiteness, while neighborhood selection penalizes each regression separately. Their tuning parameters and finite-sample errors therefore need not yield identical supports. The correction is to distinguish the shared population graph from the different finite-sample procedures, and to report the neighborhood-selection symmetrization rule whenever that method is used.