function [membership,objectives] = partitionfnd_averages(prox)

% PARTITIONFND_AVERAGES uses dynamic programming to
% construct a linearly constrained cluster analysis that
% consists of a collection of partitions with from 1 to
% n ordered classes.
%
% syntax: [membership,objectives] = partitionfnd_averages(prox)
%
% PROX is the input proximity matrix (with a zero main diagonal
% and a dissimilarity interpretation);
% MEMBERSHIP is the n x n matrix indicating cluster membership,
% where rows correspond to the number of ordered clusters, 
% and the columns are in the identity permutation input order 
% used for PROX. 
% OBJECTIVES is the vector of merit values minimized in the
% construction of the ordered partitions, each defined by the
% maximum over clusters of the average proximities within subsets.

%initializing arrays

n = size(prox,1);
xi = (1.0e+010)*ones(n,n);
nxi = zeros(n,n);
clus = zeros(n,n);

%initial single subset processing and
%homogeneity calculations

for nlb = 1:n
    for nub = nlb:n
        
                 
                subsetsize = (nub - nlb) + 1;              
               
          if(subsetsize == 1) 
              
              index = 0;
              
          end
          
          if(subsetsize >= 2)
              
              index = (sum(sum(prox(nlb:nub,nlb:nub))))/(subsetsize*(subsetsize - 1));
              
          end
          
          
                                 
     
        
        if (nlb == 1)
            xi(1,nub) = index;
            nxi(1,nub) = 1;
        end
        
        clus(nlb,nub) = index;
       
    end
end

%processing the recursion

for k = 2:n
    
    for i = (k-1):(n-1)
        for jone = (i+1):n
            
            temp = max(xi(k-1,i), clus(i+1,jone));
            
            if(temp < xi(k,jone))
                
                xi(k,jone) = temp;
                nxi(k,jone) = i + 1;
                
            end
        end
    end
end



%retriving membership

membership=[];
objectives=[];
for j = 1:n
for k = j:j
    membershipvec=zeros(1,n);
    cluster = 0;

    
    nend = n;
    nk = k;
    objective = xi(nk,nend);
    
          
        
        add=nxi(nk,nend):nend;
        cluster = cluster+1;
        a = length(add);
        for i=1:a
            membershipvec(add(i)) = cluster;
        end
      
        
      
     while(nxi(nk,nend) ~= 1)
        
        nend = nxi(nk,nend) - 1;
        nk = nk - 1;
                        
        
        add=nxi(nk,nend):nend;
        cluster = cluster+1;
        a = length(add);
        for i=1:a
            membershipvec(add(i)) = cluster;
        end

        
    end
    cluster = cluster+1;
    
end
membership = [membership; membershipvec];

objectives = [objectives;objective];
end


