7
$\begingroup$

Sorry if this is too damned long. I did what I could to abbreviate it.

The question is about Least Angle Regression (LARS).

I'm new to numerical work with matrices.

I believe I have a way to compute Least Angle without explicit matrix inversion. I'm hoping someone with more experience can check my work to see if I'm doing anything dangerous, illegal, or otherwise ignorant.

This paper: Least Angle Regression, says:

"...we compute $m−k$ inner products $c_{jk}$ of the non-active $x_j$ with the current residuals to identify the next active variable, and then invert the $k \times k$ matrix $\mathcal{G}_k = X'_kX_k$ to find the next LARS direction."

(Emphasis mine.)

They update the Cholesky decomposition in each iteration, so it seems they're doing something like this:

  • At $k$th iteration:
  • $X_k$ is a matrix formed from $k$ predictors.
  • $X'_k X_k = \mathcal{G}_{k}$
  • We already have $\mathcal{G}_{k-1}$, and $\mathcal{G}_{k} = \mathcal{G}_{k-1}$ except for a new row and column of dot products, so let's say $$\mathcal{G}_{k} = \left[ \begin{array}{cc} \mathcal{G}_{k-1} & \mathbf{p} \\ \mathbf{p}' & x \end{array} \right]$$
  • We can find the Cholesky decomposition $L_k$ from $L_{k-1}$: $$L_{k} = \left[ \begin{array}{cc} L_{k-1} & \mathbf{r}=L_{k-1}^{-1}\mathbf{p} \\ \mathbf{r}' & y=(x - \mathbf{r}\cdot\mathbf{r})^{-\frac12} \end{array} \right]$$

  • Find $L_k^{-1}$. Note this explicit inversion.

  • $\mathcal{G}^{-1}_k = {L'}_k^{-1} L_k^{-1}$
  • $A_k = (1'_k \mathcal{G}^{-1}_k 1_k)^{-\frac12}$, where $1_k$ is a vector of $k$ 1s.
  • ${\mathbf w}_k = A_k {\mathcal G}_k^{-1} 1_k$
  • Use ${\mathbf w}_k$ to select the $(k+1)$th predictor.

Now here's my approach:

  • At $k$th iteration:
  • Find $L_k$, $\mathbf{r}$, and $y$ as above
  • Define $L_k {\mathbf q}_k = 1_k$
  • Find $\mathbf{q}_k$: $${\mathbf q}_{k} = \left[ \begin{array}{c} {\mathbf q}_{k-1} \\ \text{$\scriptsize{(1 - {\mathbf q}_{k-1} \cdot \mathbf{r})/y}$} \end{array} \right]$$

  • $\begin{align}A_k &= (1'_k \mathcal{G}^{-1}_k 1_k)^{-\frac12} \\ &= (1'_k L'{_k}^{-1} L_k^{-1} 1_k)^{-\frac12} \\ &= \left[(L_k^{-1} 1_k)' (L_k^{-1} 1_k)\right]^{-\frac12} \\ &= ({\mathbf q}'_k {\mathbf q}_k)^{-\frac12} \end{align}$

  • $\begin{align}{\mathbf w}_k &= A_k {\mathcal G}_k^{-1} 1_k \\ &= A_k L'{_k}^{-1} {\mathbf q}_k \end{align}$
  • So solve for ${\mathbf w}_k$: $\frac{1}{A_k} L'_k {\mathbf w}_k = {\mathbf q}_k$

Comments? Does it even matter if I'm right?

$\endgroup$
9
  • $\begingroup$ This is a very specialized question, so it may be difficult or take a while to get a good answer. Do you seem to get the same answers with both approaches? $\endgroup$
    – Jason R
    Jun 5, 2014 at 0:44
  • 1
    $\begingroup$ As I see $L_k$ is a matrix where as $1_k$ is a column vector. Therefore, you cannot break the matrix into two vectors and then take inverse of vectors in 3rd step of decomposing $A_k$. $\endgroup$
    – learner
    Jun 5, 2014 at 2:30
  • $\begingroup$ Yes Jason, I get the same answers with both approaches. learner, I have removed the offending line. It was inconsequential. I have also rearranged the steps to match what my implementation actually does. $\endgroup$ Jun 5, 2014 at 22:54
  • $\begingroup$ Your link to the paper is broken $\endgroup$
    – Phonon
    Jun 6, 2014 at 17:50
  • $\begingroup$ That's weird. Somehow the address was abbreviated with ellipses. Fixed now. $\endgroup$ Jun 7, 2014 at 3:21

1 Answer 1

2
$\begingroup$

If you want to solve for single value of $ \lambda $ in the model:

$$ \arg \min_{x} \frac{1}{2} {\left\| A x - b \right\|}_{2}^{2} + \lambda {\left\| x \right\|}_{1} $$

Then you can use Coordinate Descent method which is the fastest and simplest and doesn't require any matrix inversion.

I have a MATLAB code for in my $ {L}_{1} $ Regularized Least Squares Analysis GitHub Repository:

enter image description here

$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.