clear all
close all
clc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% genera dati
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all
close all
clc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MU1=[0,0];
SIGMA1=[0.3 0;0 2];
Npoints=200;
X1 = mvnrnd(MU1,SIGMA1,Npoints);
X1=X1';
%plot(X1(1,:),X1(2,:),'co')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%MU2=[7,7];
MU2=[5,5];
SIGMA2=[1.5 0.7;0.7 1.5];
X2 = mvnrnd(MU2,SIGMA2,Npoints);
X2=X2';
%hold on
%plot(X2(1,:),X2(2,:),'bo')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MU3=[-5,-10];
SIGMA3=[1.5 0.7;0.7 1.5];
X3 = mvnrnd(MU3,SIGMA3,Npoints);
X3=X3';
%hold on
%plot(X3(1,:),X3(2,:),'ro')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MU4=[10,-10];
SIGMA4=[1.5 0;0 1.5];
X4 = mvnrnd(MU4,SIGMA4,Npoints);
X4=X4';
%hold on
%plot(X4(1,:),X4(2,:),'mo')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%MU5=[0,-10];
%MU5=[-5,10];
MU5=[-5,5];
SIGMA5=[1 -0.8;-0.8 1];
X5 = mvnrnd(MU5,SIGMA5,Npoints);
X5=X5';
%hold on
%plot(X5(1,:),X5(2,:),'ko')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
X=[X1 X2 X3 X4 X5];
%X=[X1 X2  X4 X3];
%X=[X1 X2];
[DIM,Ntot]=size(X);
h1=figure
plot(X(1,:),X(2,:),'co')
hold on
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
k=5; %%% numero clusters

%for k=2:10

C=-5+10*rand(2,k);
iter=10;
%%%%%%%%%%%%%%%%%%%
Ctot1=C(1,:);
Ctot2=C(2,:);
%%%%%%%%%%%%%%%%%
%%%% posiciones iniciales de los centroides
for j=1:k
  plot(C(1,j),C(2,j),'sk','MarkerEdgeColor','r','MarkerFaceColor','b','MarkerSize',10)
end

for r=1:iter
r

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
dist=@(x,z) sum((x-z).^2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:Ntot
    for j=1:k
       D(i,j)=dist(X(:,i),C(:,j));
    end
     clust_index(i)=find(D(i,:)==min(D(i,:)));     
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%  clust_index = say the index of the cluster that the sample bolongs
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for j=1:k
  pos_j=find(clust_index==j);
    if isempty(pos_j)==0
      C(:,j)=mean(X(:,pos_j),2);
        E(r,j)=sum(var(X(:,pos_j)'));
    end
      
end
%%%%%%%%%%
Ctot1(r+1,:)=C(1,:);
Ctot2(r+1,:)=C(2,:);
end %%%%% end iter of kmeans
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% plot final locations of centroid %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for j=1:k
  plot(C(1,j),C(2,j),'dk','MarkerEdgeColor','r','MarkerFaceColor','g','MarkerSize',10)
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%% Silhoutte %%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
vec_aux=1:k;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:length(X)
     %%%%%%%%%%%%%%%%%%%%%%%%%%%
     clusnow=clust_index(i);
     vec_aux_now=vec_aux;
     vec_aux_now(clusnow)=[];
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     pos_data_in_clus=find(clust_index==clust_index(i)); %%%% indices of all samples in a cluster
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     count_in_clus=length(pos_data_in_clus); %%%%% number of samples in a cluster 
     for j=1:count_in_clus
         %%%% distancia entre la i-esima muestra y las otras muestras en
         %%%% el MISMO cluster
         distaux(j)=dist(X(:,i),X(:,pos_data_in_clus(j))); 
     end
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     %%%%% media distancias entre la i-esima muestras y las otras muestras
     %%%%% el mismo cluster 
     a(i)=1/(count_in_clus-1)*sum(distaux);
      %%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
     for j=1:k-1
            pos_now=find(clust_index==vec_aux_now(j)); %%%% otro cluster
            count_now=length(pos_now); %%%% numero de muestras del otro cluster
         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
          for cc=1:count_now
              %%%% distancia entre la i-esima muestra y las otras muestras en
              %%%% el otro cluster
            psi_quasi(i,cc)= dist(X(:,i),X(:,pos_now(cc)));
          end
           psi(i,j)=mean(psi_quasi(i,:));
     end
    %%%%%%%%%%%%%%%%%%%%%%%  %%%%%%%%%%%%%%%%%%%%%%%  %%%%%%%%%%%%%%%%%%%%%%%   
     b(i)=min(psi(i,:)); %%%% minima distancia con otro cluster la muestra i-esima
     %%%% cluster alternativo la muestra i-esima
     aux_pos_now=find(psi(i,:)==b(i));
     cluster_index_alter(i)=vec_aux_now(aux_pos_now);
         %%%%%%%%%%%%%%%%%%%%%%%
    s(i)=(b(i)-a(i))/max(b(i),a(i));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure
plot(s)
disp('---------------------------------------------------------------')
disp('---------------------------------------------------------------')
disp('--- Respecto a la muestra peor agrupada ---')
posmin_s=find(s==min(s));
disp('que es:')
X(1,posmin_s)
X(2,posmin_s)
hold on
plot(posmin_s,s(posmin_s),'ro')
figure(h1)
plot(X(1,posmin_s),X(2,posmin_s),'rs','LineWidth',2,'MarkerEdgeColor','k','MarkerFaceColor','r','MarkerSize',10)
disp(['Centroide del cluster actual de la muestra peor agrupada num. ',num2str(clust_index(posmin_s)),' :'])
C(1,clust_index(posmin_s))
C(2,clust_index(posmin_s))
disp(['Centroide del posible cluster alternativo de la muestra peor agrupada num. ',num2str(cluster_index_alter(posmin_s)),' :'])
C(1,cluster_index_alter(posmin_s))
C(2,cluster_index_alter(posmin_s))
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%indice_k(k)=mean(s)
disp('---------------------------------------------------------------')
disp('---------------------------------------------------------------')
disp('Indice total de la silhouette de la conf. de centroides:')
indice_k=mean(s)
%%%%%%%%%
%end

%figure
%plot(indice_k)