iga.m

人工免疫算法基于遗传MATLAB代码很有用哦

function iga(prompt,def)


clear all
%clc
prompt={'适应度函数Func:','边界范围Bounds:','自变量个数Nvar:','种群总数N:','终止代数R:','求最大值（0）or最小值（1）:',...
        '交叉概率pcross(0.5,0.9):','变异概率pmutation(0.05,0.2):','选择参数(编码,求解精度):'};
def={'x + 10*sin(5*x)+7*cos(4*x)','[-100,100]','2','60','80','0','0.7','0.05','[0 1e-6]'};
dlgTitle='input';%'参数输入';
lineNo=1;
answer=inputdlg(prompt,dlgTitle,lineNo,def);
answer=char(answer);
FUN='0.5+(sin(sqrt(x.*x+y.*y)).* sin(sqrt(x.*x+y.*y))-0.5)./((0.001*(x.*x+y.*y)+1.0).^2)';%(answer(1,:));
bounds=str2num(answer(2,:));
nvar=str2num(answer(3,:));
popsize=str2num(answer(4,:));
gen=str2num(answer(5,:));
tag=str2num(answer(6,:));
pcross=str2num(answer(7,:));
pmutation=str2num(answer(8,:));
options=str2num(answer(9,:));

if find((bounds(1)-bounds(2))>0)
   error('数据输入错误,请重新输入(bounds(1)<bounds(2)):');
end
tic;
bits=[];
precision=options(2);%由求解精度确定二进制编码长度
bits=ceil(log2((bounds(:,2)-bounds(:,1))' ./ precision));%由设定精度划分区间
[Pop]=initpop(popsize,bits+1);        %初始化种群
[m,len]=size(Pop);                    % m行，len列
memo=round(m/5);
pop_memo=zeros(memo,len-1);

pm0=pmutation;
BestPop=zeros(gen,len-1);Trace=zeros(gen,4);%分配初始解空间
it=1;
while it<=gen
    for j=1:m
        if nvar==1
           x=b2f(Pop(j,1:len-1),bounds);                    %计算初始种群的十进制转换
         elseif nvar==2 
           x=b2f(Pop(j,1:(len-1)/2),bounds);
           y=b2f(Pop(j,1+(len-1)/2:len-1),bounds);
         else printf('the nvar>3!\n');
       end
          Pop(j,len)=eval(FUN);         %计算适应度
    end
    
    [MaxValue,Index]=max(Pop(:,len));
    BestPop(it,:)=Pop(Index,1:len-1);
    Trace(it,1)=MaxValue;                       % the best
    Trace(it,2)=mean(Pop(:,len));               % mean
    Trace(it,3)=min(Pop(:,len));                % min
    if nvar==1          % the best ch
          Trace(it,4)=b2f(BestPop(it,1:len-1),bounds);             
         elseif nvar==2 
          Trace(it,4)=b2f(BestPop(it,1:(len-1)/2),bounds) + i*(b2f(BestPop(it,1+(len-1)/2:len-1),bounds));
         else printf('the nvar>3!\n');
       end
    
   
    %群体更新
%    [qiw]=qiwang1(Pop,0.5);%计算期望值  值0.5影大
     [qiw]=qiwang(Pop,0.35);%计算期望值  值0.35影大
   
    NewPop1=zeros(m,len);
    NewPop1(:,1:len-1)=Pop(:,1:len-1);
    NewPop1(:,len)=qiw';  
    % 在matlab中对矩阵的值按某一列升序排序。
     result = sortrows(NewPop1,len);
     result=flipud(result);
     Pop1=result(1:(m-memo),1:len-1);
     pop_memo=result(1:memo,1:len-1);
        
    [selectpop]=SelectChrom(FUN,Pop1,bounds,nvar);% 选择
    [CrossOverPop]=CrossOver(selectpop,pcross); % 交叉
    [NewPop]=Mutation(CrossOverPop,pmutation);  % 变异
    Pop(:,1:len-1)=[NewPop;pop_memo];
       
    pmutation=pm0+(it^4)*(pcross/2-pm0)/(gen^4); %随着种群向前进化，逐步增大变异率
    p(it)=pmutation;
    it=it+1;
end
t=1:gen;
plot(t,Trace(:,1)');
title('函数优化的免疫遗传算法');xlabel('进化世代数(gen)');ylabel('每一代最优适应度(maxfitness)');
[MaxFval,I]=max(Trace(:,1));
X=Trace(I,4);
hold on;  plot(I,MaxFval,'*');
text(I+5,MaxFval,['FMAX=' num2str(MaxFval)]);
str1=sprintf('进化到 %d 代 ,自变量为 %s 时,得本次求解的最优值 %f\n对应染色体是：%s\n use times: %5.2f',...
             I,num2str(X),MaxFval,num2str(BestPop(I,:)),toc);
disp(str1);



%*****************************************
%                子函数
%*****************************************

%初始化种群，采用二进制编码
function [pop]=initpop(popsize,bits)
len=sum(bits);
pop(1,:)=zeros(1,len);%The whole zero encoding
for i=2:popsize-1
    pop(i,:)=round(rand(1,len));
end
pop(popsize,:)=ones(1,len);%The whole one encoding 

%解码
function [fval] = b2f(bval,bounds,nvar)
% fval   - 表征各变量的十进制数
% bval   - 表征各变量的二进制编码串
% bounds - 各变量的取值范围
% bits   - 各变量的二进制编码长度
if nargin == 2,
   nvar=1;
end
bits=size(bval,2)/nvar;
scale=(bounds(:,2)-bounds(:,1))'./(2.^bits-1); %The range of the variables
numV=size(bounds,1);
cs=[0 cumsum(bits)]; 
for i=1:numV
  a=bval((cs(i)+1):cs(i+1));
  fval(i)=sum(2.^(size(a,2)-1:-1:0).*a)*scale(i)+bounds(i,1);
end 

%选择操作
function [selectpop]=SelectChrom(FUN,pop,bounds,nvar)%计算各个体的适应度并采用轮盘赌进行选择
[m,n]=size(pop);
for i=1:m
    if nvar==1
           x=b2f(pop(i,1:n),bounds);                    %计算初始种群的十进制转换
         elseif nvar==2 
           x=b2f(pop(i,1:n/2),bounds);
           y=b2f(pop(i,1+n/2:n),bounds);
         else printf('the nvar>3!\n');
    end
    fit(i)=eval(FUN);         %计算适应度
end
selectprob=fit/sum(fit);%选择概率
prob=cumsum(selectprob);%累计选择概率
sumprob=[0 prob];
for i=1:m
    selectpop(i,:)=pop(length(find(rand>=sumprob)),:);
end    
    
%交叉操作
function [NewPop]=CrossOver(OldPop,pcross)%OldPop为父代种群，pcross为交叉概率
[m,n]=size(OldPop);
r=rand(1,m);
y1=find(r<pcross);
y2=find(r>=pcross);
len=length(y1);
if len>2&mod(len,2)==1%如果用来进行交叉的染色体的条数为奇数，将其调整为偶数
    y2(length(y2)+1)=y1(len);
    y1(len)=[];
end
if length(y1)>=2
   for i=0:2:length(y1)-2
       [NewPop(y1(i+1),:),NewPop(y1(i+2),:)]=EqualCrossOver(OldPop(y1(i+1),:),OldPop(y1(i+2),:));
   end     
end
NewPop(y2,:)=OldPop(y2,:); 

function [children1,children2]=EqualCrossOver(parent1,parent2)
%采用均匀交叉 例：
%父1：0 1 1 1 0 0 1 1 0 1 0
%父2：1 0 1 0 1 1 0 0 1 0 1
%掩码：0 1 1 0 0 0 1 1 0 1 0
%交叉后新个体：
%子1：1 1 1 0 1 1 1 1 1 1 1 
%子2：0 0 1 1 0 0 0 0 0 0 0
L=length(parent1);
hidecode=round(rand(1,L));%随机生成掩码，如hidecode=[0 1 1 0 0 0 1 1 0 1 0];
children1=zeros(1,L);
children2=zeros(1,L);
children1(find(hidecode==1))=parent1(find(hidecode==1));%掩码为1，父1为子1提供基因
children1(find(hidecode==0))=parent2(find(hidecode==0));%掩码为0，父2为子1提供基因
children2(find(hidecode==1))=parent2(find(hidecode==1));%掩码为1，父2为子2提供基因
children2(find(hidecode==0))=parent1(find(hidecode==0));%掩码为0，父1为子2提供基因 

%变异操作
function [NewPop]=Mutation(OldPop,pmutation)
[m,n]=size(OldPop);
r=rand(1,m);
position=find(r<=pmutation);
len=length(position);
if len>=1
   for i=1:len
       k=unidrnd(n,1,1); %设置变异点数，一般设置1点
       for j=1:length(k)
           if OldPop(position(i),k(j))==1
              OldPop(position(i),k(j))=0;
           else
              OldPop(position(i),k(j))=1;
           end
       end
   end
end
NewPop=OldPop;

%计算抗体浓度及期望值—距离法
function [qiw]=qiwang(pop,det)
  [nm,len]=size(pop);
 det=1/(1+sqrt(det*nm));
for t=1:nm
     cv=0;
    for j=1:nm
        d=0;
        for k=1:(len-1)
           if(pop(t,k)~=pop(j,k))
               d=d+1;
           end
        end
        hn=1/(1+sqrt(d));
%      Aff=1/(1+hn);  
      if (hn>=det)
          cv=cv+1;
      end
  end
  % or qiw(t)=pop(t,len)*exp(-cv/nm)
 qiw(t)=pop(t,len)*nm/cv;
end