8
$\begingroup$

Is there formula for the inverse of a matrix which is diagonal plus a sum of rank one matrices?

$$S=\alpha I + \sum_1^N u_iu_i^T$$

$$S^{-1} = ?$$

Is there any decomposition or special trick that I can solve this inverse faster?

$\endgroup$
5
  • 1
    $\begingroup$ Any matrix can be written in this form for $N$ high enough. $\endgroup$
    – tom
    Sep 18, 2013 at 19:20
  • $\begingroup$ well, symmetric matrix $\endgroup$
    – Will Jagy
    Sep 18, 2013 at 20:10
  • $\begingroup$ Jing, what happens in the 2 by 2 and 3 by 3 cases? $\endgroup$
    – Will Jagy
    Sep 18, 2013 at 20:16
  • $\begingroup$ cause I am gonna do this at large-scale, say N = 10 billion. I am wondering how can I efficiently invert this matrix $\endgroup$
    – Jing
    Sep 19, 2013 at 3:52
  • $\begingroup$ Miller (1981) $\endgroup$
    – DJJ
    Apr 15, 2015 at 21:55

1 Answer 1

14
$\begingroup$

Use the formula of Mr Woodbury. You can express the sum of the rank one matrices as $UU^T$ where $U=[u_1,...,u_N]$.

EDIT: If $$ S_k = \alpha I + \sum_{i=1}^ku_iu_i^T, $$ you can use it as well (or more precisely the Sherman-Morrison formula) to update the inverse of $S_k$ from that of $S_{k-1}$ (because $S_k=S_{k-1}+u_ku_k^T$): $$ S_k^{-1} = S_{k-1}^{-1} - \frac{w_kw_k^T}{1+u_k^Tw_k}, \quad w_k=S_{k-1}^{-1}u_k. $$ If $\alpha>0$, then there should not be any problem with the denominator since $S_k$ would be positive definite for all $k$.

EDIT: Another solution: It would probably be better to instead of updating the inverse of $S$ to update its factorization. If $\alpha>0$, then $S$ (or any $S_k$) is SPD so that $S$ can be factorized using the Cholesky factorization $S=LL^T$, where $L$ is lower triangular. Now assume that $\tilde{S}=S+uu^T$. You can write $\tilde{S}$ in the form $$ \tilde{S}=\begin{bmatrix}L^T \\ u^T\end{bmatrix}^T\begin{bmatrix}L^T \\ u^T\end{bmatrix}. $$ This matrix $[L,u]^T$ is "almost" upper triangular except the last row. You can apply the product of Givens rotations $Q=Q_1\cdots Q_n$ to zero out the nonzero entries in the row containing $u$ so that you obtain a QR factorization $$ Q^T\begin{bmatrix}L^T\\u^T\end{bmatrix}=\begin{bmatrix}R\\0\end{bmatrix}. $$ You can check that $$ \tilde{S}=\begin{bmatrix}L^T \\ u^T\end{bmatrix}^T\begin{bmatrix}L^T \\ u^T\end{bmatrix}=R^TR, $$ that is, $R^T$ is the Cholesky factor after the rank one update. Note that you need only to update the matrix $L$ and you actually do not need to form the orthogonal matrix $Q$ (just to apply the Givens rotations on the rows of $L^T$ (or equivalently, the columns of $L$) and $u$).

This approach could look like an overkill. It of course depends on the dimensions and the number of vectors. Anyway, this is numerically more stable than updating the inverse using the Sherman-Morrison formula.

$\endgroup$
4
  • $\begingroup$ Thank you, but I am gonna do this at large-scale say N = 10 billion. I can't use this. I think I am looking for some kind of sequential formula. Something like $S^{-1}_k = S^{-1}_{k-1} + ....$$ $\endgroup$
    – Jing
    Sep 19, 2013 at 3:54
  • 2
    $\begingroup$ Well, that was not specified. Of course the choice depends on the dimension of $S$ and the number of vectors $N$. $\endgroup$ Sep 19, 2013 at 11:47
  • $\begingroup$ Your comments and answer is very helpful. My problem has something like 1 billion by 1 million matrix. It sparse, with density about 0.00001. I think your first suggestion fits my need more, cuz rank one update in sparse matrix can be done very fast. cholesky seems more risky for me. Thank you $\endgroup$
    – Jing
    Sep 19, 2013 at 14:52
  • 1
    $\begingroup$ Actually if the vectors $u$ are sparse, then the approach based on updating the Cholesky factorization using Givens rotations could be more feasible. Usually the inverse of a sparse matrix is dense. On the other hand, you could maintain some sparsity of its (Cholesky) factors and Givens rotations are normally used in sparse QR factorizations. $\endgroup$ Sep 19, 2013 at 15:31

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .