Why the PCA image doesnt resemble the original image at all?

I am trying to implement PCA without any library for image dimension reduction. I tried the code in the O'Reilly Computer Vision book and implement it on a sample lenna picture:

    from PIL import Image
    from numpy import *

    def pca(X):
        num_data, dim = X.shape

        mean_X = X.mean(axis=0)
        X = X - mean_X

        if dim > num_data:
            # PCA compact trick
            M = np.dot(X, X.T) # covariance matrix
            e, U = np.linalg.eigh(M) # calculate eigenvalues an deigenvectors
            tmp = np.dot(X.T, U).T
            V = tmp[::-1] # reverse since the last eigenvectors are the ones we want
            S = np.sqrt(e)[::-1] #reverse since the last eigenvalues are in increasing order
            for i in range(V.shape[1]):
                V[:,i] /= S
        else:
            # normal PCA, SVD method
            U,S,V = np.linalg.svd(X)
            V = V[:num_data] # only makes sense to return the first num_data
        return V, S, mean_X
img=color.rgb2gray(io.imread('D:\lenna.png'))
x,y,z=pca(img)
plt.imshow(x)

but the image plot of the pca doesnt look like the original image like at all. As far as i know PCA kinda reduce the image dimension but it will still somehow resemble the original image but in lower detail. Whats wrong with the code?

Well, nothing is wrong per se in your code, but you're not displaying the right thing if I do understand what you actually want to do!

What I would write for your problem is the following:

def pca(X, number_of_pcs):
    num_data, dim = X.shape

    mean_X = X.mean(axis=0)
    X = X - mean_X

    if dim > num_data:
        # PCA compact trick
        M = np.dot(X, X.T) # covariance matrix
        e, U = np.linalg.eigh(M) # calculate eigenvalues an deigenvectors
        tmp = np.dot(X.T, U).T
        V = tmp[::-1] # reverse since the last eigenvectors are the ones we want
        S = np.sqrt(e)[::-1] #reverse since the last eigenvalues are in increasing order
        for i in range(V.shape[1]):
            V[:,i] /= S

        return V, S, mean_X

    else:
        # normal PCA, SVD method
        U, S, V = np.linalg.svd(X, full_matrices=False)

        # reconstruct the image using U, S and V
        # otherwise you're just outputting the eigenvectors of X*X^T
        V = V.T
        S = np.diag(S)
        X_hat = np.dot(U[:, :number_of_pcs], np.dot(S[:number_of_pcs, :number_of_pcs], V[:,:number_of_pcs].T))      

        return X_hat, S, mean_X

The change here lies in the fact that we want to reconstruct the image using a given number of eigenvectors (determined by number_of_pcs).

The thing to remember is that in np.linalg.svd, the columns of U are the eigenvectors of X.X^T.

When doing that, we obtain the following results (displayed here using 1 and 10 principal components):

---

X_hat, S, mean_X = pca(img, 1)
plt.imshow(X_hat)

---

X_hat, S, mean_X = pca(img, 10)
plt.imshow(X_hat)

PS: note that the picture aren't displayed in grayscale because of matplotlib.pyplot, but this is a very minor issue here.