5
$\begingroup$

What is most efficient way to updating cholesky factorization of for removing a column from the matrix?

T = Chol(X' *X)

If I remove a column from X, how to effienciently update cholesky. It is trivial to remove the last column of X as new cholesky will be removing last row and column form T. What about removing other column?

$\endgroup$
3
  • 1
    $\begingroup$ I don't think the question is Bayesian in any way. The problem you're applying this to may be, but the question contains no trace of it, and the problem you do ask about comes up in frequentist statistics too. $\endgroup$
    – Glen_b
    Commented Feb 8, 2014 at 10:01
  • $\begingroup$ Here the algorithm description: normalsplines.blogspot.com/2019/02/… $\endgroup$
    – NetTvor
    Commented Feb 11, 2019 at 10:04
  • $\begingroup$ See also: math.stackexchange.com/questions/1896467/… $\endgroup$
    – Damien
    Commented Aug 16, 2019 at 5:41

1 Answer 1

12
$\begingroup$

The usual approach involves interchanging the position of the column you want removed and the last one, and returning the Cholesky to triangular form via either Givens rotations or Householder transforms, whereupon you can simply drop the last one.

Similarly you can update an additional column into any position by adding to the end and then swapping into the desired position followed by the Givens/Householder transformations to fix it back up to triangular form.

When you do the swap of columns, you destroy the triangular form you need it to have. Each of the Givens rotations (or Householder transforms) then fixes up one of the values that are 'sticking out', until they are all fixed. This is quite a lot faster than recomputing from scratch.

I'll talk about rotations (but the Householder transforms do the same job, so much of the discussion carries over directly). For each non-zero value that's 'sticking out', you apply a Givens rotation (a 2x2 matrix of $\sin\theta$ and $\cos\theta$ values representing a rotation through $\theta$), so that it 'zeroes' that one value. It changes other values, but retains the $A=RR'$ property. By starting with the value that's sticking out the most and working in toward the diagonal of the triangle, you end up not losing the values you choose to zero out - that is, if you do things in the right order, the zeros you make at one step are still zeros the next step.

The Householder approach is a bit more efficient.


References. Hmm.

Matrix Computations (Golub and Van Loan) has a section on updating calculations - 12.5.2 "Adding or deleting a column" - it mostly talks about QR but much of the approaches and ideas carry over.

The documentation of the LINPACK routine SCHEX (*) and some of the other routines was actually where I first saw (and eventually understood, many, many years ago now) how to do it (if you have access to a library an old users guide may be available). I found its explanation relatively straightforward at the time.

* that function simply permutes the order of the variables in a Cholesky; this uses (I think) Givens rotations, which is relatively easy to find information on.

In fact, the documentation right in the code you can see online** is fairly informative. I am not sure the documentation said very much more than that. The description is clearly using Givens rotations.

** I mean the block of comments at the top

I'm sure I have something else on it somewhere; I'll have to see what other references I have.

$\endgroup$
2
  • $\begingroup$ Thanks Glen_b! If you can give reference so that I can implement this it will be great. $\endgroup$
    – avi
    Commented Feb 7, 2014 at 17:04
  • $\begingroup$ Thanks Glen. This is very informative. I am surprised there is not much information online and references considering it is useful everywhere in optimization and sampling. Please post here whenever you have other references. $\endgroup$
    – avi
    Commented Feb 8, 2014 at 17:05

Your Answer

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

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