15 Pre-Class Assignment: Singular Value Decomposition
Contents
15 Pre-Class Assignment: Singular Value Decomposition#
Goals for today’s pre-class assignment#
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
sym.init_printing()
1. Video On Signular Value Decomposition#
Calculating the eigenvaluse and eigenvectors of a matrix requires that the matrix is square (\(n \times n\)). Think of SVD as a way to get similar information from a non-square matrix.
This video is a little longer than most videos. I recommend reading the questions first and then watching the video and following along with the provided code examples.
from IPython.display import YouTubeVideo
YouTubeVideo("EfZsEFhHcNM",width=640,height=360)
2. Diagonalizing a Square Matrix#
As we saw in previous assignments, a square matrix \(n\times n\), whose eigenvectors form a basis of \(\mathbb R^n\) is diagonalizable. Consider the following matrix provided in the above video (~1:50 into the video)
A = np.matrix([[1,1,3], [-3, -5, -3], [3, 3, 1]])
sym.Matrix(A)
✅ DO THIS: Calculate the Eigenvalues and Eigenvectors for the above matrix \(A\). name the eigenvalues (vals) and the eigenvectors (vecs)
##Put your code here
The following code sorts the eigenvectors based on the ordering of the eigenvalues. It does this using the argsort algorithm which puts the values in ascending order. The [::-1] notation reverses the order of the indexes. We will put the sorted vectors into a matrix called \(P\)
idx = vals.argsort()[::-1]
vals = vals[idx]
vecs = vecs[:,idx]
print(idx)
V = vecs
sym.Matrix(V)
print(V)
P = vecs
sym.Matrix(P)
The following code generates a matrix of the same size as \(A\) and puts the eigenvalues on the diagonals.
D = np.zeros(A.shape)
for i in range(len(vals)):
D[i,i] = vals[i]
sym.Matrix(D)
✅ DO THIS: Show that \(A=PDP^{-1}\)
#Put your answer here
###Answer###
P*D*np.linalg.inv(P)-A
###Answer###
✅ DO THIS: Show that \(A^{10}=PD^{10}P^{-1}\). I.e. multiply A by itself 10 times and then taking \(\sigma^{10}\) for each of the eigenvalues in the diagonal of \(D\).
Put your answer to the above question here
✅ QUESTION 1: What is the estimated complexity of multiplying an arbitrary \(n\times n\) matrix by itself \(m\) times? (Hint: both \(m\) and \(n\) should be included in your Big-O notation answer).
Put your answer to the above question here
Note:#
Diagonalization of an \(n\times n\) matrix takes \(O(n^3)\) operations.
✅ QUESTION 2: What is the estimated complexity (Big-O) of multiplying a diagonlized \(n \times n\) matrix by itself \(m\) times? i.e. you already have \(A=PDP^{-1}\) calculated what does it cost to calculate \(A^m\)?
Put your answer to the above question here
Inverse of an orthonormal matrix \(U\) is \(U^\top\)#
Caluclate the \(P\) matrix for the following matrix \(A\) such that \(A=PDP^{-1}\)
A = np.matrix([[6,-2,-1],[-2,6,-1], [-1,-1,5]])
sym.Matrix(A)
✅ DO THIS: Check the columns of \(P\). Are they orthogonal?
Put your answer here
✅ DO THIS: Show that in this special case (Where \(A\) is symmetric) \(A=PDP^\top\) with a special \(P\).
Put your answer here
✅ QUESTION 3: What do you need to do to invert an arbitrary \(n \times n\) matrix?
Put your answer here.
✅ QUESTION 4: What do you need to do to invert a \(n \times n\) orthonormal matrix?
Put your answer here.
3. Singular Value Decomposition (SVD)#
SVD is not restricted to diagonalizable matrices, for example it can be applied to square non-diagonalizable matrices, as well as to non-square matrices.
Now consider the non-square \(n \times m\) matrix \(A\)
A = np.matrix([[4, 11, 14], [8, 7, -2]])
sym.Matrix(A)
The following code calculates \(A^\top A=VDV^\top\):
sym.Matrix(A.T*A)
eigvals, eigvec = np.linalg.eig(A.T*A)
idx = eigvals.argsort()[::-1]
eigvals = eigvals[idx]
eigvec = eigvec[:,idx]
V = eigvec
sym.Matrix(V)
✅ DO THIS: Calculate \(AA^\top=UDU^\top\):
#Put your answer here
The following code calculates \(\Sigma\) by putting the singular values on the diagonal of an \(n \times m\) zero matrix:
np.sqrt(np.abs(eigvals))*np.identity(1)
E = np.zeros(A.shape)
for i in range(len(eigvals)):
E[i,i] = np.sqrt(eigvals[i])
sym.Matrix(E)
✅ DO THIS: Verify that \(A=U \Sigma V^\top\) or \(A=-U\Sigma V^\top\).
Note: \(\sigma_i^2=\lambda_i\) is an eigenvalue of \(A^TA\) and also \(AA^T\). When we put the singular values in descending order, \(\sigma_1\geq \sigma_2 \geq \dots \sigma_r>0\), they are uniquely determined.
#Put your answer here