diff --git a/RBM/UFLDL/pca/pca_gen.m b/RBM/UFLDL/pca/pca_gen.m
index f528a56..443bedd 100644
--- a/RBM/UFLDL/pca/pca_gen.m
+++ b/RBM/UFLDL/pca/pca_gen.m
@@ -113,7 +113,7 @@ xPCAWhite = diag(sqrt(1./(diag(S) + epsilon))) * xrot;
 %  becoming smaller.
 
 %%% YOUR CODE HERE %%%
-xpca_nonreg = diag(sqrt(1./(diag(S) + 1e-6))) * xrot;
+xpca_nonreg = diag(sqrt(1./(diag(S)))) * xrot;
 sigma = xpca_nonreg * xpca_nonreg' / size(xpca_nonreg, 2);
 
 % Visualise the covariance matrix. You should see a red line across the
diff --git a/RBM/mnist_gen_pca.m b/RBM/mnist_gen_pca.m
new file mode 100644
index 0000000..0964713
--- /dev/null
+++ b/RBM/mnist_gen_pca.m
@@ -0,0 +1,64 @@
+function mnist_gen_pca(retained_variance_target)
+
+close all;
+
+addpath(genpath('UFLDL/common'))
+fprintf('Load mnist raw data\n');
+
+x = loadMNISTImages('UFLDL/common/train-images-idx3-ubyte');
+
+rbmWrite(x, 'mnist.trainingStates.dat')
+figure('name','Raw images');
+randsel = randi(size(x,2),64,1); % A random selection of samples for visualization
+display_network(x(:,randsel));
+
+fprintf('Zero mean raw data\n');
+avg = mean(x, 1);     % Compute the mean pixel intensity value separately for each patch. 
+x = x - repmat(avg, size(x, 1), 1);
+
+fprintf('Do the PCA whitening\n');
+[xHat, k, xZCAWhite] = pca(x, retained_variance_target);
+
+figure('name',['PCA processed images ',sprintf('(%d / %d dimensions)', k, size(x, 1)),'']);
+display_network(xHat(:,randsel));
+
+fprintf('Normalize the whitened data\n');
+minZ = min(xZCAWhite(:))
+maxZ = max(xZCAWhite(:))
+xZCAWhite_norm = (xZCAWhite - minZ)/(maxZ-minZ);
+
+minZ = min(xHat(:))
+maxZ = max(xHat(:))
+xHat_norm = (xHat - minZ)/(maxZ-minZ);
+
+figure('name','ZCA whitened images');
+display_network(xZCAWhite_norm(:,randsel));
+
+figure('name','xHat');
+display_network(xHat_norm(:,randsel));
+
+fprintf('Save data\n');
+rbmWrite(xZCAWhite_norm, 'mnist_zca.trainingStates.dat')
+rbmWrite(xHat_norm, 'mnist_xhat.trainingStates.dat')
+
+
+function rbmWrite(data, name)
+[numPixel, numTrain] = size(data)
+nx = sqrt(numPixel)
+ny = nx
+
+data = reshape(data, nx, ny, numTrain);
+
+fid = fopen(name, 'w');
+fprintf(fid, '%d\n', numTrain);
+fprintf(fid, '%d\n', numPixel);
+
+for m=1:numTrain,
+	d = data(:,:,m)';
+	d = d(:);
+	for n=1:numPixel,
+		fprintf(fid, '%f\n', d(n));
+	end
+end
+
+fclose(fid);
diff --git a/RBM/pca.m b/RBM/pca.m
new file mode 100644
index 0000000..b29ea2b
--- /dev/null
+++ b/RBM/pca.m
@@ -0,0 +1,75 @@
+%%================================================================
+%% Step 0a: Load data
+%  Here we provide the code to load natural image data into x.
+%  x will be a 784 * 600000 matrix, where the kth column x(:, k) corresponds to
+%  the raw image data from the kth 12x12 image patch sampled.
+%  You do not need to change the code below.
+
+function [xHat, k, xZCAWhite] = pca(x, retained_variance_target)
+ 
+%%================================================================
+%% Step 1a: Implement PCA to obtain xRot
+%  Implement PCA to obtain xRot, the matrix in which the data is expressed
+%  with respect to the eigenbasis of sigma, which is the matrix U.
+
+fprintf('Calculate covariance matrix of size %d x %d\n', size(x, 1), size(x, 1));
+sigma = x * x' / size(x, 2);
+
+fprintf('Perform Singular Value Decomposition\n');
+[U,S,V] = svd(sigma);
+
+fprintf('Perform PCA transformation\n');
+xrot = U' * x;
+
+%%================================================================
+%% Step 2: Find k, the number of components to retain
+%  Write code to determine k, the number of components to retain in order
+%  to retain at least 99% of the variance.
+
+fprintf('Find K for variance target of %f\n', 100*retained_variance_target);
+retained_variance = 1.0;
+lambda = diag(S);
+k = length(lambda);
+while(retained_variance > retained_variance_target)
+	retained_variance = sum(lambda(1:k))./sum(lambda);
+	k = k - 1;
+end
+k = k
+
+%%================================================================
+%% Step 3: Implement PCA with dimension reduction
+%  Now that you have found k, you can reduce the dimension of the data by
+%  discarding the remaining dimensions. In this way, you can represent the
+%  data in k dimensions instead of the original 144, which will save you
+%  computational time when running learning algorithms on the reduced
+%  representation.
+% 
+%  Following the dimension reduction, invert the PCA transformation to produce 
+%  the matrix xHat, the dimension-reduced data with respect to the original basis.
+%  Visualise the data and compare it to the raw data. You will observe that
+%  there is little loss due to throwing away the principal components that
+%  correspond to dimensions with low variation.
+
+fprintf('Calculate xHat\n');
+xHat = U(:, 1:k) * xrot(1:k, :);
+
+%%================================================================
+%% Step 4a: Implement PCA with whitening and regularisation
+%  Implement PCA with whitening and regularisation to produce the matrix
+%  xPCAWhite. 
+
+epsilon = 1e-1; 
+%%% YOUR CODE HERE %%%
+fprintf('Calculate xPCAWhite\n');
+xPCAWhite = diag(sqrt(1./(diag(S) + epsilon))) * xrot;
+
+%%================================================================
+%% Step 5: Implement ZCA whitening
+%  Now implement ZCA whitening to produce the matrix xZCAWhite. 
+%  Visualise the data and compare it to the raw data. You should observe
+%  that whitening results in, among other things, enhanced edges.
+
+%%% YOUR CODE HERE %%%
+fprintf('Calculate xZCAWhite\n');
+xZCAWhite = U * xPCAWhite;
+