Instantly share code, notes, and snippets.

anonymous /pca_animation.m
Created Jan 26, 2017

Embed
What would you like to do?
Matlab code to produce PCA animations
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab code to produce PCA animations shown here:
% http://stats.stackexchange.com/questions/2691
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Static image
clear all
rng(42)
X = randn(100,2);
X = X*chol([1 0.6; 0.6 0.6]);
X = bsxfun(@minus, X, mean(X));
[a,b] = eig(cov(X));
fig = figure('Position', [100 100 1000 400]);
set(gcf,'color','w');
axis([-3.5 3.5 -3.5 3.5])
axis square
hold on
scatter(X(:,1), X(:,2), 'b', 'filled')
%% Rotating animation
fig = figure('Position', [100 100 1000 400]);
set(gcf,'color','w');
axis([-3.5 3.5 -3.5 3.5])
axis square
hold on
for alpha = 1:1:179
w = [cosd(alpha) sind(alpha)];
z = X*w'*w;
cla
for i=1:100
plot([X(i,1) z(i,1)], [X(i,2) z(i,2)], 'r')
end
plot(w(1)*3.5*[-1 1], w(2)*3.5*[-1 1], 'k')
plot(-w(2)*2*[-1 1], w(1)*2*[-1 1], 'Color', [.6 .6 .6])
scatter(z(:,1), z(:,2), 'r', 'filled')
scatter(X(:,1), X(:,2), 'b', 'filled')
scatter(0,0,65,'k','filled', 'LineWidth', 2, 'MarkerFaceColor', [1 1 1], 'MarkerEdgeColor', [0 0 0])
a1 = 3.5;
a2 = 4.5;
plot(a(1,2)*[-a2 -a1], a(2,2)*[-a2 -a1], 'm', 'LineWidth', 2)
plot(a(1,2)*[ a1 a2], a(2,2)*[ a1 a2], 'm', 'LineWidth', 2)
drawnow
frame = getframe(fig);
if alpha == 1
[imind,map] = rgb2ind(frame.cdata,16,'nodither');
else
imind(:,:,1,alpha) = rgb2ind(frame.cdata,map,'nodither');
end
end
imwrite(imind,map, 'animation_pca.gif', 'DelayTime', 0, 'LoopCount', inf)
%% Pendulum animation
fig = figure('Position', [100 100 1000 400]);
set(gcf,'color','w');
axis([-3.5 3.5 -3.5 3.5])
axis square
hold on
alpha = -45;
omega = 0;
for t = 1:1:1000
w = [cosd(alpha) sind(alpha)];
z = X*w'*w;
M = sum(sum((z * [0 1; -1 0]) .* (X-z), 2));
omega = omega + M;
omega = omega * 0.93;
alpha = alpha + omega/40;
if(abs(omega)<1 && abs(alpha-abs(atand(a(2,2)/a(1,2)))) < 1)
break
end
cla
for i=1:100
plot([X(i,1) z(i,1)], [X(i,2) z(i,2)], 'r')
end
plot(w(1)*3.5*[-1 1], w(2)*3.5*[-1 1], 'k')
plot(-w(2)*2*[-1 1], w(1)*2*[-1 1], 'Color', [.6 .6 .6])
scatter(z(:,1), z(:,2), 'r', 'filled')
scatter(X(:,1), X(:,2), 'b', 'filled')
scatter(0,0,65,'k','filled', 'LineWidth', 2, 'MarkerFaceColor', [1 1 1], 'MarkerEdgeColor', [0 0 0])
a1 = 3.5;
a2 = 4.5;
plot(a(1,2)*[-a2 -a1], a(2,2)*[-a2 -a1], 'm', 'LineWidth', 2)
plot(a(1,2)*[ a1 a2], a(2,2)*[ a1 a2], 'm', 'LineWidth', 2)
pause(0.01)
drawnow
frame = getframe(fig);
if t == 1
[imind,map] = rgb2ind(frame.cdata,16,'nodither');
else
imind(:,:,1,t) = rgb2ind(frame.cdata,map,'nodither');
end
end
imwrite(imind,map, 'animation_pca_pendulum.gif', 'DelayTime', 0, 'LoopCount', inf)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment