# coding: utf-8 # # Pseudoinverse and Least Squares # In[1]: import numpy as np import numpy.linalg as la np.set_printoptions(precision=4, linewidth=100) # In[2]: A = np.random.randn(5, 3) # Now compute the SVD of \$A\$. Note that `numpy.linalg.svd` returns \$V^T\$: # In[3]: U, singval, VT = la.svd(A) V = VT.T # Let's first understand the shapes of these arrays: # In[4]: print(U.shape) print(singval.shape) print(V.shape) # Check the orthogonality of \$U\$ and \$V\$: # In[5]: U.T.dot(U) # In[6]: V.T.dot(V) # Now build the matrix \$\Sigma\$: # In[7]: Sigma = np.zeros(A.shape) Sigma[:3, :3] = np.diag(singval) Sigma # Now piece \$A\$ back together from the factorization: # In[8]: U.dot(Sigma).dot(V.T) - A # ---------------- # Next, compute the pseudoinverse: # In[9]: SigmaInv = np.zeros((3,5)) SigmaInv[:3, :3] = np.diag(1/singval) SigmaInv # In[10]: A_pinv = V.dot(SigmaInv).dot(U.T) # ------------- # Now use the pseudoinverse to "solve" \$Ax=b\$ for our tall-and-skinny \$A\$: # In[11]: b = np.random.randn(5) # In[12]: A_pinv.dot(b) # --------------- # Compare with the solution from QR-based Least Squares: # In[13]: Q, R = la.qr(A, "complete") la.solve(R[:3], Q.T.dot(b)[:3])