Using sparse inverse covariance matrix in estimating least squares coefficients

Question

I am reading the paper introducing the graphical lasso, which is a way to estimate a sparse inverse covariance matrix.

http://www-stat.stanford.edu/~tibs/ftp/graph.pdf

Finding a sparse inverse covariance matrix can be useful because 0's in the inverse of the covariance matrix imply that two variables are independent. In this paper, it is mentioned that it is related to the least squares solution because $$B = (X^T X)^{-1} X^T Y$$. In this solution $(X^TX)^{-1}$ can be thought of as the inverse covariance matrix.

I am curious about two things:

1) Say I have some data where I suspect that variables of $X$ are independent of some other variables in $X$. Perhaps computing the empirical covariance matrix of $X$ doesn't reflect that idea, so instead I try to find a sparser inverse covariance matrix using the graphical lasso. Now, I introduce a response variable $Y$ and want to model a linear relationship that takes into consideration the idea that some variables in $X$ are completely independent of other variables in $X$. Does this mean I can solve a "graphical" least squares problem by taking a response matrix $Y$ and a design matrix $X$, and computing the sparse inverse covariance matrix of $X$ (call this $W$) then multiply by the crossproduct $X^T Y$, ie $$B = W X^T Y$$?

When computing the crossproduct $X^T Y$, should I still use the original design matrix $X$, or should it be some function of $W$ now?

2) I am curious about the case where the response is a tuple, maybe the response $Y$ is a three column matrix. Suppose I want to measure a linear model for $Y$ with a design matrix $X$ again, but this time I just want sparsity in the regression coefficients, taking us back to the standard lasso problem. I could iterate over each column of $Y$ (call these $j_1, j_2, j_3$), and run a lasso with $j_i$ and $X$ which would yield 3 coefficient vectors, each vector having sparsity. Maybe instead I would like to estimate an entire coefficient matrix, where the matrix as a whole is sparse but each column individually need not be sparse. The paper has this to say: (does this mean solving for a sparse B matrix is more appropriately done using a graphical lasso?)

Blockquote There is a simple, conceptually appealing way to view this procedure. Given a data matrix X and outcome vector y, we can think of the linear least squares regression estimates $(X^TX)^{-1}X^TY$ as functions not of the raw data, but instead the inner products $X^TX$ and $X^TY$ . Similarly, one can show that the lasso estimates are functions of these inner products as well. Hence in the current problem, we can think of the lasso estimates for the pth variable on the others as having the functional form lasso(S11; s12; ): (9) But application of the lasso to each variable does not solve problem (1); to solve this via the graphical lasso we instead use the inner products W11 and s12. That is, we replace (9) by lasso(W11; s12; ): (10) The point is that problem (1) is not equivalent to p separate regularized regression problems, but to p coupled lasso problems that share the same W and $\theta = W^{-1}$ . The use of W11 in place of S11 shares the information between the problems in an appropriate fashion. Blockquote

Answer 1

Bias-Variance Tradeoff

Assuming the coefficents' estimator has the form

$$\hat{\beta} = WX^TY$$

where $W$ is some arbitrary matrix, we have

$$E[ WX^TY ] = E[ WX^T( X\beta + \epsilon) ] = W(X^TX)\beta$$

Hence we only have an unbiased estimator if $W = (X^TX)^{-1}$, ie. $\hat{\beta} = \beta_{ols}$. This is ok, as unbiased-ness is overrated and we typically can acheive greater variance reductions by ignoring it. But do we reduce the variance of $\hat{\beta}$ by using the $p$-penalised sparse matrix, $\Theta_p$, from the Glasso output?

Considering the optimization problem of the Glasso with penality $p$, we know

$$|| \Theta_p ||_1 \le \frac{1}{p} \le || (X^TX)^{-1} ||_1 = s$$

WLOG, assume that $Var( \epsilon) = 1$. Set $W = \Theta_p$ in $\hat{\beta}$. Then $$Var( \beta_{ols}) = (X^TX)^{-1}$$ $$\Rightarrow || Var( \beta_{ols} ) ||_1 = s$$

More generally,

$$Var( \hat{\beta} ) = W( X^TX )W^T $$

Note the relation that as $p \rightarrow 0$, $\Theta_p \rightarrow (X^TX)^{-1} $. On the other hand, $p \rightarrow \infty, \; \Theta_p \rightarrow 0$. This second fact implies that

$$||Var(\hat{\beta})||_1 \rightarrow 0, \; \text{ as } p\rightarrow \infty $$

It can be shown that $||Var(\hat{\beta})||_1$ is monotonic in $p$, hence we have that the variance is less than the OLS variance for all $p>0$, where less is under the $L1$ norm.

Furthermore, we can bound the rate at which it goes to zero:

$$ ||Var(\hat{\beta})||_1 = ||W( X^TX )W^T||_1 \le \frac{|| X^TX ||_1}{p^2} $$

This shows that most of the variance reduction occurs very quickly, and then does not reduce by much thereafter. Keep in mind, as we increase $p$, our bias increases too.

So you have the bias-variance tradeoff: increase bias but decrease variance by using $\Theta_p$ instead of $(X^TX)^{-1}$.

$X^TY$ or $f(W)Y$

You make a good point:

When computing the crossproduct $X^TY$, should I still use the original design matrix $X$, or should it be some function of $W$ now?

• There is no intuitive reason why the crossproduct matrix should not be $X^T$,
• If you were to use a function of $W$, the most natural candidate would be to use the Cholesky decomposition of $W^{-1}$.

Answer 2

I'm curious about what maximization problem this new estimator, with the Glasso output used instead of $(X^TX)^{-1}$, is trying to solve. This may give us insight about the use of it. Assume the coefficients $\beta$ have the following form: $$ \beta = M X^T Y$$ for some unknown $M$. It can be shown that $M = (X^TX)^{-1}$ for the standard least squares problem. In fact, this solution is the solution to the following optimization problem:

$$\min_M \; (Y - X\beta)^T(Y - X\beta)$$ $$ = \min_M \; (Y - X(MX^TY))^{T}(Y - X(MX^TY)) $$

assuming $\beta$ has the form above. If you are careful with your matrix calculus, like I said you can show that the solution is $M = (X^TX)^{-1}$.

I think what you are asking is equivilant to, again assuming $\beta$ has the same form as above, what is the solution to

$$ \min_M \; (Y - XMX^TY)^{T}(Y - XMX^TY) + \lambda||M||_1 $$

Is this is the same solution as the Graphical Lasso? It is possible. One could argue against this by saying that this problem has a $Y$ vector, while the Glasso algorithm does not include any response vector. But, when performing the matrix calculus for the first problem, the $Y$ vector drops out (hence $M$ is ind. of $Y$), and that could happen here too.

To prove this you would have to play around with the above opt. program and show it is the same as

$$ \max_ { W} \;\; \log{ \det{ W} } - Tr(SW) - \alpha||W||_1$$

where $M = W^{-1}$. I don't think this is true. So perhaps using the Glasso output is not an optimal choice, but may be good nonetheless.