Updating the Cholesky factorization of $A'A=LL'$ when one or more columns are added to or removed from matrix $A$ can be done very efficiently obviating the re-factorization from scratch.

Removing columns

Assume we know the Cholesky factorisation of $A'A=LL'$ and we remove a column from matrix $A$ which can be written as follows

$A = \begin{bmatrix} A_{1:\kappa} & c & A_{\kappa+1:n} \end{bmatrix}$

Then,

\begin{aligned} A'A &= \begin{bmatrix} A_{1:\kappa} & c & A_{\kappa+1:n} \end{bmatrix}'\begin{bmatrix} A_{1:\kappa} & c & A_{\kappa+1:n} \end{bmatrix}\\ &=\begin{bmatrix} A_{1:\kappa}'A_{1:\kappa} & A_{1:\kappa}'c & A_{1:\kappa}'A_{\kappa+1:n}\\ c'A_{1:\kappa} & c'c & c'A_{\kappa+1:n}\\ A_{\kappa+1:n}'A_{1:\kappa} & A_{1:\kappa}c & A_{\kappa+1:n}'A_{\kappa+1:n} \end{bmatrix}\\ &=LL'\\ &=\begin{bmatrix} L_{11}\\ \ell_{12} & \ell_{22}\\ L_{31} & \ell_{32} & L_{33} \end{bmatrix} \begin{bmatrix} L_{11}\\ \ell_{12} & \ell_{22}\\ L_{31} & \ell_{32} & L_{33} \end{bmatrix}' \end{aligned}

Once we delete the column $c$ at $\kappa+1$ we have

$\begin{bmatrix} A_{1:\kappa}'A_{1:\kappa} & A_{1:\kappa}'A_{\kappa+1:n}\\ A_{\kappa+1:n}'A_{1:\kappa} & A_{\kappa+1:n}'A_{\kappa+1:n} \end{bmatrix}= \begin{bmatrix} \bar{L}_{11}\\ \bar{L}_{31} & \bar{L}_{33} \end{bmatrix} \begin{bmatrix} \bar{L}_{11}\\ \bar{L}_{31} & \bar{L}_{33} \end{bmatrix}'$

from which we have that

\begin{aligned} \bar{L}_{11} &= L_{11}\\ \bar{L}_{31} &= L_{31} \end{aligned}

and

$\bar{L}_{33}\bar{L}_{33}' = L_{33}L_{33}' + \ell_{32}\ell_{32}'.$

The last equation is a rank-1 update (see Golub and Van Loan) and it comes at the
cost of $2(n-\kappa)^2+4(n-\kappa)$.

Appending columns using permutation matrix

Having computed a Cholesky factorisation with permutation for the matrix $A'A$, that is $A'A = PLL'P'$, we need to compute the Cholesky factorisation of $\bar{A}'\bar{A}$, where $\bar{A}=[A\ c]$.

We have

\begin{aligned} \bar{A}'\bar{A} &= \begin{bmatrix}A & c\end{bmatrix}'\begin{bmatrix}A & c\end{bmatrix}\\ &= \begin{bmatrix}A'A & A'c\\c'A & c'c\end{bmatrix}\\ &= \begin{bmatrix}P \\ & 1\end{bmatrix} \begin{bmatrix}L_1\\\ell_1' & \ell_2\end{bmatrix} \begin{bmatrix}L_1\\\ell_1' & \ell_2\end{bmatrix}' \begin{bmatrix}P \\ & 1\end{bmatrix}' \end{aligned}

From which we have that

$A'A = PL_1L_1'P' \Rightarrow L_1 = L,$

and $\ell_1$ is computed from the following nice linear system

$PL\ell_1 = A'c \Leftrightarrow L\ell_1 = P'A'c,$

and, provided that $c'c - \ell_1'\ell_1>0$, $\bar{A}'\bar{A}$ is positive definite and

$\ell_2^2 = c'c - \ell_1'\ell_1.$

Remark. to solve $A'Ax=b$ we do $PLL'P'x= b$ and we set $x=Py$ so we then need to solve

$PLL'y = b \Leftrightarrow LL'y = P'b.$

The permuted Cholesky factor of $\bar{A}'\bar{A}=\bar{P}\bar{L}\bar{L}'\bar{P}'$ is

$\bar{L}=\begin{bmatrix}L_1\\\ell_1' & \ell_2\end{bmatrix},$

and the corresponding permutation matrix is

$\bar{P}=\begin{bmatrix}P \\ & 1\end{bmatrix}.$

The computational cost of this update is $n^2+3n$.

Inserting columns

In this example we shall insert a column in $A$, so we shall define the matrix

$\tilde{A} = \begin{bmatrix} A_{1:\kappa} & c & A_{\kappa+1:n} \end{bmatrix},$

where $\kappa\in\mathbb{N}_{[1,n]}$.

There is then a permutation matrix $\tilde{P}$ so that $\tilde{A}\tilde{P}' = [A\ c]$.

It is then easy to update the factorisation
\begin{aligned} &(\tilde{A}\tilde{P})'(\tilde{A}\tilde{P})=\tilde{L}\tilde{L}', \\ \Leftrightarrow &\ \tilde{A}'\tilde{A} = \tilde{P}\tilde{L}\tilde{L}'\tilde{P}' \end{aligned}

which is the updated Cholesky factorisation of $\tilde{A}'\tilde{A}$ with permutation matrix $\tilde{P}$.

Inserting many columns

Now we are going to insert many columns, recursively, at various positions in $\tilde{A}$. Essentially, we are going to update a permuted Cholesky factorisation by inserting a column in any position in $A$, so the updated matrix becomes

$\bar{A} = \begin{bmatrix} A_{1:\kappa} & c & A_{\kappa+1:n} \end{bmatrix}$

There is a permutation matrix $\tilde{P}$ so that

$\bar{A}\tilde{P} = \begin{bmatrix} A & c \end{bmatrix} \triangleq \tilde{A}$

Say we have $A'A=PLL'P'$. We can then compute a permuted Cholesky factorisation for $\tilde{A}'\tilde{A}$,

\begin{aligned} &\tilde{A}'\tilde{A}= \bar{P}\tilde{L}\tilde{L}'\bar{P}'\\ \Leftrightarrow &\ \tilde{P}'\bar{A}'\bar{A}\tilde{P} = \bar{P}\tilde{L}\tilde{L}'\bar{P}'\\ \Leftrightarrow &\ \bar{A}'\bar{A} = \tilde{P}\bar{P}\tilde{L}\tilde{L}'\bar{P}'\tilde{P}', \end{aligned}
which is the permuted Cholesky factorisation of $\bar{A}'\bar{A}$ with permutation matrix $\tilde{P}\bar{P}$.

Let now $A\in\Re^{m\times n}$ be a given matrix and $\alpha$ be a collection of column indexes of $A$, and $A_\alpha\in\Re^{m\times |\alpha|}$ and we know the Cholesky factorisation of $A_{\alpha}'A_{\alpha}=L_{\alpha}L_{\alpha}'$.

Let $\bar{\alpha} = \alpha\cup \{\alpha^*\}$, i.e., the $\alpha^*$-th column of $A$ is added to $A_{\alpha}$ to form the new matrix

$A_{\bar\alpha} = \begin{bmatrix} A_{\alpha} & A_{\alpha^*} \end{bmatrix}.$

Having augmented $A_{\alpha}$ we can now update the factorisation of $A_{\alpha}'A_{\alpha}$ and compute an $L_{\bar\alpha}$ such that

$A_{\bar\alpha}'A_{\bar\alpha} = L_{\bar\alpha}L_{\bar\alpha}'.$

mathbabe

Exploring and venting about quantitative issues

Look at the corners!

The math blog of Dmitry Ostrovsky

The Unapologetic Mathematician

Mathematics for the interested outsider

Almost Sure

A random mathematical blog

Mathematix

Mathematix is the mathematician of the village of Asterix