Daniel@0: function mlr_plot(X, Y, W, D)
Daniel@0: 
Daniel@0:     [d, n] = size(X);
Daniel@0:     %%%
Daniel@0:     % Color-coding
Daniel@0:     CODES = {   'b.', 'r.', 'g.', 'c.', 'k.', 'm.', 'y.', ...
Daniel@0:                 'b+', 'r+', 'g+', 'c+', 'k+', 'm+', 'y+', ...
Daniel@0:                 'b^', 'r^', 'g^', 'c^', 'k^', 'm^', 'y^', ...
Daniel@0:                 'bx', 'rx', 'gx', 'cx', 'kx', 'mx', 'yx', ...
Daniel@0:                 'bo', 'ro', 'go', 'co', 'ko', 'mo', 'yo'};
Daniel@0: 
Daniel@0:     %%%
Daniel@0:     % First, PCA-plot of X
Daniel@0: 
Daniel@0:     figure;
Daniel@0: 
Daniel@0:     z = sum(X,3);
Daniel@0:     subplot(3,2,[1 3]), pcaplot(z, eye(d), Y, CODES), title('Native');
Daniel@0:     subplot(3,2,[2 4]), pcaplot(X, W, Y, CODES), title('Learned');
Daniel@0: 
Daniel@0:     [vecs, vals] = eig(z * z');
Daniel@0:     subplot(3,2,5), bar(sort(real(diag(vals)), 'descend')), title('X*X'' spectrum'), axis tight;
Daniel@0: 
Daniel@0:     if size(X,3) == 1
Daniel@0:         [vecs, vals] = eig(W);
Daniel@0:         vals = real(diag(vals));
Daniel@0:     else
Daniel@0:         if size(W,3) == 1
Daniel@0:             vals = real(W(:));
Daniel@0:         else
Daniel@0:             vals = [];
Daniel@0:             for i = 1:size(W,3)
Daniel@0:                 [vecs, vals2] = eig(W(:,:,i));
Daniel@0:                 vals = [vals ; real(diag(vals2))];
Daniel@0:             end
Daniel@0:         end
Daniel@0:     end
Daniel@0:     subplot(3,2,6), bar(sort(vals, 'descend')), title('W spectrum'), axis tight;
Daniel@0: 
Daniel@0:     if nargin < 4
Daniel@0:         return;
Daniel@0:     end
Daniel@0: 
Daniel@0:     %%%
Daniel@0:     % Now show some diagnostics
Daniel@0: 
Daniel@0:     figure;
Daniel@0:     subplot(2,1,1), plot(D.f), title('Objective');
Daniel@0:     subplot(2,1,2), barh([D.time_SO, D.time_solver D.time_total]), ...
Daniel@0:                     title('Time% in SO/solver/total');
Daniel@0: 
Daniel@0: 
Daniel@0: function pcaplot(X, W, Y, CODES)
Daniel@0:     if size(X,3) == 1
Daniel@0:         A = kernelCenter(X' * W * X);
Daniel@0:     else
Daniel@0: %         A = zeros(size(X,1));
Daniel@0:         A = 0;
Daniel@0:         if size(W,3) == 1
Daniel@0:             for i = 1:size(X,3)
Daniel@0:                 A = A + X(:,:,i)' * bsxfun(@times, W(:,i), X(:,:,i));
Daniel@0:             end
Daniel@0:         else
Daniel@0:             if size(W,1) == size(W,2)
Daniel@0:                 for i = 1:size(W,1)
Daniel@0:                     for j = i:size(W,2)
Daniel@0:                         A = A + (X(:,:,i) + X(:,:,j)') * bsxfun(@times, squeeze(W(i,j,:)), X(:,:,i) + X(:,:,j));
Daniel@0:                     end
Daniel@0:                 end
Daniel@0: 
Daniel@0:             else
Daniel@0:                 for i = 1:size(X,3)
Daniel@0:                     A = A + X(:,:,i)' * W(:,:,i) * X(:,:,i);
Daniel@0:                 end
Daniel@0:             end
Daniel@0:         end
Daniel@0:     end
Daniel@0:     [v,d] = eigs(A, 3);
Daniel@0:     X2 = d.^0.5 * v';
Daniel@0:     hold on;
Daniel@0:     for y = 1:max(Y)
Daniel@0:         z = Y == y;
Daniel@0:         scatter3(X2(1, z), X2(2, z), X2(3,z), CODES{y});
Daniel@0:     end
Daniel@0: 
Daniel@0:     axis equal;