%SELO-sigma-KLDA.m
%SELO algorithm %Consider 10 families with 2 parents %and 2 kids in each
%family is %considered. %Training vectors subjected to KLDA. 
pk1=0.3;
pk2=0.7; 
X1=randn(3,4)*10+ones(3,4)*1;
X2=randn(3,4)*10+ones(3,4)*2;
X3=randn(3,4)*10+ones(3,4)*3; 
%Construction of Gram-matrix 
M=[X1 X2 X3]; 
%Initialize population (with 10 families, %2 parents and 2 kids in each family)
%Initially the ranges are identical for two parents %and two kids 
FL=[0.1 1.5 2.1 3.1 4.1 5.1 10 20 30 40];
FU=[1 2 3 4 5 10 20 30 40 50];
FL=[FL FL FL FL]; 
FU=[FU FU FU FU];
BESTVAL=[];
BESTSOL=[];
V=[]; 
F=[]; 
for i=1:1:40
    F(i)=0.1*rand*(FU(i)-FL(i))+FL(i);
    [PD,value]=seloalgo(M,F(i));
    V=[V value];
end
%Identifying the best member 
[P,Q]=min(V); 
[PD,value]=seloalgo(M,F(Q));
BESTVAL=value; 
BESTSOL=F(Q);
figure
plot(PD(1:1:4),zeros(1,4),'r*')
hold on 
plot(PD(5:1:8),zeros(1,4),'b*') 
plot(PD(9:1:12),zeros(1,4),'k*')

for iteration=1:1:20 
    V=[]; 
    F=[];
    for i=1:1:40
        F(i)=FL(i)+0.1*rand*(FU(i)-FL(i)); 
        [PD,value]=seloalgo(M,F(i));
        V=[V value];
    end
%Identifying the best member 
[P,Q]=min(V); 
if(V(Q)<BESTVAL(length(BESTVAL)))
    [PD,value]=seloalgo(M,F(Q));
    BESTVAL=[BESTVAL value]; 
    BESTSOL=[BESTSOL F(Q)]; 
end
%Adjusting the range for the parents 
for i=1:1:20 
    if(mod(Q,10)~=i) 
        temp=1./V(1:1:20); 
        R=cumsum(temp); 
        [temp]=cumsum(R)-rand; 
        [g,h]=find(temp>0); 
        FL(i)=F(i)-0.1*rand*norm([FL(h) FU(h)]);
        FU(i)=F(i)+0.1*rand*norm([FL(h) FU(h)]); 
    else if(Q<=20) FL(Q)=F(Q)-0.1*rand*norm([FL(Q) FU(Q)]); 
            FU(Q)=F(Q)+0.1*rand*norm([FL(Q) FU(Q)]); 
        end
    end
end
%%Update for kids 
for i=1:1:20 
    if(mod(Q,10)~=i) 
        if(rand>pk1)
temp=1./V(21:1:40); 
R=cumsum(temp);
[temp]=cumsum(R)-rand;
[g,h]=find(temp>0);
FL(20+i)=F(i+20)-0.1*rand*norm([FL(h+20) FU(h+20)]);
FU(20+i)=F(i+20)+0.1*rand*norm([FL(h+20) FU(h+20)]); 
        else
            temp=1./V(1:1:20);
R=cumsum(temp);
[temp]=cumsum(R)-rand;
[g,h]=find(temp>0);
FL(20+i)=F(i+20)-0.1*rand*norm([FL(h) FU(h)]);
FU(20+i)=F(i+20)+0.1*rand*norm([FL(h) FU(h)]);
        end
    else
        if(Q>20)
        if(rand>pk2) 
            if(rand>0.5) 
                FL(Q)=F(Q)-0.1*rand*norm([FL(i+20) FU(i+20)]);
                FU(Q)=F(Q)+0.1*rand*norm([FL(i+20) FU(i+20)]);
            else
                FL(Q)=F(Q)-0.1*rand*norm([FL(i+30) FU(i+30)]);
                FU(Q)=F(Q)+0.1*rand*norm([FL(i+30) FU(i+30)]);
            end
        else
            if(rand>0.5) 
                FL(Q)=F(Q)-0.1*rand*norm([FL(i) FU(i)]); 
                FU(Q)=F(Q)+0.1*rand*norm([FL(i) FU(i)]); 
            else
                FL(Q)=F(Q)-0.1*rand*norm([FL(i+10) FU(i+10)]); 
                FU(Q)=F(Q)+0.1*rand*norm([FL(i+10) FU(i+10)]); 
            end
        end
        end
    end
end
end
[PD,value]=seloalgo(M,BESTSOL(length(BESTSOL))); 
figure
plot(PD(1:1:4),zeros(1,4),'r*')
hold on, plot(PD(5:1:8),zeros(1,4),'b*') 
plot(PD(9:1:12),zeros(1,4),'k*')
figure
plot(BESTVAL) 
title('Convergence graph')
