function [monproxpermut, vaf, diff] = proxmon(proxpermut, fitted)

%PROXMON produces a monotonically transformed proximity matrix (MONPROXPERMUT) 
%  from the order constraints obtained from each pair of entries in the input 
%  proximity matrix PROXPERMUT (symmetric with a zero main diagonal and a dissimilarity
%  interpretation).
%  MONPROXPERMUT is close to the $n \times n$ matrix FITTED in the least-squares sense;
%  The variance accounted for (VAF) is how much variance in MONPROXPERMUT can be accounted for by 
%  FITTED; DIFF is the value of the least-squares criterion.

n = size(proxpermut,1);
work = zeros(n*(n-1)*n*(n-1),1);
targ = proxpermut;
fit = fitted;
cr = 1.0;

while (cr >= 1.0e-006)
   
  
    
   cr = 0.0;
   indexll = 0;
   
   for jone = 1:(n-1)
      for jtwo = (jone+1):n
         for jthree = 1:(n-1)
            for jfour = (jthree+1):n
               
               if((jone ~= jthree) | (jtwo ~= jfour))
                  
                  p1 = fit(jone,jtwo);
                  p2 = fit(jthree,jfour);
                  fit(jone,jtwo) = fit(jone,jtwo) - work(indexll+1);
                  fit(jthree,jfour) = fit(jthree,jfour) - work(indexll+2);
                                    
                  
                  if((abs(targ(jone,jtwo) - targ(jthree,jfour)) ...
                        > 1.0e-006) & (targ(jone,jtwo) < ...
                        targ(jthree,jfour)))
                     
                     if(fit(jone,jtwo) <= fit(jthree,jfour))
                        
                        work(indexll+1) = 0;
                        work(indexll+2) = 0;
                        
                       
                        
                     elseif(fit(jone,jtwo) > fit(jthree,jfour))
                        
                        ave = (fit(jone,jtwo) + fit(jthree,jfour))/2.0;
                        work(indexll+1) = ave - fit(jone,jtwo);
                        work(indexll+2) = ave - fit(jthree,jfour);
                                                
                        fit(jone,jtwo) = ave;
                        fit(jthree,jfour) = ave;
                        
                    
                        
                     end
                     
                  elseif((abs(targ(jone,jtwo) - targ(jthree,jfour)) ...
                        > 1.0e-006) & (targ(jone,jtwo) > ...
                        targ(jthree,jfour)))
                     
                     if(fit(jone,jtwo) >= fit(jthree,jfour))
                        
                        work(indexll+1) = 0;
                        work(indexll+2) = 0;
                        
                       
                        
                     elseif(fit(jone,jtwo) < fit(jthree,jfour))
                        
                        ave = (fit(jone,jtwo) + fit(jthree,jfour))/2.0;
                        work(indexll+1) = ave - fit(jone,jtwo);
                        work(indexll+2) = ave - fit(jthree,jfour);
                        fit(jone,jtwo) = ave;
                        fit(jthree,jfour) = ave;
                        
                     
                        
                     end
                  end
                             
               cr = cr + abs(p1-fit(jone,jtwo)) + ...
                  abs(p2-fit(jthree,jfour));
            end
                      
               indexll = indexll + 2;
               
            end
         end
      end
   end
end

   
      
                  
                     
   for jone = 1:(n-1)
      for jtwo = (jone+1):n
                  
         fit(jtwo,jone) = fit(jone,jtwo);
         
      end
   end
   
  avefit = sum(sum(fit))/(n*(n-1));

for i = 1:n
   for j = 1:n
      if( i ~= j)
         proxave(i,j) = avefit;
      else
         proxave(i,j) = 0;
      end
   end
end

diff = sum(sum((fit - fitted).^2));

denom = sum(sum((fit - proxave).^2));

vaf = 1 - (diff/denom);

monproxpermut = fit;

diff = (.5)*diff;
               
      
               
   
   
   
   
