Instructions
Requirements and Specifications
Source Code
!pip install otter-grader
# Initialize Otter
import otter
grader = otter.Notebook("lab6.ipynb")
# Lab 6: Image Compression and Matrix Factorization
Matrix factorization is a way to find a set of basis vectors that describe a given dataset. Depending on the factorization used, the set of bases are different.
In this notebook, we use singular value decomposition (SVD) and nonnegative matrix factorization (NMF)
## Faces dataset
We use "[labeled faces in the wild](http://vis-www.cs.umass.edu/lfw/)" dataset.
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
from sklearn.datasets import fetch_lfw_people
lfw_people = fetch_lfw_people(min_faces_per_person=70, resize=0.4)
img_count, img_height, img_width = lfw_people.images.shape
print('number of faces in dataset:', img_count)
print('image width in pixels :', img_width)
print('image height in pixels :', img_height)
Each face is vectorized into a row in data matrix `X`
## Question 1a: Data transformation
Inspecting `lfw_people.images.shape` shows images are stored as a 3-dimensional array of size (1288, 50, 37). Use `numpy.ndarray.reshape()` to tranform matrix `lfw_people.images` to a 2-dimensional array of size `image_count` by `image_width` * `image_height` with name `X`. Take first image, `X[0]`, and `numpy.ndarray.reshape()` it back to a 2-dimensional array of size `image_width` * `image_height` with name `X0_img`.
X = lfw_people.images.reshape((img_count, img_width*img_height))
X0_img = X[0].reshape((img_height, img_width))
grader.check("q1a")
If everything went correctly, you should see a gray scale image below
# each row of X is a vectorized image
plt.imshow(X0_img, cmap=plt.cm.gray);
## Question 1b: Visualization
To make plotting easier, create a plotting function.
def draw_img(img_vector, h=img_height, w=img_width):
"""
1. takes img_vector,
2. reshapes into right dimensions,
3. draws the resulting image
"""
plt.imshow( img_vector.reshape((img_height, img_width)), cmap=plt.cm.gray)
plt.xticks(())
plt.yticks(())
_Cell Intentionally Blank_
# check draw_img function
draw_img(X[49])
## Question 1c: Standardization
Since SVD looks for singular values (related to eigenvalues) and eigenvectors of `X`, center each column (pixel) of `X` so that mean of each column is zero. Otherwise the result can be strange. This is a detail that is not critical for understanding the conceptual aspect of SVD. The variance of each pixel can be left alone. Use [`sklearn.preprocessing.StandardScaler`](https://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.StandardScaler.html) to standardize `X` (hint: use `.fit_transform`).
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler(with_std=False)
Xstd = scaler.fit_transform(X)
grader.check("q1c")
The mean centered images look unnatural, but all the pertinent data is retained
# standardization transforms image data
draw_img(Xstd[49])
### Inverse Transformation
We can recover the original data by putting the means back
# inverse tranformation recovers original image units
Xorig = scaler.inverse_transform(Xstd)
draw_img(Xorig[49])
## Question 2a: Singular Value Decomposition (SVD)
Numpy package has SVD decomposition function [`numpy.linalg.svd`](https://numpy.org/doc/stable/reference/generated/numpy.linalg.svd.html). Decompose `Xstd` into `U`, `S`, `VT`, i.e., $X_{std} = U S V^T$.
from numpy.linalg import svd
U, S, VT = np.linalg.svd(Xstd)
grader.check("q2a")
## Question 2b: Matrix Factors from SVD
Since `S` is a diagonal matrix, multiplying `U` by `S` is like scaling each column of `U` by the corresponding diagonal entry. Check that we can recover the original data from the matrix factors by inverse transforming the resulting `Xhat`
# Compute real S
Strue = np.zeros((U.shape[1], VT.shape[0]))
Strue[:S.size, :S.size] = np.diag(S)
US = U.dot(Strue)
# reconstruct standardized images from matrix factors
Xhat = US.dot(VT)
# inverse transform Xhat to remove standardization
Xhat_orig = scaler.inverse_transform(Xhat).astype('float32')
draw_img(Xhat_orig[49])
grader.check("q2b")
### Dimensionality reduction
We can describe each face using smaller portions of matrix factors. Because of how `US` and `VT` are ordered, the first portions of `US` and `VT` retain the most information. So, to keep the most relevant parts, keep the first columns of `US` and first rows of `VT`. Below illustrate why this is called a dimensionality reduction method.
In the following, we keep 500 columns of `US` and 500 rows of `VT` out of 1288.
# reconstruct Xhat with less information: i.e. dimensionality is reduced
Xhat_500 = US[:, 0:500] @ VT[0:500, :]
# inverse transforms Xhat to remove standardization
Xhat_500_orig = scaler.inverse_transform(Xhat_500)
# draw recovered image
draw_img(Xhat_500_orig[49])
Using even smaller number degrades the reconstruction; however, still the reconstructed image captures the "gist" of the original data.
# reconstruct Xhat with less information: i.e. dimensionality is reduced
Xhat_100 = US[:, 0:100] @ VT[0:100, :]
# inverse transforms Xhat to remove standardization
Xhat_100_orig = scaler.inverse_transform(Xhat_100)
# draw recovered image
draw_img(Xhat_100_orig[49])
To make the dimensionality reduction and inverse transform easier, write a function:
def dim_reduce(US_, VT_, dim=100):
Xhat_ = US_[:, 0:dim] @ VT_[0:dim, :]
return scaler.inverse_transform(Xhat_)
We can see how the increasing the rows and columns of matrix factors used increases fidelity of reconstructions
dim_vec = [50, 100, 200, 400, 800]
plt.figure(figsize=(1.8 * len(dim_vec), 2.4))
for i, d in enumerate(dim_vec):
plt.subplot(1, len(dim_vec), i + 1)
draw_img(dim_reduce(US, VT, d)[49])
### Matrix Factors and "Eigenfaces"
What are in these matrix factors `US` and `VT`?
`VT` contains the set of "basis" vectors. In this setting rows of `VT` are called eigenfaces
# each row of VT is an "eigenface"
draw_img(VT[0])
Plotting more eigenfaces show how the information in each eigenfaces are highlighting different components of photos
num_faces = 10
# each row of VT is an "eigenface"
plt.figure(figsize=(1.8 * 5, 2.4 * 2))
for i in range(0, 10):
one_face = VT[i]
plt.subplot(2, 5, i + 1)
draw_img(one_face)
### Each face is a linear combination of eigenfaces
Reconstructing a face can be thought of as combining these eigenfaces according to some row of `US` (coefficients).
# each face is a linear combination of eigenfaces VT
face_num = 3 # which face to reconstruct?
dim = 300 # higher dim is more accurate fit
draw_img(scaler.inverse_transform(US[face_num, 0:dim].reshape((1,-1)) @ VT[:dim,:]))
## Question 3: Nonnegative Matrix Factorization
NMF is a matrix factorization method that require nonnegative data matrix. Images are represented as light intensities between 0 and 255: i.e. nonnegative numbers.
NMF decomposes `X` into `W` and `H` such that $X \approx WH$. NMF is slower than SVD. So, we only choose a very small number of basis here: 200. Obtain `W` and `H`. Use [`sklearn.decomposition.NMF`](https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.NMF.html).
from sklearn.decomposition import NMF
model = NMF(n_components=200, init='nndsvd', random_state=0)
W = model.fit_transform(X)
H = model.components_
grader.check("q3")
### Matrix factor H
NMF matrix factor `H` contain the set of basis faces.
draw_img(H[0])
Following shows that the set of basis vectors are very different than what SVD chose
num_faces = 20
plt.figure(figsize=(1.8 * 5, 2.4 * 4))
for i in range(0, num_faces):
one_face = VT[i]
plt.subplot(4, 5, i + 1)
draw_img(H[i])
However, each face is still a linear combination of matrix `H`
draw_img(W[30]@H) # using 200 NMF basis vectors
draw_img(dim_reduce(US, VT, 200)[30]) # using 200 SVD basis vectors
draw_img(X[30]) # original image
---
To double-check your work, the cell below will rerun all of the autograder tests.
grader.check_all()
## Submission
Make sure you have run all cells in your notebook in order before running the cell below, so that all images/graphs appear in the output. The cell below will generate a zip file for you to submit. **Please save before exporting!**
# Save your notebook first, then run this cell to export your submission.
grader.export()