Import the aligned face dataset and build a list of Image objects from it


In [ ]:
import scipy.io as sio
import numpy as np
from menpo.image import MaskedImage
import menpo.io as pio

im_db = sio.loadmat(pio.data_path_to('alignedbwfaces.mat'))
imagedata = im_db['images']
mask = im_db['mask']
images = []
for i in range(imagedata.shape[-1]):
    images.append(MaskedImage(imagedata[...,i], mask=mask))

images is a Python list of our menpo.image.Image class. We can grab one of them and view it


In [ ]:
%matplotlib inline
image = images[0]
image.view()

the view() method is a quick and easy visualization on all our types. For images, it's roughly equivalent to this


In [ ]:
import matplotlib.pyplot as plt
plt.gray()
plt.imshow(image.pixels[...,0])

but does some legwork to give a sensible result for B+W images too. See the menpo.images.Images.ipynb notebook for more examples of how our Image class works.

Set aside the last 10 images for testing, use the rest for constructing the model


In [ ]:
training_images = images[:-10]
test_images = images[-10:]

from menpo.model import PCAModel
pca = PCAModel(training_images)

Firstly, lets take a look at the mean


In [ ]:
pca.mean.view()

Note that the method above return Image instances!


In [ ]:
type(pca.mean)

This makes working with our models really convienient. Of course statistical models generally only work with vectors of data, so how is this possible? It's all done through the Vectorizable interface (see menpo.vectorizable.ipynb for details). If you ever want to see the low level vectors, just call the mean_vector property.


In [ ]:
pca.mean_vector

What else can we see? How about the proportion of variance captured in each of component


In [ ]:
import matplotlib.pyplot as plt
plt.plot(pca.eigenvalues_ratio)

or cumulative variance captured in the first 40 components


In [ ]:
plt.plot(np.cumsum(pca.eigenvalues_ratio[:40]))

Or just the components themselves


In [ ]:
print 'shape of all the components: {}'.format(pca.components.shape)
print 'pca.components[:5] grabs the first 5 components: {}'.format(pca.components[:5].shape)

Now lets use our model. We can project an in-sample image onto the PCA and check the reconstruction is perfect


In [ ]:
insample_image = training_images[9]
weightings = pca.project(insample_image)
rebuilt_image = pca.instance(weightings)
print 'if the reconstruction is close, an image of the reconstruction is displayed'
if np.allclose(insample_image.pixels, rebuilt_image.pixels):
    rebuilt_image.view()

In [ ]:
print "Examining the numerical error"
plt.imshow((insample_image.pixels - rebuilt_image.pixels)[...,0])

Now lets limit the number of components used in the reconstruction a little...


In [ ]:
rebuilt_image_50_comps = pca.instance(weightings[:50])
rebuilt_image_50_comps.view()

Now lets try an out of sample image


In [ ]:
test_image = test_images[0]
test_image.view()

In [ ]:
weights = pca.project(test_image)
reconstructed_image = pca.instance(weights)
reconstructed_image.view()

This operation of projecting to learn weights followed by reconstruction of an instance is fairly common, so we have some syntactic sugar to make things a little more natural.


In [ ]:
reconstructed_image_alt = pca.reconstruct(test_image)
print "Does the easier reconstruct syntax yield the same result as before?: {}".format(
np.all(reconstructed_image_alt.pixels == reconstructed_image.pixels))

Note that you can also pass in a number of Components to limit the reconstruction power


In [ ]:
reconstructed_image_limited = pca.reconstruct(test_image, n_components=50)
reconstructed_image_limited.view()