function lyngby_ikm_plot_w(Within, Clusters)

% lyngby_ikm_plot_w    - Plot Within variance for iterative K-means
%
%       function lyngby_ikm_plot_w(Within, Cluster)
%
%       Make a plot of the intracluster variance as a function of the
%       number of clusters together with the curvature of the
%       function. 
%
%       Similar to a plot in Goutte et al: "On Clustering fMRI Time
%       series", NeuroImage 9(3):298-310, 1999 March 
%
%       See also LYNGBY, LYNGBY_IKM_MAIN, LYNGBY_KM_MAIN,
%                LYNGBY_KM_WITHIN. 
%
% $Id: lyngby_ikm_plot_w.m,v 1.2 2003/02/03 17:47:49 fnielsen Exp $ 



    minWithin = min(Within, [], 2);

    for n = 2:length(Clusters)-1
      y = [ minWithin(n-1) minWithin(n) minWithin(n+1) ]';
      X = [ Clusters(n-1)^2 Clusters(n-1) 1 ; ...
	Clusters(n)^2 Clusters(n) 1 ; ...
	Clusters(n+1)^2 Clusters(n+1) 1 ];
      theta = inv(X) * y;
      curvature(n) = 2*theta(1);
    end
    
    scale = floor(max(Within(:,1))/max(curvature));
    
    h1 = bar(Clusters, minWithin, 'y');
    hold on
    h2 = plot(Clusters, Within, 'rx');
    h3 = plot(Clusters(2:end-1), scale*curvature(2:end));
    set(h3, 'LineWidth', 3, 'Color', [0 0 0]);
    hold off
    
    % axis([0 max(Clusters)+1 -1 max(Within(:,1))*1.1])
    
    xlabel('Number of clusters');
    ylabel('Within variance');
    legend([ h1 h2(1) h3], 'Minimum Within variance', 'Within variance', ...
	sprintf('Curvature (x%d)', scale)) 
    title('Within variance');


