function [supAlph,supVector,bias,className,ite,maxE,gap,E]=mySMO(X,Y,className,C,tol,kerF,maxIte)
% Syntax: [supAlph,supVector,bias,className,ite,maxE,gap,E]=mySMO5(X,Y,className,C,tol,kerF,maxIte).
rng(1); % For reproducibility
fc=find(Y(:,1)==className(1,1));
sc=find(Y(:,1)==className(1,2));
Xlength=size(fc,1)+size(sc,1);
supAlph=[];
supVector=[];
maxE=[];
gap=[];
IsParaTrue=1;
if(Xlength~=size(X,1) || className(1,1)==className(1,2) )
IsParaTrue=0;
end
Alph=zeros(Xlength,1);
Yclass=ones(Xlength,1);
Yclass(sc)=-1;
bias=0;
E=ones(Xlength,1);
E(fc)=-1;
Con=zeros(Xlength,1);
if size(C,1)==1 || size(C,1)==size(X,1)
Con=C;
if size(C,1)==1
Con=C*ones(Xlength,1);
end
else
IsParaTrue=0;
end
isAll=1;
ite=[0;0];
incre=0;
% i1=0;
% i2=0;
% vbias=0;
% AlphOld1=0;
% AlphOld2=0;
% delta1=0;
% delta2=0;
while (incre>0 || isAll==1 ) && ite(1,1)+ite(2,1)<maxIte && IsParaTrue==1
incre=0;
i1=0;
i2=0;
vbias=0;
AlphOld1=0;
AlphOld2=0;
delta1=0;
delta2=0;
if (isAll~=1)
index=find(Alph>0);
sizeIndex=size(index,1);
randIndex=randperm(sizeIndex);
for indexi1=1:sizeIndex
if(incre>0)
break;
end
i1=index(randIndex(1,indexi1));
delta1=Yclass(i1,1)*E(i1,1);
vbias=0;
i2=0;
AlphOld1=Alph(i1,1);
AlphOld2=0;
if (delta1>tol && AlphOld1~=0) || (AlphOld1~=Con(i1,1)&& delta1<-tol)
randi2=randperm(Xlength);
for i3=1:Xlength
i2=randi2(1,i3);
if i1~=i2
delta2=Yclass(i2,1)*E(i2,1);
if (delta2>tol && Alph(i2,1)~=0) || (Alph(i2,1)~=Con(i2,1)&& delta2<-tol)
AlphOld2=Alph(i2,1);
kc=feval(kerF,X(i1,:),X(i1,:))+feval(kerF,X(i2,:),X(i2,:))-2*feval(kerF,X(i1,:),X(i2,:));
[incre,Alph(i1,1),Alph(i2,1)]=takestep(Yclass(i1,1),Yclass(i2,1),Alph(i1,1),Alph(i2,1),Con(i1,1),Con(i2,1),E(i1,1),E(i2,1),kc);
if incre>0
break;
end
end
end
end
end
end
end
if (isAll==1)
randi=randperm(Xlength);
for i=1:Xlength
if(incre>0)
break;
end
i1=randi(1,i);
delta1=Yclass(i1,1)*E(i1,1);
vbias=0;
i2=0;
AlphOld1=Alph(i1,1);
AlphOld2=0;
if (delta1>tol && AlphOld1~=0) || (AlphOld1~=Con(i1,1)&& delta1<-tol)
for j=i+1:Xlength
i2=randi(1,j);
delta2=Yclass(i2,1)*E(i2,1);
if (delta2>tol && Alph(i2,1)~=0) || (Alph(i2,1)~=Con(i2,1)&& delta2<-tol)
AlphOld2=Alph(i2,1);
kc=feval(kerF,X(i1,:),X(i1,:))+feval(kerF,X(i2,:),X(i2,:))-2*feval(kerF,X(i1,:),X(i2,:));
[incre,Alph(i1,1),Alph(i2,1)]=takestep(Yclass(i1,1),Yclass(i2,1),Alph(i1,1),Alph(i2,1),Con(i1,1),Con(i2,1),E(i1,1),E(i2,1),kc);
end
if(incre>0)
break;
end
end
end
end
if (incre<1)
randi=randperm(Xlength);
for i=1:Xlength
if(incre>0)
break;
end
i1=randi(1,i);
delta1=Yclass(i1,1)*E(i1,1);
vbias=0;
i2=0;
AlphOld1=Alph(i1,1);
AlphOld2=0;
if (delta1>tol && AlphOld1~=0) || (AlphOld1~=Con(i1,1)&& delta1<-tol)
for j=1:Xlength
if j~=i
i2=randi(1,j);
AlphOld2=Alph(i2,1);
kc=feval(kerF,X(i1,:),X(i1,:))+feval(kerF,X(i2,:),X(i2,:))-2*feval(kerF,X(i1,:),X(i2,:));
[incre,Alph(i1,1),Alph(i2,1)]=takestep(Yclass(i1,1),Yclass(i2,1),Alph(i1,1),Alph(i2,1),Con(i1,1),Con(i2,1),E(i1,1),E(i2,1),kc);
end
if(incre>0)
break;
end
end
end
end
end
end
if incre>0
% [i1,i2,Yclass(i1,1),Yclass(i2,1),AlphOld1,AlphOld2,Alph(i1,1), Alph(i2,1),delta1,delta2,isAll]
if Alph(i1,1)~=Con(i1,1)&& Alph(i1,1)~=0
vbias=(AlphOld1-Alph(i1,1))*Yclass(i1,1)*feval(kerF,X(i1,:),X(i1,:))+...
(AlphOld2-Alph(i2,1))*Yclass(i2,1)*feval(kerF,X(i2,:),X(i1,:))-E(i1,1);
else if Alph(i2,1)~=Con(i2,1)&& Alph(i2,1)~=0
vbias=(AlphOld1-Alph(i1,1))*Yclass(i1,1)*feval(kerF,X(i1,:),X(i2,:))+...
(AlphOld2-Alph(i2,1))*Yclass(i2,1)*feval(kerF,X(i2,:),X(i2,:))-E(i2,1);
end
end
bias=bias+vbias;
index=find(Alph>0);
supAlph=Alph(index).*Yclass(index);
supVector=X(index,:);
E=zeros(Xlength,1);
for i=1:Xlength
for j=1:size(supAlph,1)
E(i,1)=E(i,1)+supAlph(j,1)*feval(kerF,supVector(j,:),X(i,:));
end
E(i,1)= E(i,1)+bias-Yclass(i,1);
end
% scores=myPredict(supAlph,supVector,bias,kerF,X);
% E=scores-Yclass;
isAll=0;
ite(1,1)=ite(1,1)+1;
else if isAll==0
isAll=1;
else
isAll=0;
end
ite(2,1)=ite(2,1)+1;
end
if mod(ite(1,1)+ite(2,1),500)==0
[ite(1,1),ite(2,1),ite(1,1)+ite(2,1)]
end
end
if IsParaTrue==1
supAlph=[];
supVector=[];
for i=1:size(E,1)
if Alph(i,1)>0
supAlph=[supAlph;Alph(i,1)*Yclass(i,1),i];
supVector=[supVector;X(i,:)];
end
end
scores=myPredict(supAlph,supVector,bias,kerF,X);
eror=(scores-Yclass).*Yclass;
gE=-(scores-Yclass).*Yclass;
sE=-(scores-Yclass).*Yclass;
for i=1:size(gE,1)
if sE(i,1)<0
sE(i,1)=0;
end
if (Alph(i,1)==0 && gE(i,1)<0)
gE(i,1)= 0;
end
if Alph(i,1)==Con(i,1) && gE(i,1)>0
gE(i,1)=0;
end
end
%%% sE��vector of slack variables;
gap=abs(Con'*sE+Alph'*eror);
maxE=max(abs(gE));
E=[Alph(:),Yclass(:),scores,E(:),gE(:),sE(:)];
end
end
function [Inc, AlphNew1, AlphNew2]= takestep(Yclass1,Yclass2,Alph1,Alph2,Con1,Con2,E1,E2,kc)
Inc=0;
AlphNew1=Alph1;
AlphNew2=Alph2;
if Yclass1==Yclass2
H=max(0,Alph1+Alph2-Con2);
L=min(Con1,Alph1+Alph2);
else
H=max(0,Alph1-Alph2);
L=min(Con1,Alph1-Alph2+Con2);
end
if H<L
if kc~=0
AlphNew1=Alph1+Yclass1*(E2-E1)/kc;
if AlphNew1>L
