function[c,sigma,W_output]=SOFNN(X,d,Kd)%SOFNNSelf-OrganizingFuzzyNeuralNetworks%InputParameters%X(r,n)-rthtraningdatafromnthobservation%d(n)-thedesiredoutputofthenetwork(mustbearowvector)%Kd(r)-predefineddistancethresholdfortherthinput%OutputParameters%c(IndexInputVariable,IndexNeuron)%sigma(IndexInputVariable,IndexNeuron)%W_outputisavector%SettingupParametersforSOFNNSigmaZero=4;delta=0.12;threshold=0.1354;k_sigma=1.12;%Formoreaccurateresultsuncommentthefollowing%formatlong;%ImplementationofaSOFNNmodel[size_R,size_N]=size(X);%size_R-thenumberofinputvariablesc=[];sigma=[];W_output=[];u=0;%thenumberofneuronsinthestructureQ=[];O=[];Psi=[];forn=1:size_Nx=X(:,n);ifu==0%Noneuroninthestructure?c=x;sigma=SigmaZero*ones(size_R,1);u=1;Psi=GetMePsi(X,c,sigma);[Q,O]=UpdateStructure(X,Psi,d);pT_n=GetMeGreatPsi(x,Psi(n,:))';else[Q,O,pT_n]=UpdateStructureRecursively(X,Psi,Q,O,d,n);end;KeepSpinning=true;whileKeepSpinning%Calculatetheerrorandif-partcriteriaae=abs(d(n)-pT_n*O);%approximationerror[phi,~]=GetMePhi(x,c,sigma);[maxphi,maxindex]=max(phi);%maxindexreferstotheneuron'sindexifaedeltaifmaxphithreshold%enlargewidth[minsigma,minindex]=min(sigma(:,maxindex));sigma(minindex,maxindex)=k_sigma*minsigma;Psi=GetMePsi(X,c,sigma);[Q,O]=UpdateStructure(X,Psi,d);pT_n=GetMeGreatPsi(x,Psi(n,:))';else%Addanewneuronandupdatestructurectemp=[];sigmatemp=[];dist=0;forr=1:size_Rdist=abs(x(r)-c(r,1));distIndex=1;forj=2:uifabs(x(r)-c(r,j))distdistIndex=j;dist=abs(x(r)-c(r,j));end;end;ifdist=Kd(r)ctemp=[ctemp;c(r,distIndex)];sigmatemp=[sigmatemp;sigma(r,distIndex)];elsectemp=[ctemp;x(r)];sigmatemp=[sigmatemp;dist];end;end;c=[cctemp];sigma=[sigmasigmatemp];Psi=GetMePsi(X,c,sigma);[Q,O]=UpdateStructure(X,Psi,d);KeepSpinning=false;u=u+1;end;elseifmaxphithreshold%enlargewidth[minsigma,minindex]=min(sigma(:,maxindex));sigma(minindex,maxindex)=k_sigma*minsigma;Psi=GetMePsi(X,c,sigma);[Q,O]=UpdateStructure(X,Psi,d);pT_n=GetMeGreatPsi(x,Psi(n,:))';else%DonothingandexitthewhileKeepSpinning=false;end;end;end;end;W_output=O;endfunction[Q_next,O_next,pT_n]=UpdateStructureRecursively(X,Psi,Q,O,d,n)%O=O(t-1)O_next=O(t)p_n=GetMeGreatPsi(X(:,n),Psi(n,:));pT_n=p_n';ee=abs(d(n)-pT_n*O);%|e(t)|temp=1+pT_n*Q*p_n;ae=abs(ee/temp);ifee=aeL=Q*p_n*(temp)^(-1);Q_next=(eye(length(Q))-L*pT_n)*Q;O_next=O+L*ee;elseQ_next=eye(length(Q))*Q;O_next=O;end;endfunction[Q,O]=UpdateStructure(X,Psi,d)GreatPsiBig=GetMeGreatPsi(X,Psi);%M=u*(r+1)%n-thenumberofobservations[M,~]=size(GreatPsiBig);%OthersWaysofgettingQ=[P^T(t)*P(t)]^-1%**************************************************************************%opts.SYM=true;%Q=linsolve(GreatPsiBig*GreatPsiBig',eye(M),opts);%%Q=inv(GreatPsiBig*GreatPsiBig');%Q=pinv(GreatPsiBig*GreatPsiBig');%**************************************************************************Y=GreatPsiBig\eye(M);Q=GreatPsiBig'\Y;O=Q*GreatPsiBig*d';end%Thisfunctionworkstoowithx%(X=XandPsiisaMatrix)-GetsyouthewholeGreatPsi%(X=xandPsiistherowrelatedtox)-Getsyoujustthecolumnrelatedwiththeobservationfunction[GreatPsi]=GetMeGreatPsi(X,Psi)%Psi-Inarowyougothroughtheneuronsandinacolumnyougothroughnumberof%observations****Psi(#obs,IndexNeuron)****GreatPsi=[];[N,U]=size(Psi);forn=1:Nx=X(:,n);GreatPsiCol=[];foru=1:UGreatPsiCol=[GreatPsiCol;Psi(n,u)*[1;x]];end;GreatPsi=[GreatPsiGreatPsiCol];end;endfunction[phi,SumPhi]=GetMePhi(x,c,sigma)[r,u]=size(c);%u-thenumberofneuronsinthestructure%r-thenumberofinputvariablesphi=[];SumPhi=0;forj=1:u%movingthroughtheneuronsS=0;fori=1:r%movingthroughtheinputvariablesS=S+((x(i)-c(i,j))^2)/(2*sigma(i,j)^2);end;phi=[phiexp(-S)];SumPhi=SumPhi+phi(j);%phi(u)=exp(-S)end;end%Thisfunctionworkstoowithx,itwillgiveyoutherowrelatedtoxfunction[Psi]=GetMePsi(X,c,sigma)[~,u]=size(c);[~,size_N]=size(X);%u-thenumberofneuronsinthestructure%size_N-thenumberofobservationsPsi=[];forn=1:size_N[phi,SumPhi]=GetMePhi(X(:,n),c,sigma);PsiTemp=[];forj=1:u%PsiTempisarowvectorex:[123]PsiTemp(j)=phi(j)/SumPhi;end;Psi=[Psi;PsiTemp];%Psi-Inarowyougothroughtheneuronsandinacolumnyougothroughnumberof%observations****Psi(#obs,IndexNeuron)****end;end