# Pseudoinverse of perturbed matrix

How does the pseudo inverse of a full column rank matrix change if I rescale a single row?

In more detail the problem is the following:

We have a fixed matrix $V$ with linear independent columns and lots of matrices $D_i$ of the form

$$D_i = \begin{pmatrix} d_i & 0 \\ 0 & E_n \end{pmatrix}$$ where $E_n$ is the identity matrix. All matrices are real and have nice condition numbers. So far so sweet.

We know the pseudoinverse $V^+=(V^*V)^{-1} V^*$ of $V$ and the job is to compute all the pseudo inverses

$$(D_i \cdot V)^+$$

in a fast way. Any ideas how to exploit the structure of the problem here?

For example if the problem was the other way round, algebra of pseudo inveses would allow to break it down to inverting a single number and multiplying a single row:

$$(V \cdot D_i )^+ = D_i^{-1} \cdot V^+$$

In my applications $V$ is of small size, like $4 \times 3$ for example.

• How big is the size of $E_n$ compared to $[d_i]$ ? – Tolga Birdal Feb 5 '16 at 11:18
• For me $d_i$ a single number and $E_n$ a $3 \times 3$ matrix is the most important case. – Jan Weidner Feb 6 '16 at 9:19

So, since $V\in\mathbb R^{m\times n}$ has full column rank $n$, implying also $m\geq n$, we might obtain the pseudoinverse using $$V^+=(V^TV)^{-1}V^T$$ since $V^TV$ is invertible and $V^*=V^T$ in $\mathbb R$.

Now, let's denote by $i$ the index of the row which is scaled and by $d_i$ - the scaling factor. Thus, a perturbed matrix $V^\prime$ can be expressed as $$V^\prime = D_i V= \text{diag}(\underbrace{1,\ldots,1}_{i-1},d_i,\underbrace{1,\ldots,1}_{m-i})V$$

Now, how can we express it as a low-rank update?

Consider $$u=e_i=\left[ \begin{array}{c} 0\\ \vdots\\ 1 \\ \vdots\\ 0 \end{array}\right], \quad\quad v=(V_{i,*})^T$$

where $u$ is the $i$th elementary unit vector in $m$ dimensional space and $v$ is the transposed $i$th row of the original matrix $V$. Then: $$V^\prime = V+(d_i-1)uv^T$$

Now, let's express the update to the "correlation" matrix $V^T V$: \begin{aligned} {V^\prime}^T{V^\prime}&=(V+(d_i-1)uv^T)^T(V+(d_i-1)uv^T)\\ &=V^TV+(d_i-1)v\underbrace{u^TV}_{v^T}+(d_i-1)\underbrace{V^Tu}_{v}v^T+(d_i-1)(d_i-1)v\underbrace{u^Tu}_{1}v^T\\ &=V^TV+\left((d_i-1)+(d_i-1)+(d_i-1)^2\right)vv^T\\ &=V^TV+(d_i^2-1)vv^T \end{aligned}

With that, we expressed a rank-1 update to the correlation matrix $V^T V$ for which we already know the inverse. For compactness, let's denote

$$\gamma=d_i^2-1$$.

Now, let's use Sherman-Morrison to find an update to the inverse of the correlation matrix:

\begin{aligned} ({V^\prime}^T{V^\prime})^{-1}&=(V^TV+\gamma v v^T)^{-1}\\ &=(V^TV)^{-1}-\frac{\gamma (V^TV)^{-1}vv^T(V^TV)^{-1}}{1+\gamma v^T (V^TV)^{-1}v} \end{aligned}

Let's denote the following matrix-vector product by $$w = (V^TV)^{-1}v$$ Then, since $V^T V$ is symmetric

\begin{aligned} ({V^\prime}^T{V^\prime})^{-1}&=(V^TV)^{-1}-\frac{\gamma ww^T}{1+\gamma v^T w}\\ &=(V^TV)^{-1}-\tau ww^T, \quad \tau=\frac{d_i^2-1}{1+(d_i^2-1)v^Tw} \end{aligned}

So, in order to update the inverse of the correlation matrix (provided you have the computation for the original matrix done):

1. Compute matrix-vector product $w = (V^TV)^{-1}v$ in $n^2$ operations.
2. Calculate dot-product of $w$ and original row $v$: $v^T w$ in $n$ operations.
3. Calculate $\tau$ in constant time.
4. Calculate the scaled by $\tau$ outer product $\tau ww^T$ in $n^2$ operations
5. Calculate the update to the inverse of the correlation matrix by subtracting the outer product in $n^2$ operations.

So, you achieved $({V^\prime}^T{V^\prime})^{-1}$ in $3n^2+n$ operations. Probably, can be done in $2n^2+n$ if fused flops are used combining steps 4 and 5.

Now, you can calculate

$${V^\prime}^+ = ({V^\prime}^T{V^\prime})^{-1}V^TD_i$$

or leave it in the unassembled form depending on your future usage.

I am not sure if that will lead to significant speed-up for the requested $\mathbb R^{4\times 3}$ case but is worth trying.

Note, this particular technique relies on $V$ having full-column rank, otherwise, the pseudoinverse cannot be found this way and SVD decomposition has to be used -> where updates are also possible, but much trickier.