2

I have difficulties to handle singular values close to zero. My SVD implementation $A = U \Sigma V {}^{T}$ performs first an Eigen decomposition of the matrix $A {}^{T} A$. That is done with the QR-Algorithm, where the QR decomposition is computed with Householder reflections.

From the resulting Eigenvalues $\lambda_i$ and Eigenvectors $v_i$, I compute the singular values $$ \Sigma {}{_{i,i}} = \sigma {}_i= \sqrt{\lambda_i} $$ the right singular vectors $$ V_i = v_i $$ and finally the left singular vectors $$ U_i = \frac{Av_i}{\sigma {}_i}. $$

In case $\sigma_i = 0$, the corresponding left singular vector is determined as discussed here. As I work with double-percision, I do get small numbers instead of zero, which lead to an incorrect $U_i$. The situation is actually worse, as the QR algorithm sometimes computes even very small negative Eigenvalues.

So far, I use a static thresholds to decide what is zero and also checks, such as testing in the end if $U$ is orthogonal. This appeared to be a lousy approach, because I have to adapt it regularly.

Does there exists a good concept or strategy to overcome these issues? Should the problem be tackled in the Eigen decomposition?

If it of interest: The whole code is here and the latest SVD problem arose when implementing this paper.

Adrian Schneider
  • 153
  • 3
    Mathematically, it is true that one can work with $A^T A$ to get the SVD, but numerically this is a bad idea. There are better algorithms for the SVD that don't construct $A^T A$ or $A A^T$ but instead work with $A$ directly. You should look into them. – Ian Aug 30 '18 at 12:41
  • The method proposed by Golub and Kahan seems to be numerically more stable, even though it does not work on A directly. Would that be a good alternative or do you have a specific method in mind? Thanks for your help. – Adrian Schneider Aug 31 '18 at 08:14
  • Did you end up solving this issue? I'd think that perhaps iterative methods would not run into these issues. – Nathaniel Bubis Jan 22 '21 at 07:37
  • 1
    e.g. here: https://math.stackexchange.com/a/3990583/28743 – Nathaniel Bubis Jan 22 '21 at 08:23
  • Thank you for the hint. I implemented the Golub and Kahan approach. The floating-point zero problem is there explicitly stated. From an engineer's perspective, it is well explained in the book "Matrix Computations, Gene H. Golub". – Adrian Schneider Jan 26 '21 at 07:31

0 Answers0