0
$\begingroup$

I need to compute the inverse of an ill-conditioned matrix. Since condition number is ratio of high/low singular values. I am approximating the matrix by removing small singular values. But the conditioned number of obtained matrix is even higher. The python code is:

UU,SS,VV=scipy.linalg.svd(A) # A is 100 x 150 
Sigma = numpy.zeros((70, A.shape[1]))
Sigma[:70, :70] = numpy.diag(SS[0:70])
numpy.linalg.cond(numpy.matmul(numpy.matmul(UU[:,0:70],Sigma),numpy.transpose(VV)))

The condition number of A is 3391639000000000.0 but after removing singular values it becomes 1.712286461590629e+23

$\endgroup$
5
  • 1
    $\begingroup$ Why would one remove the small singular values? The resulting matrix would be singular. The new condition number would be $\sigma_\max/0 = \infty$. By the way, if your initial matrix has such a large condition number, I'm not optimistic that you will be able to invert it. Maybe take a step back and find a way to avoid the necessity of inverting this nearly singular matrix. $\endgroup$
    – littleO
    Commented May 26, 2019 at 11:16
  • $\begingroup$ Why have you put the smallest singular value as 0. The logic I am using is: let max is 1000 and min earlier is 10 : giving ratio 100. But if we remove some smaller singular values: now minimum let is 500 which gives ratio 2 and hence the condition number is reduced. $\endgroup$ Commented May 26, 2019 at 11:20
  • $\begingroup$ Like in my code, I have removed 30 smaller singular values for instance. This number can be decided empirically. Please correct if I my understanding is wrong. $\endgroup$ Commented May 26, 2019 at 11:27
  • $\begingroup$ What you're doing is equivalent to setting all but the top 70 singular values equal to 0. $\endgroup$
    – littleO
    Commented May 26, 2019 at 11:28
  • $\begingroup$ No, in my code I have reduced the size of Sigma (matrix with singular values) as 70 x 150 instead of 100 x 150. $\endgroup$ Commented May 26, 2019 at 11:31

1 Answer 1

1
$\begingroup$

What you're doing is equivalent to setting all but the top 70 singular values equal to $0$, as the code below illustrates. The resulting matrix is singular.

import numpy as np
import scipy as sp

numRows = 100
numCols = 150
numKeep = 70
A = np.random.randn(numRows,numCols)

U,S,VT = sp.linalg.svd(A)  

# This is how we recover A from U, S, and VT
Sigma = np.zeros((numRows,numCols))
for i in range(numRows):
    Sigma[i,i] = S[i]

A_check = U @ Sigma @ VT
print("Is A equal to A_check?")
print(np.allclose(A,A_check))

#Now compute B using method 1, "removing smaller singular values"
Sigma = np.zeros((numKeep, numCols))
for i in range(numKeep):
    Sigma[i,i] = S[i]

B = U[:,0:numKeep] @ Sigma @ VT

# Now compute B using method 2, 
# by setting all but top 70 singular values set equal to 0.
Sigma = np.zeros((numRows,numCols))
for i in range(numKeep):
    Sigma[i,i] = S[i]

B_check = U @ Sigma @ VT
print("Is B equal to B_check?")
print(numpy.allclose(B,B_check))
$\endgroup$

You must log in to answer this question.

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