Both other answers are wrong and misunderstand the Meinshausen–Bühlmann approach. The method is implemented correctly, but does not return what you might expect. It is also not necessarily faster.
A quick answer
If you use approx = TRUE
, proceed as follows:
A <- glasso(S, rho = 0, approx = TRUE)$wi
Pcor <- 1/2 * (A + t(A))
That's it. You now have a matrix of approximate partial correlations.
The covariance matrix is not estimated by this method, which is why $w
is all zero.
Comparison to an exact solution
Let's use some example data and $\lambda = 0$ to see that this does indeed do what we want.
Example data:
set.seed(2024)
n <- 20
p <- 3
X <- matrix(rnorm(n * p), nrow = n, ncol = p)
Xc <- scale(X, center = TRUE, scale = FALSE)
S <- 1/n * crossprod(Xc)
Meinshausen–Bühlmann approach:
library("glasso")
MB <- glasso(S, rho = 0, approx = TRUE)
(Pcor <- 1/2 * (MB$wi + t(MB$wi)))
# [,1] [,2] [,3]
# [1,] 0.0000000 0.2150574 -0.2212086
# [2,] 0.2150574 0.0000000 -0.1252797
# [3,] -0.2212086 -0.1252797 0.0000000
Exact solution:
P <- solve(S)
P <- -P
diag(P) <- -diag(P)
(Pcor <- cov2cor(P))
# [,1] [,2] [,3]
# [1,] 1.0000000 0.2057225 -0.2164845
# [2,] 0.2057225 1.0000000 -0.1247561
# [3,] -0.2164845 -0.1247561 1.0000000
Some background
Inverting a large matrix can be slow, even with the coordinate descent algorithm used in glasso
. The option approx = TRUE
uses an alternative based on neighborhood selection (Meinshausen & Bühlmann, 2006). Importantly, it does not estimate the empirical covariance matrix $\mathbf{S}$, nor the precision matrix $\mathbf\Theta$. Instead, it relies on the observation that partial correlations are closely related to regression coefficients (1), and attempts to estimate non-zero off-diagonal elements directly using (LASSO-penalized) regression.
You can demonstrate that this is what happens by implementing your own version:
library("glmnet")
A <- matrix(0, nrow = p, ncol = p)
for(i in 1:p){
y_current <- X[, i]
X_current <- X[, -i]
LASSO <- glmnet(X_current, y_current, alpha = 1, lambda = 0)
A[-i, i] <- coef(LASSO)[-1]
}
(Pcor <- 1/2 * (A + t(A)))
# [,1] [,2] [,3]
# [1,] 0.0000000 0.2150574 -0.2212086
# [2,] 0.2150574 0.0000000 -0.1252797
# [3,] -0.2212086 -0.1252797 0.0000000
(Compare this to the output of glasso(..., approx = TRUE)
above.)
The asymmetry in the solution is a consequence of the fact that regressing y ~ x
$\neq$ x ~ y
. Some papers try using an OR rule or AND rule to decide which elements should be non-zero, like in Banerjee et al. (2008), but the symmetrization trick shown here provides a much more reasonable approximation in my opinion.
When should you use Meinshausen–Bühlmann (approx = TRUE
)?
The main advantage over the block-coordinate descent algorithm in glasso is that each element $\theta_{ij}$ is estimated separately, so this is very easy to parallelize. However, Witten et al. (2011) showed that there is a simple way to check whether $\mathbf\Theta$ is block-diagonal for a given value of $\lambda$. If it is, then we can separately estimate each block $k \in K$ of $$\mathbf\Theta = \begin{pmatrix} \mathbf\Theta_1 & & & \\ & \mathbf\Theta_2 & & \\ & & \ddots & \\ & & & \mathbf\Theta_K\end{pmatrix}.$$
The authors conclude:
Strikingly, when the graphical lasso is performed with a large tuning parameter value, then the algorithms proposed in this article lead to such massive speed improvements that computing the exact graphical lasso solution is much faster even than computing the Meinshausen and Buhlmann (2006) approximate solution.
From glasso
1.7 on, if you expect $\mathbf\Theta$ to be sparse, approx = FALSE
can be orders of magnitude faster.
References
Meinshausen, N., & Bühlmann, P. (2006). High-dimensional graphs and variable selection with the Lasso. Annals of Statistics, 34(3), 1436–1462. https://doi.org/10.1214/009053606000000281
Banerjee, O., El Ghaoui, L., & d’Aspremont, A. (2008). Model Selection Through Sparse Maximum Likelihood Estimation for Multivariate Gaussian or Binary Data. J. Mach. Learn. Res., 9, 485–516. https://www.jmlr.org/papers/volume9/banerjee08a/banerjee08a.pdf
Witten, D. M., Friedman, J. H., & Simon, N. (2011). New Insights and Faster Computations for the Graphical Lasso. Journal of Computational and Graphical Statistics, 20(4), 892–900. http://www.jstor.org/stable/23248939
glasso(var(X), rho = 0.01, zero = info[-c(1000:1950000),], approx=TRUE)
, the estimated cov matrix is still zero!! $\endgroup$