search

j2kun / svd

Python code implementing the power method for Singular Value Decomposition
https://jeremykun.com/2016/05/16/singular-value-decomposition-part-2-theorem-proof-algorithm/
103 stars 41 forks source link

Fails on singularity. #4

Open congxinUSC opened 6 years ago

congxinUSC commented 6 years ago

If the matrix is not full-ranked (one at least one of the dimensions, the SVD gives wrong results for the eigenvalues and eigenvectors that corresponds to the eigenvalues should be 0.

j2kun commented 6 years ago

I can confirm that the output of this is nondeterministic.

When computing the third singular value for this matrix [[1,2,3], [4,5,6], [7,8,9]] I get (run a few times, to show it's nondeterministic):

5.00408346
1.87754565
5.07695215
11.89394783
16.00065575
1.68317765

This is due to numerical instability of the part where I subtract the outer product. For instance, when computing the SVD and printing the matrix in each 1d step, on the last one the matrix is reduced to:

[[ -4.44089210e-15  -5.78703752e-15  -6.43929354e-15]
 [ -1.36002321e-15  -1.59247615e-15  -1.85962357e-15]
 [  1.55431223e-15   1.78329573e-15   2.05391260e-15]]

So one fix can be to threshold each entry by epsilon, or threshold the entire matrix based on the max magnitude value. IIRC, if at any step the remaining matrix is zero, then all remaining singular values and vectors are zero.

However, I think this method is inherently not numerically stable. I'm not sure if thresholding is enough to fix the problem. A literature search on numerical stability for the power method might clarify this.