daniele@193: function SMALL_DL_test
daniele@193: clear, clc, close all
daniele@193: % Create a 2-dimensional dataset of points that are oriented in 3
daniele@194: % directions on a x/y plane
daniele@194: 
daniele@194: %
daniele@194: %   Centre for Digital Music, Queen Mary, University of London.
daniele@194: %   This file copyright 2012 Daniele Barchiesi.
daniele@194: %
daniele@194: %   This program is free software; you can redistribute it and/or
daniele@194: %   modify it under the terms of the GNU General Public License as
daniele@194: %   published by the Free Software Foundation; either version 2 of the
daniele@194: %   License, or (at your option) any later version.  See the file
daniele@194: %   COPYING included with this distribution for more information.
daniele@194: 
daniele@193: nData = 10000;						%number of data
daniele@193: theta = [pi/6 pi/3 4*pi/6];			%angles
daniele@194: nAngles = length(theta);			%number of angles
daniele@193: Q	  = [cos(theta); sin(theta)];	%rotation matrix
daniele@194: X	  = Q*randmog(nAngles,nData);	%training data
daniele@193: 
daniele@194: % find principal directions using PCA
daniele@194: XXt = X*X';							%cross correlation matrix
daniele@194: [U ~] = svd(XXt);					%svd of XXt
daniele@194: 
daniele@194: scale = 3;							%scale factor for plots
daniele@193: subplot(1,2,1), hold on
daniele@193: title('Principal Component Analysis')
daniele@194: scatter(X(1,:), X(2,:),'.');		%scatter training data
daniele@194: O = zeros(size(U));					%origin
daniele@194: quiver(O(1,1:2),O(2,1:2),scale*U(1,:),scale*U(2,:),...
daniele@194: 	'LineWidth',2,'Color','k')		%plot atoms
daniele@194: axis equal							%scale axis
daniele@193: 
daniele@193: subplot(1,2,2), hold on
daniele@193: title('K-SVD Dictionary')
daniele@193: scatter(X(1,:), X(2,:),'.');
daniele@193: axis equal
daniele@194: 
daniele@194: nAtoms  = 3;						%number of atoms in the dictionary
daniele@194: nIter   = 1;						%number of dictionary learning iterations
daniele@194: initDict = normc(randn(2,nAtoms));	%random initial dictionary
daniele@194: O = zeros(size(initDict));			%origin
daniele@194: 
daniele@193: % apply dictionary learning algorithm
daniele@194: ksvd_params = struct('data',X,...			%training data
daniele@194: 					 'Tdata',1,...			%sparsity level
daniele@193: 					 'dictsize',nAtoms,...	%number of atoms
daniele@194: 					 'initdict',initDict,...%initial dictionary
daniele@194: 					 'iternum',10);			%number of iterations
daniele@194: DL = SMALL_init_DL('ksvd','ksvd',ksvd_params); %dictionary learning structure
daniele@194: DL.D = initDict;							%copy initial dictionary in solution variable
daniele@194: problem = struct('b',X);					%copy training data in problem structure
daniele@194: 
daniele@193: xdata = DL.D(1,:);
daniele@193: ydata = DL.D(2,:);
daniele@193: qPlot = quiver(O(1,:),O(2,:),scale*initDict(1,:),scale*initDict(2,:),...
daniele@194: 	'LineWidth',2,'Color','k','UDataSource','xdata','VDataSource','ydata');
daniele@193: 
daniele@193: for iIter=1:nIter
daniele@193: 	DL.ksvd_params.initdict = DL.D;
daniele@194: 	pause
daniele@193: 	DL = SMALL_learn(problem,DL);			%learn dictionary
daniele@194: 	xdata = scale*DL.D(1,:);	
daniele@194: 	ydata = scale*DL.D(2,:);
daniele@193: 	refreshdata(gcf,'caller');
daniele@193: end
daniele@193: 
daniele@193: 
daniele@193: function X = randmog(m, n)
daniele@193: % RANDMOG - Generate mixture of Gaussians
daniele@193: s = [0.2 2];
daniele@193: % Choose which Gaussian
daniele@193: G1 = (rand(m, n) < 0.9);
daniele@193: % Make them
daniele@193: X = (G1.*s(1) + (1-G1).*s(2)) .* randn(m,n);