We don't know why there is a error. How do you expect us to know? Either you should provide the input file needed to run the given code or you should mention what error you get.
Why there is an error in line 132? ( [Fit(1,j), tetaDE]=selectmodel2(Pop(:,j),NP,D,mphi,start,nodata,y,u);)
1 vue (au cours des 30 derniers jours)
Afficher commentaires plus anciens
% %this program for comparison between MGA & DE
clear;
%load input-output data
tt=cputime;
load ui.dat; %load input
u=ui(:,1);
mx=2;
y(1:mx+1)=0;
for t=mx+1:500
% for i=1:500
% e(i)=(rand-0.5)*0.002;
% end
% y(t)=0.797*y(t-1)-0.237*y(t-3)+0.441*u(t-1)+0.105*y(t-4)*u(t-4)...
% +0.333*u(t-3)*u(t-5);%+e(t);
y(t)=0.3*u(t-1)*u(t-2)-0.8*y(t-2)*u(t-1)*u(t-1)-0.1*u(t-1)*u(t-1)*u(t-2);%+e(t);
end
y=y';%load output
nodata=500;
noy=2; nou=2;
nl=3; start=max(nou,noy)+1;
m=nou+noy;
mphi=[]; %form regressor for ARMAX model
for i=start:nodata
phi=[];
if nl==1
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
end
if nl==2
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
for j=1:noy
for l=j:noy
phi=[phi y(i-j)*y(i-l)];
end
for k=1:nou
phi=[phi y(i-j)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi u(i-k)*u(i-n)];
end
end
end
if nl==3
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
for j=1:noy
for l=j:noy
phi=[phi y(i-j)*y(i-l)];
end
for k=1:nou
phi=[phi y(i-j)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi u(i-k)*u(i-n)];
end
end
for j=1:noy
for l=j:noy
for p=l:noy
phi=[phi y(i-j)*y(i-l)*y(i-p)];
end
for k=1:nou
phi=[phi y(i-j)*y(i-l)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi y(i-j)*u(i-k)*u(i-n)];
end
end
end
for k=1:nou
for n=k:nou
for q=n:nou
phi=[phi u(i-k)*u(i-n)*u(i-q)];
end
end
end
end
mphi=[mphi;phi]; %end regressor
end
[temp m]=size(mphi);
%DIFFERENTIAL EVOLUTION
%CONTROL PARAMETERS
D=m; %dimension of problem
NP=30; %Population size
F=0.9; %differentiation constant
CR=0.5; %crossover rate
GEN=100; %generations
L=0.0; %low boundary constraint 1
H=1.0; %high boundary constraint 1
% L2=4.1; %low boundary constraint 2
% H2=5.8; %high boundary constraint 2
%**********************************%
% Algorithm Variables
%**********************************%
X=zeros(D,1); %trial vctor X1
% X2=zeros(D,1); %trial vctor X2
Pop=zeros(D,NP); %population X1
% Pop2=zeros(D,NP); %population X1
Fit=zeros(1,NP); %fitness population
iBest=1; %index of best solution
r=zeros(3,1); %randomly selected indices
%**********************************%
% Create Population
%**********************************%
%Initialize random number of generator
%rand('state',sum(100*clock));
for j=1:NP %Initialize each individual
Pop(:,j)=L+(H-L)*rand(D,1); %within b.constraint 1
% Pop2(:,j)=L2+(H2-L2)*rand(D,1); %within b.constraint 2
[Fit(1,j), tetaDE]=selectmodel2(Pop(:,j),NP,D,mphi,start,nodata,y,u); %evaluate fitness 1
end
%**********************************%
% Optimization
%**********************************%
fbestDE=[];
for g=1:GEN %for each generation
for j=1:NP %for each individual
%choose 3 random individuals from population
%mutually different and different from j
r(1)=floor(rand()*NP)+1;
while r(1)==j
r(1)=floor(rand()*NP)+1;
end
r(2)=floor(rand()*NP)+1;
while (r(2)==r(1))||(r(2)==j)
r(2)=floor(rand()*NP)+1;
end
r(3)=floor(rand()*NP)+1;
while (r(3)==r(2))||(r(3)==r(1))||(r(3)==j)
r(3)=floor(rand()*NP)+1;
end
%create trial individual
%in which at least one parameter is changed
Rnd=floor(rand()*D)+1;
for i=1:D
if (rand()<CR)||(Rnd==i)
X(i)=Pop(i,r(3))+F*(Pop(i,r(1))-Pop(i,r(2)));
% X2(i)=Pop2(i,r(3))+F*(Pop2(i,r(1))-Pop2(i,r(2)));
else
X(i)=Pop(i,j);
% X2(i)=Pop2(i,j);
end
end
%verify b.constraint
for i=1:D
if (X(i)<L)||(X(i)>H)
X(i)=L+(H-L)*rand();
end
% if (X2(i)<L2)||(X2(i)>H2)
% X2(i)=L2+(H2-L2)*rand();
% end
end
%select the best individual
%between trial and current ones
%calculate fitness of trial vector
[f, tetaDE]=selectmodel2(X,NP,D,mphi,start,nodata,y,u);
%if trial is better or equal than current
% if f<=Fit(j)
% Pop1(:,j)=X1; %replace current by trial
% Pop2(:,j)=X2;
% Fit(j)=f;
% %if trial is better than the best
% if f<=Fit(iBest)
% iBest=j; %update the best's index
% end
% end
if f<=Fit(j)
Pop(:,j)=X; %replace current by trial
% Pop2(:,j)=X2;
Fit(j)=f;
%if trial is better than the best
if f<=Fit(iBest)
iBest=j; %update the best's index
end
end
end
fbestDE=[fbestDE Fit(iBest)];
numgenDE(g)=g;
end
% %MGA
%CONTROL PARAMETERS
evol=100;%input('# of evolution >:');
popsize=30;%input('population size >:');
pcr=0.8;%input('crossover probability=');
pmut=0.01;%input('mutation probability=');
chr=initbp(popsize,m); %initialize population
chrom=chr;
mseGA=[];
for t=1:evol %start evolution
gen1=[];
for i=1:popsize %model structure chrom's=1
mmphi=[];
krom=chrom(i,:);
for b=1:m
index=krom(:,b);
if index==1
mmphi=[mmphi mphi(:,b)];
end
end
[N, g]=size(mmphi);
yt=y(start:nodata);ut=u(start:nodata);
y2=sum(yt.*yt);
teta=pinv(mmphi'*mmphi)*mmphi'*yt; %find teta by LSE
num=0;
for c=1:g
tetahat=teta(c,:);
if tetahat<0.001
num=num+0.1; %penalty to added noise
end
end
yhat=mmphi*teta;
error=sum((yt-yhat).*(yt-yhat))+num;
gen=[krom error];
gen1=[gen1;gen];
pm=max(gen1(:,m+1));
fitness=pm-gen1(:,m+1);
% bestftn=pm-fitness-num;
gene1=[gen1 fitness];
end
generation=sortrows(gene1,[m+1]);
% mseGA=[mseGA gene1(:,m+1)];
sse=generation(:,m+1);
bestchr=0.0; secondchr=0.0;
for i=1:popsize
sume2=sse(i,:);
if sume2<1.0
bestchr=bestchr+1; %good chromosomes
else
secondchr=secondchr+1; %2nd class chromosomes
end
end
sumerror2(t)=sum(gene1(:,m+1));
fitnv=gene1(:,m+2); %assign fitness
sumfit(t)=sum(fitnv);
%newchr=selsus(fitnv,popsize);
newchr=seltour(fitnv,popsize,2);
selchrom=chrom(newchr,:);
%selchrom=generation(:,1:m);
mutchrom=mutbint(selchrom(1:bestchr,:),pmut); %mutate good chrom
xxchrom=xovsp(selchrom(bestchr+1:popsize-1,:),pcr); %xover 2nd chrom
xmchrom=mutbint(xxchrom,pmut);
newnew=initbp(1,m); %worst chrom,reproduce
chrom=[mutchrom;xmchrom;newnew];
numgen(t)=t;
end %end of evolution
% newchromosome=generation(:,1:m);
% [r f]=size(newchromosome);
% Out=[];MSE1=[];CX=[];CX2=[];MSE21=[];YHT2=[];
% for i=1:r
% selectedchrom=newchromosome(i,:);
% numbits=sum(newchromosome(i,1:m));
% cx=nou+noy+nl+numbits;
% CX=[CX;cx];
% pi=[];
% for r=1:m
% ind=selectedchrom(:,r);
% if ind==1
% pi=[pi mphi(:,r)];
% else
% pi=[pi mphi(:,r)*0];
% end
% end
% tet=pinv(pi'*pi)*pi'*yt;
% noterm=newchromosome(i,1:m);
% nterm=noterm';
% yht=pi*tet;
% msedata=mse(yt-yht);
% at=fix(log10(msedata))-1;a1=round(100*(msedata/10^at))/100;
% msedata=a1*10^at;
% MSE=msedata;
% MSE1=[MSE1;MSE];
% Out=[nterm tet];
% end
% teta=[tetaDE tetaGA]
BestMSE=[Fit(iBest) gene1(1,m+1)]
% BestTerm=Pop(:,iBest)
plot(numgenDE,fbestDE,'m-',numgen,sumfit,'g--')
xlabel('Generation');ylabel('Fitness Value');
RUN_TIME=cputime-tt;
4 commentaires
Réponse acceptée
Walter Roberson
le 23 Juin 2022
Modifié(e) : Walter Roberson
le 23 Juin 2022
The following fails because you did not provide initbp()
ui = [ -3.7287749e-001
1.5462074e-001
3.6396491e-001
-2.2540157e-001
3.4016394e-001
-4.2924476e-001
-1.2120792e-001
-2.3183057e-001
-3.4707882e-001
1.3099958e-001
-1.8362518e-001
4.5911185e-001
-1.3260108e-003
2.3860852e-001
-4.8724433e-001
1.0535319e-001
7.6451321e-002
3.0737882e-001
1.5496806e-001
3.7822840e-001
4.0237293e-001
-3.4776714e-001
-3.0741991e-001
2.9097565e-001
-4.3929490e-001
-1.1017281e-001
-2.0003415e-001
2.3417991e-001
-3.9579075e-001
2.9257549e-001
2.8272894e-001
3.2397693e-002
-2.4664766e-001
-4.2904541e-001
1.2580333e-001
-4.7531853e-001
-4.3795751e-001
-3.7038802e-001
-4.9385928e-002
1.7233584e-001
3.5611083e-001
-1.5545757e-003
-4.5121550e-001
-1.8616764e-001
1.4163117e-001
2.8638664e-001
-2.1084951e-001
-2.1319893e-003
3.1843462e-001
9.5128293e-002
3.6424609e-002
-1.6912686e-001
-8.8309774e-002
2.9400624e-001
-1.5679178e-001
-3.7393264e-002
-1.3217633e-001
1.7957062e-001
6.7779032e-002
1.5177685e-001
-8.8862517e-003
-1.0154177e-001
-2.2515791e-002
-4.3341238e-001
-8.8971513e-002
4.6909196e-001
2.8072103e-001
2.2901789e-001
2.6565062e-001
2.5658080e-001
3.4327032e-001
2.7015685e-001
4.7866093e-001
-3.8864056e-001
-1.0392306e-001
-7.9393492e-003
-2.4190538e-001
-4.6303378e-001
4.7438132e-001
2.2642598e-001
-3.5203517e-001
-3.5210544e-001
2.0482410e-001
-1.1900648e-001
-4.2358856e-001
-8.9159657e-002
-3.5700945e-001
2.9892012e-001
4.3024816e-001
-4.9528052e-001
1.5004058e-001
1.7853246e-001
-2.4637716e-001
3.4318091e-001
-2.0604123e-001
-4.7313962e-001
-4.0669284e-001
2.9788497e-001
2.1140263e-001
2.8340678e-001
1.2392264e-001
3.2541546e-001
-4.6497741e-001
-9.4525028e-002
-2.5033124e-001
-1.9103158e-002
3.8083913e-001
-2.1931518e-001
9.9142509e-002
-4.7377386e-001
-3.4480133e-001
3.3391087e-001
-3.0510663e-001
3.2978804e-001
-1.6192514e-001
1.7111655e-001
-4.4763368e-001
2.3430907e-001
-5.1822218e-004
4.4329994e-001
-2.1022758e-001
-1.2343957e-001
-3.8624464e-001
4.6485071e-001
-6.7488433e-002
-4.1543680e-001
2.1667641e-001
6.7873795e-003
-1.7191226e-001
2.5352414e-001
3.3600582e-001
-2.4628365e-001
3.4425017e-002
-6.4834767e-002
-3.4229585e-001
1.0048120e-001
4.3745076e-001
-3.9224142e-001
3.9998143e-001
5.0464991e-002
-7.2643313e-002
-3.4761917e-001
-2.5245400e-001
-5.2631206e-002
3.2783438e-002
-1.4534887e-001
2.7311481e-001
3.8168114e-001
2.3409209e-001
-9.3567794e-002
1.0417875e-001
1.4110505e-001
-3.7253316e-001
-3.8080031e-003
-1.8953064e-001
7.8574464e-002
4.4360894e-001
-7.3060391e-002
-4.6687156e-001
4.2943388e-001
4.2498225e-001
-1.4169308e-001
-2.4001088e-001
2.8685965e-001
1.1582243e-002
6.2527153e-002
1.8479376e-001
-4.0760340e-001
3.7257897e-001
4.4293788e-001
-4.0340576e-001
3.4588251e-001
4.0939561e-001
-4.8865886e-001
2.3682704e-002
1.5034527e-001
-1.1485465e-001
1.4930205e-001
2.6285348e-001
7.5685981e-002
1.3191915e-001
-2.2179600e-001
3.3983642e-001
-7.3165044e-002
1.3162238e-001
3.3346698e-001
-2.2981454e-001
-9.9199264e-002
5.4254592e-002
-5.6134643e-002
-4.0961588e-001
2.4438131e-001
-4.6738455e-001
-7.0257457e-002
-4.6273772e-001
4.7579775e-001
2.2339908e-002
4.0962998e-001
-1.1675208e-001
3.8445224e-001
-2.4498239e-001
4.0904984e-001
3.9456044e-001
-1.0148284e-001
1.2502000e-001
6.7597286e-002
3.9451189e-001
-2.8583420e-001
-4.9614084e-001
3.8058146e-001
-2.6487836e-001
-2.5513538e-001
1.4091704e-001
-1.9547552e-001
3.2562081e-001
3.8368871e-001
4.4537296e-001
-1.0920466e-001
3.0132006e-001
-3.4288671e-001
1.2516540e-001
1.9898535e-001
-4.1413107e-001
3.1179573e-002
3.8856652e-001
-2.3633399e-001
-2.6521447e-001
3.3965723e-001
-4.4599703e-003
-3.4763363e-001
-2.6923184e-001
1.5795396e-001
6.2948381e-002
-2.0817126e-001
1.2230498e-001
2.1590541e-001
-2.1926897e-001
-8.7726841e-002
-1.3779429e-001
2.8139202e-001
-3.6451285e-001
4.0206951e-001
-2.1036391e-001
-4.4207389e-004
2.8359039e-001
1.7706206e-001
-3.5018827e-001
1.9661948e-001
-3.7098696e-001
4.4594629e-001
3.8641185e-001
1.5004613e-002
1.7940875e-001
4.7679336e-001
-3.7453964e-001
2.5224305e-001
3.2705191e-001
2.8143033e-001
-3.0911733e-001
-7.1358659e-002
-4.8554325e-001
-1.7471451e-001
-3.6529797e-001
-4.9482538e-002
7.2274930e-002
2.9202328e-001
-8.0263390e-002
3.2536544e-002
4.2570439e-001
3.9908184e-001
4.4833398e-002
4.0112395e-001
-4.4817329e-001
3.0860252e-001
-1.6509542e-001
-2.7131688e-001
3.2240181e-001
-1.5176542e-001
-3.3452940e-001
-4.7186643e-001
4.5536977e-001
1.8029071e-001
3.6056211e-001
4.3909352e-001
1.8019396e-001
4.1742375e-001
-2.4330833e-001
3.8561813e-001
4.2004279e-001
-1.9993656e-001
-4.2660905e-001
2.6739873e-001
-4.1504778e-001
2.2876408e-001
-5.2110356e-002
1.5124866e-001
-3.3049813e-001
3.1449241e-002
1.3380086e-001
-4.8590400e-001
-2.9628665e-002
3.8632597e-001
-3.8597156e-001
-5.7458848e-002
1.5954789e-001
-2.0522703e-001
4.5036840e-001
1.9428625e-001
-2.9319347e-001
5.4762326e-002
3.7927848e-001
5.7857691e-002
2.5233406e-001
3.9490079e-001
3.4184417e-001
-3.6914332e-001
-3.1084607e-001
-3.4636017e-001
-4.7109793e-001
-4.9091510e-001
9.6451742e-002
1.0904932e-001
4.1892343e-001
2.3357444e-001
-1.9885305e-001
-4.4255210e-003
-2.4183698e-001
2.3285436e-001
-3.8323925e-001
2.4604156e-001
3.0978911e-001
2.4523368e-001
-1.6285669e-001
8.4324980e-002
-3.1047895e-002
-4.1273497e-001
3.2871743e-001
1.8594471e-001
-2.3267494e-001
4.6948356e-001
-3.1622580e-001
-2.0005908e-001
-8.8815109e-002
-2.6351081e-001
-3.0494468e-001
2.0537928e-001
-3.1945175e-001
2.2332918e-002
-2.0382818e-001
-3.7218112e-002
4.2523242e-001
-2.8411086e-001
-4.9899112e-001
4.0660624e-001
1.8004195e-001
1.4951805e-002
2.2068881e-002
-3.9708142e-001
4.9688049e-001
-1.4103118e-001
1.2524260e-001
-1.0663348e-001
-4.9233755e-001
4.5284934e-002
9.1103744e-003
-2.5322108e-001
-4.5460473e-001
3.4172974e-001
-4.5175008e-001
-1.8367950e-001
2.8341943e-001
4.7240006e-001
8.6464081e-002
2.7804385e-001
2.2769744e-001
1.5098680e-001
1.6461502e-001
4.3877980e-001
3.5081266e-002
-1.0155988e-001
1.7045766e-001
-5.9465528e-002
-3.6712606e-001
-6.0796252e-002
4.7642742e-002
-1.0486362e-001
-1.0172801e-001
2.5134865e-001
2.2350346e-002
-9.5669663e-003
-4.1132107e-001
-2.4914837e-001
-5.2440977e-002
1.3796090e-001
2.0944529e-001
4.9261981e-001
4.3219487e-001
-4.0777050e-001
4.5353997e-001
-3.3720434e-001
4.7045509e-001
9.7006665e-002
-2.5977288e-001
-4.2970484e-001
-1.9995878e-001
3.1354482e-001
-4.2329013e-001
-1.4552690e-001
-3.6798913e-001
-3.4182037e-001
-4.3785295e-001
2.0184348e-001
-4.1351830e-001
1.1678724e-001
-3.2622842e-001
1.5140106e-001
-1.3035940e-003
-2.1548931e-001
3.3055974e-001
3.1835988e-001
4.3816940e-001
-4.9967390e-001
1.4038922e-001
-4.9264431e-001
-3.9357889e-001
-3.9320559e-001
-1.3289096e-001
-2.6039230e-001
-1.5385980e-001
-2.5038016e-001
-1.1293569e-001
-7.8961630e-002
1.4007729e-001
2.8755297e-001
-2.3000558e-001
3.4398236e-001
2.4046838e-001
3.2610195e-001
-3.1780775e-001
-4.3456382e-001
1.1035026e-001
2.0155281e-001
-3.8838229e-001
-4.0417627e-001
9.7833772e-002
3.1223261e-001
3.1457817e-001
-4.1056303e-001
2.3127871e-001
4.0385658e-001
-4.7768213e-002
-4.2931232e-001
-2.5872184e-001
2.3186506e-001
-4.5950838e-001
-7.5477237e-002
4.0215087e-002
4.5382784e-001
-2.9109426e-001
-3.8366790e-001
1.4622007e-001
-3.9158891e-001
4.8349692e-001
-2.5165596e-001
1.0635572e-001
3.1669527e-001
3.3005421e-001
-1.0960668e-002
2.6072826e-001
4.1510794e-001
4.0097496e-001
-2.8576249e-001
4.7061106e-002
2.8470948e-001
-3.0555891e-001
2.4688915e-001
-2.4441776e-002
8.3258581e-002
-2.3945192e-001
-4.1517754e-001
-2.0186715e-001
4.1712936e-001
-2.9481807e-002
-2.3053245e-001
2.6297048e-001
2.7217211e-001
-4.7870038e-001
3.7998619e-001
2.9817018e-001
-1.7583476e-001
1.6904442e-001
-2.0370614e-001
4.2995204e-001
-2.1803957e-001
-3.3112003e-001
2.4516660e-001
-2.2865648e-002
1.5344454e-001
4.6657593e-001
-1.8697251e-001];
save('ui.dat', 'ui', '-ascii');
% %this program for comparison between MGA & DE
clear;
%load input-output data
tt=cputime;
load ui.dat; %load input
u=ui(:,1);
mx=2;
y(1:mx+1)=0;
for t=mx+1:500
% for i=1:500
% e(i)=(rand-0.5)*0.002;
% end
% y(t)=0.797*y(t-1)-0.237*y(t-3)+0.441*u(t-1)+0.105*y(t-4)*u(t-4)...
% +0.333*u(t-3)*u(t-5);%+e(t);
y(t)=0.3*u(t-1)*u(t-2)-0.8*y(t-2)*u(t-1)*u(t-1)-0.1*u(t-1)*u(t-1)*u(t-2);%+e(t);
end
y=y';%load output
nodata=500;
noy=2; nou=2;
nl=3; start=max(nou,noy)+1;
m=nou+noy;
mphi=[]; %form regressor for ARMAX model
for i=start:nodata
phi=[];
if nl==1
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
end
if nl==2
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
for j=1:noy
for l=j:noy
phi=[phi y(i-j)*y(i-l)];
end
for k=1:nou
phi=[phi y(i-j)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi u(i-k)*u(i-n)];
end
end
end
if nl==3
for j=1:noy
phi=[phi y(i-j)];
end
for k=0:nou
phi=[phi u(i-k)];
end
for j=1:noy
for l=j:noy
phi=[phi y(i-j)*y(i-l)];
end
for k=1:nou
phi=[phi y(i-j)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi u(i-k)*u(i-n)];
end
end
for j=1:noy
for l=j:noy
for p=l:noy
phi=[phi y(i-j)*y(i-l)*y(i-p)];
end
for k=1:nou
phi=[phi y(i-j)*y(i-l)*u(i-k)];
end
end
for k=1:nou
for n=k:nou
phi=[phi y(i-j)*u(i-k)*u(i-n)];
end
end
end
for k=1:nou
for n=k:nou
for q=n:nou
phi=[phi u(i-k)*u(i-n)*u(i-q)];
end
end
end
end
mphi=[mphi;phi]; %end regressor
end
[temp m]=size(mphi);
%DIFFERENTIAL EVOLUTION
%CONTROL PARAMETERS
D=m; %dimension of problem
NP=30; %Population size
F=0.9; %differentiation constant
CR=0.5; %crossover rate
GEN=100; %generations
L=0.0; %low boundary constraint 1
H=1.0; %high boundary constraint 1
% L2=4.1; %low boundary constraint 2
% H2=5.8; %high boundary constraint 2
%**********************************%
% Algorithm Variables
%**********************************%
X=zeros(D,1); %trial vctor X1
% X2=zeros(D,1); %trial vctor X2
Pop=zeros(D,NP); %population X1
% Pop2=zeros(D,NP); %population X1
Fit=zeros(1,NP); %fitness population
iBest=1; %index of best solution
r=zeros(3,1); %randomly selected indices
%**********************************%
% Create Population
%**********************************%
%Initialize random number of generator
%rand('state',sum(100*clock));
for j=1:NP %Initialize each individual
Pop(:,j)=L+(H-L)*rand(D,1); %within b.constraint 1
% Pop2(:,j)=L2+(H2-L2)*rand(D,1); %within b.constraint 2
[Fit(1,j), tetaDE]=selectmodel2(Pop(:,j),NP,D,mphi,start,nodata,y,u); %evaluate fitness 1
end
%**********************************%
% Optimization
%**********************************%
fbestDE=[];
for g=1:GEN %for each generation
for j=1:NP %for each individual
%choose 3 random individuals from population
%mutually different and different from j
r(1)=floor(rand()*NP)+1;
while r(1)==j
r(1)=floor(rand()*NP)+1;
end
r(2)=floor(rand()*NP)+1;
while (r(2)==r(1))||(r(2)==j)
r(2)=floor(rand()*NP)+1;
end
r(3)=floor(rand()*NP)+1;
while (r(3)==r(2))||(r(3)==r(1))||(r(3)==j)
r(3)=floor(rand()*NP)+1;
end
%create trial individual
%in which at least one parameter is changed
Rnd=floor(rand()*D)+1;
for i=1:D
if (rand()<CR)||(Rnd==i)
X(i)=Pop(i,r(3))+F*(Pop(i,r(1))-Pop(i,r(2)));
% X2(i)=Pop2(i,r(3))+F*(Pop2(i,r(1))-Pop2(i,r(2)));
else
X(i)=Pop(i,j);
% X2(i)=Pop2(i,j);
end
end
%verify b.constraint
for i=1:D
if (X(i)<L)||(X(i)>H)
X(i)=L+(H-L)*rand();
end
% if (X2(i)<L2)||(X2(i)>H2)
% X2(i)=L2+(H2-L2)*rand();
% end
end
%select the best individual
%between trial and current ones
%calculate fitness of trial vector
[f, tetaDE]=selectmodel2(X,NP,D,mphi,start,nodata,y,u);
%if trial is better or equal than current
% if f<=Fit(j)
% Pop1(:,j)=X1; %replace current by trial
% Pop2(:,j)=X2;
% Fit(j)=f;
% %if trial is better than the best
% if f<=Fit(iBest)
% iBest=j; %update the best's index
% end
% end
if f<=Fit(j)
Pop(:,j)=X; %replace current by trial
% Pop2(:,j)=X2;
Fit(j)=f;
%if trial is better than the best
if f<=Fit(iBest)
iBest=j; %update the best's index
end
end
end
fbestDE=[fbestDE Fit(iBest)];
numgenDE(g)=g;
end
% %MGA
%CONTROL PARAMETERS
evol=100;%input('# of evolution >:');
popsize=30;%input('population size >:');
pcr=0.8;%input('crossover probability=');
pmut=0.01;%input('mutation probability=');
chr=initbp(popsize,m); %initialize population
chrom=chr;
mseGA=[];
for t=1:evol %start evolution
gen1=[];
for i=1:popsize %model structure chrom's=1
mmphi=[];
krom=chrom(i,:);
for b=1:m
index=krom(:,b);
if index==1
mmphi=[mmphi mphi(:,b)];
end
end
[N, g]=size(mmphi);
yt=y(start:nodata);ut=u(start:nodata);
y2=sum(yt.*yt);
teta=pinv(mmphi'*mmphi)*mmphi'*yt; %find teta by LSE
num=0;
for c=1:g
tetahat=teta(c,:);
if tetahat<0.001
num=num+0.1; %penalty to added noise
end
end
yhat=mmphi*teta;
error=sum((yt-yhat).*(yt-yhat))+num;
gen=[krom error];
gen1=[gen1;gen];
pm=max(gen1(:,m+1));
fitness=pm-gen1(:,m+1);
% bestftn=pm-fitness-num;
gene1=[gen1 fitness];
end
generation=sortrows(gene1,[m+1]);
% mseGA=[mseGA gene1(:,m+1)];
sse=generation(:,m+1);
bestchr=0.0; secondchr=0.0;
for i=1:popsize
sume2=sse(i,:);
if sume2<1.0
bestchr=bestchr+1; %good chromosomes
else
secondchr=secondchr+1; %2nd class chromosomes
end
end
sumerror2(t)=sum(gene1(:,m+1));
fitnv=gene1(:,m+2); %assign fitness
sumfit(t)=sum(fitnv);
%newchr=selsus(fitnv,popsize);
newchr=seltour(fitnv,popsize,2);
selchrom=chrom(newchr,:);
%selchrom=generation(:,1:m);
mutchrom=mutbint(selchrom(1:bestchr,:),pmut); %mutate good chrom
xxchrom=xovsp(selchrom(bestchr+1:popsize-1,:),pcr); %xover 2nd chrom
xmchrom=mutbint(xxchrom,pmut);
newnew=initbp(1,m); %worst chrom,reproduce
chrom=[mutchrom;xmchrom;newnew];
numgen(t)=t;
end %end of evolution
% newchromosome=generation(:,1:m);
% [r f]=size(newchromosome);
% Out=[];MSE1=[];CX=[];CX2=[];MSE21=[];YHT2=[];
% for i=1:r
% selectedchrom=newchromosome(i,:);
% numbits=sum(newchromosome(i,1:m));
% cx=nou+noy+nl+numbits;
% CX=[CX;cx];
% pi=[];
% for r=1:m
% ind=selectedchrom(:,r);
% if ind==1
% pi=[pi mphi(:,r)];
% else
% pi=[pi mphi(:,r)*0];
% end
% end
% tet=pinv(pi'*pi)*pi'*yt;
% noterm=newchromosome(i,1:m);
% nterm=noterm';
% yht=pi*tet;
% msedata=mse(yt-yht);
% at=fix(log10(msedata))-1;a1=round(100*(msedata/10^at))/100;
% msedata=a1*10^at;
% MSE=msedata;
% MSE1=[MSE1;MSE];
% Out=[nterm tet];
% end
% teta=[tetaDE tetaGA]
BestMSE=[Fit(iBest) gene1(1,m+1)]
% BestTerm=Pop(:,iBest)
plot(numgenDE,fbestDE,'m-',numgen,sumfit,'g--')
xlabel('Generation');ylabel('Fitness Value');
RUN_TIME=cputime-tt;
function [f teta]=selectmodel2(X,popsize,m,mphi,start,nodata,y,u)
gen1=[];
for i=1:popsize %model structure chrom's=1
mmphi=[];
krom=X;
for b=1:m
index=krom(b,:);
if index>0.5
mmphi=[mmphi mphi(:,b)];
end
end
[N g]=size(mmphi);
yt=y(start:nodata);ut=u(start:nodata);
% y2=sum(yt.*yt);
teta=inv(mmphi'*mmphi)*mmphi'*yt; %find teta by LSE
num=0;
for c=1:g
tetahat=teta(c,:);
if tetahat<0.001;
num=num+0.1; %penalty to added noise
end
end
yhat=mmphi*teta;
mse=sum((yt-yhat).*(yt-yhat))+num;
%whos
gen=[krom.' mse];
gen1=[gen1, gen];
pm=max(mse);
% f=mse;
f=pm-gen1(:,m+1);
% gene1=[gen1 f];
end
end
0 commentaires
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Matrices and Arrays dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
Vous pouvez également sélectionner un site web dans la liste suivante :
Amériques
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
Asie-Pacifique
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)