% Spearman Rank Order Correlation
% Written by Wonho Song, October 2002
% E-mail: whsong@cau.ac.kr, whsong73@hotmail.com

function [Rs,rr]=spearman(a)

% Rs is the Spearman rank correlation
% rr is the ranks matrix

[n,p]=size(a);


rr=[];

for i=1:p
rr=[rr rankindx(a(:,i))];
end;

  Rs=zeros(p,p);

for i=1:p
  for j=1:p
    r1=rr(:,i); r2=rr(:,j);
    d=r1-r2;
    d2=d'*d;
    Rs(i,j)=1-6*d2/(n^3-n);
  end;
end;



function [indexs]=rankindx(a)

[n,p]=size(a);

[x,I]=sort(a);

index=zeros(n,1);

index=1:n;
index=index';

i=1;

while i<n   % terminate at n-1

   if x(i)==x(i+1); 
      j=i;   % start

      if i<n-1; 
           while x(i+1)==x(i+2); i=i+1;
                       if i==n-1; break; end;
           end;  
      end;     
      index(j:i+1)=(j+i+1)/2;


   end;  % end of if
 i=i+1;
end;  % end of while

indexs=zeros(n,1);

for i=1:n
indexs(I(i))=index(i);
end;