1 function sheet05_02
   2 % Principal Component Analysis
   3 
   4 % simple data set
   5 X = randn(500, 2) * diag([1, 2]) * rotmat(0.3*pi);
   6 X = X + repmat([3, -4], 500, 1);
   7 
   8 plot(X(:,1), X(:,2), '.');
   9 axis equal
  10 grid
  11 
  12 % first PCA without centering
  13 figure(1)
  14 p = pca(X, 0)
  15 plot_pca(X, p);
  16 title('PCA without centering')
  17 
  18 % now, with proper centering
  19 figure(2)
  20 p = pca(X)
  21 plot_pca(X, p);
  22 title('PCA with centering')
  23 
  24 % and now, we add further, far away point
  25 figure(3)
  26 X = [X; -50, -50];
  27 p = pca(X);
  28 plot_pca(X, p);
  29 title('PCA is not very stable with respect to outliers')
  30 
  31 % how about other data?
  32 figure(4)
  33 phi = rand(500, 1) * 4 * pi;
  34 X = [phi, cos(phi) + (phi.^3)/1e2] + 0.5*randn(500, 2);
  35 
  36 p = pca(X);
  37 plot_pca(X, p);
  38 
  39 function R = rotmat(phi)
  40 R = [cos(phi), sin(phi); -sin(phi), cos(phi)];
  41 
  42 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  43 % Fill in your solutions below
  44 
  45 % 1. Compute the principal components. If the second argument is given and
  46 % non-zero, properly center the data. Otherwise, don't center.  Return a
  47 % struct with fields "U", "D", and "M" with
  48 %
  49 %    U is the matrix whose columns are the principal components
  50 %    D is a vector with the variances
  51 %    M is the mean vector
  52 function result = pca(X, center)
  53 if nargin == 1
  54   center = 1;
  55 end
  56 % ...
  57 result.U = U;
  58 result.D = D;
  59 result.M = M;
  60 
  61 % 2. Plot the principal components for the PCA. First, plot the data as
  62 % points. Then, draw an ellipsis whose radii are given as 2*sqrt(D) such
  63 % that the ellipsis contains the data more or less. Add a grid and make sure
  64 % axes are sized equally.
  65 function plot_pca(X, p)
  66 % ...