Precision of SVD in solving Least Square goes against theory

137 Views Asked by At

Let $A\in\mathbb{R}^{100\times100}$ (with $\operatorname{rank}(A)=100$) and $b\in\mathbb{R}^{100}$, the goal is to solve for $x$ :

$$ \|b-Ax\|_{2}=\min_{y\in\mathbb{R}^{n}}\|b-Ay\|_{2} $$using SVD.

SVD Approach

The SVD is known to be the most stable in finding least-squares as compared to QR and Cholesky. The approach goes by defining the economic version of SVD that is $A=U_{r}S_{r}V_{r}^{T}$ where $U_{r},S_{r},V_{r}\in\mathbb{R}^{100\times 100}$, then the least-square problem can be reduced to solving the following normal system : $$ S_{r}V_{r}^{T}x=U_{r}^{T}b $$ and this can be done as follow :

$1.$ Compute $c=U_{r}^{T}b$

$2.$ Find $w$ in $S_{r}w=c$.

$3.$ Find $x$ in $V_{r}^{T}x=w$.

Experimenting goes Against Theory

After writing a MATLAB algorithm to compute the least-square of $A$ and $b$ with both $A$ and $b$ being randomly generated, then executing it, I received the vector $x\in\mathbb{R}^{100}$. However the issue began to arise when I computed the norm $\|b-Ax\|_{2}$ to determine the stability and shockingly it appeared to be less precise ($\sim 6\times 10^{-16}$) than when I computed the least-square of the same $A$ and $b$ using QR which turned out to be $\sim 3\times10^{-16}$ which is about twice as precise as SVD approach which made me wonder why this happened.

The Question

How could this be interpreted knowing the following information I gathered:

  • $\operatorname{cond(A)}\sim 10^{3}$
  • $\|A\|_{2}=\sigma_{\max}(A)\approx 100$
  • The least three singular values of $A$ where $\sigma_{1}=0.1885, \sigma_{2}=0.0580,$ and $\sigma_{3}=0.0128$.

I was wondering if singular values has to do anything with the precision whether the scale of their values can affect the accuracy of $x$ since it was supposed to be that SVD is most stable in solving least-squares problems.

Side Note : The algorithm that I wrote involves using MATLAB's built-in SVD function using its economic version. Therefore, everything in my algorithm is optimized.