通俗解释matlab之遗传算法程序部分(二)

(1)程序怎么开始

从哪里开始程序比较好了？直接先主函数吧，然后再分着说：

%-------------函数说明----------------

%             主函数      

%---------------------------------------

function main()

clear

clc

popsize = 100;     %种群大小

chromlength = 10;  %二进制编码长度

pc = 0.6;          %交叉概率

pm = 0.001;        %变异概率

pop = initpop(popsize,chromlength);   %初始种群

 

  for i=1:100

     [objvalue] = cal_objvalue(pop);     %计算适应度值（函数值）

      fitvalue = objvalue;

     [newpop] = selection(pop,fitvalue);  %选择操作

     [newpop] = crossover(newpop,pc);     %交叉操作

     [newpop] = mutation(newpop,pm);      %变异操作

      pop = newpop;                %更新种群

      [bestindividual,bestfit]=best(pop,fitvalue);%寻找最优解

      x2 = binary2decimal(bestindividual);

   

      x1 = binary2decimal(newpop);

      [y1] = cal_objvalue(newpop);    

      if  mod(i,10)==0

      figure;   

      fplot('10*sin(5*x)+7*abs(x-5)+10',[0 10]);

      hold on;

      title(['迭代次数为 n=' num2str(i)]);

      plot(x1,y1,'*');

      end

  end

      fprintf('the best X is  --->>%5.2f\n',x2);

      fprintf('the best Y is  --->>%5.2f\n',bestfit);

好了，这就是主函数，最后运行这个基本上什么都有了，包括图行什么的都画了，下面在分着说吧。。。

（2）关于二进制种群怎么生成

        从上面的主程序可以看到，我们设置了100个个体（也就是100个初始化x值），二进制编码长度为10位，那么先开始自然是生成这100个随机个体了，注意是随机的，不能使相同的，也就是说生成100个不同的二进制编码组合，像11000 10101这样的，生成100组。怎么办，直接上程序：

%-------------函数说明----------------

%    初始化种群大小

%       输入变量：

%               popsize:种群大小

%               chromlength：染色体长度--》转化的二进制长度

%       输出变量：

%               pop：种群

%---------------------------------------

function pop = initpop(popsize,chromlength)

pop = round(rand(popsize,chromlength));

很简单，一句话搞定，说一下，关于rand（m,n）用法，就是生成m行n列的0~1之间随机数，比如rand(3,4)为：

>> rand(3,4)

ans =

    0.6028    0.1174    0.4242    0.2625

    0.7112    0.2967    0.5079    0.8010

0.2217    0.3188    0.0855    0.0292

round什么意思？四舍五入，简单，这样上面就变成了什么了，round（ans）：

>> round(ans)

ans =

     1     0     0     0

     1     0     1     1

     0     0     0     0

这样当popsize=100，chromlength=10时就生成了对应的种群个体了。

（3）如何在把二进制返回对应的十进制

知道了二进制怎么返回对应范围内的x值呢？前面讲了算法怎么构造了，特别要注意的是二进制的位数以及转化完后的x值范围，程序如下：

%-------------函数说明----------------

%     二进制转化十进制函数

%       输入变量：

%               二进制种群

%       输出变量：

%               十进制数值

%---------------------------------------

function pop2 = binary2decimal(pop)

[px,py]=size(pop);

for i=1:py

    pop1(:,i) = 2.^(py-i).*pop(:,i);

end

%sum(.,2)对行求和，得到列的向量

temp = sum(pop1,2);

pop2 = temp*10/1023;

输入的为100组0、1编码的二进制，输出的是x值。开始取一下种群大小size(pop)，显然这里py=10了，接着对每一位求和，就是pop1(:,i) = 2.^(py-i).*pop(:,i);这里省略用了冒号，什么意思了，就是对所有的行都是这个操作，冒号意思就是从1到100了。那么就其中某个个体比如第一个吧，假设为11001 10010，那么先进行这么一个操作后就是什么了？是不是就是对应为的0或1乘以2的对应次幂了，若果是1就管用，是0就不管用，那么这个值就是：（2^0）*1+（2^1）*1+0+0+(2^4)*1+...,这样就算出了一个值了，因为是10位编码，所以这个数介于0~2^9即0~1023。那么最大为多少？1023吧。temp = sum(pop1,2);对行求和吧，2表示对行，1表示对列，像

ans =  

     1     0     0     0

     1     0     1     1

     0     0     0     0

则>> sum(ans,1)

ans =

     2     0     1     1    

而>> sum(ans,2)

ans =

     3

     2

     2

明白了吧，这样temp就变成了有100行1列的值，且每个值都在0~1023之间变化，最后一行不解释了，就是把它转化为100组0~10之间的数值了。

（4）计算适应度函数值

这里也就是根据100组x计算对应的100组y的值，简单，根据上面的x值带入就可以了：

%-------------函数说明----------------

%     计算函数目标值

%       输入变量：

%              二进制数值 

%       输出变量：

%               目标函数值

%---------------------------------------

function [objvalue]=cal_objvalue(pop)

x = binary2decimal(pop);

%转化二进制数为x变量的变化域范围的数值

objvalue = 10*sin(5*x)+7*abs(x-5)+10;

（5）如何选择新的个体

       上面所有个体的函数值都计算出来了，存在objvalue中，此时它是不是也是100组y值呀，恩。那么对于现有的随机生成的100组x，怎么来在选择100组新的更好的x呢？这里我们把选择放在了交叉与变异之间了，都可以。如何选择，就要构造概率的那个轮盘了，谁的概率大，是不是选择的个体就会多一些？也就是现在的选择就是100中选100个，最后出现的结果就是以前的100个中最优的x有一个的话，选择完后，可能就变成了5个这个x了，多余的4个是不是相当于顶替了以前的不好的4个x值，这样才能达到x总数100不变呀。

%-------------函数说明----------------  

%       输入变量：

%                pop  ：  二进制种群

%                fitvalue  ：  适应度值

%       输出变量：

%                newpop：选择以后的二进制种群

%---------------------------------------

function [newpop] = selection(pop,fitvalue)

%构造轮盘

[px,py]=size(pop);   

totalfit = sum(fitvalue);

p_fitvalue = fitvalue/totalfit; 

p_fitvalue = cumsum(p_fitvalue);%概率求和排序

%-------

ms = sort(rand(px,1));%从小到大排列

fitin = 1;

newin = 1;

while newin<=px

     if (ms(newin))<p_fitvalue(fitin)

         newpop(newin,:)=pop(fitin,:);

         newin=newin+1;

     else fitin=fitin+1;

     end

 end

这一部分可能不太好理解。自己好好想想吧。前三句，求出每个个体被选择的概率吧，第四句，求和排序是干什么的呢？比如现在假设一个概率可能是：

p_fitvalue =     

    0.0816

    0.0408

    0.1020

    0.1224

    0.0408

    0.1429

    0.1837

    0.0204

    0.0816

    0.1837

那么执行这一句p_fitvalue = cumsum(p_fitvalue)后就变为：

p_fitvalue =

    0.0816

    0.1224

    0.2245

    0.3469

    0.3878

    0.5306

    0.7143

    0.7347

    0.8163

    1.0000

也就是后面的值加到前面去在替换，这样做有什么意义呢？看后面吧。ms = sort(rand(px,1))这一句，生成100个0~1的随机概率数，然后在从小到大排序，100次，好理解，轮盘赌不得赌100次吗？为什么排序呢？这和前面的cumsum求和了有关。好了真正的选择开始了，while里面，先对小的概率进行选择，

if (ms(newin))<p_fitvalue(fitin)，就选择一次，个体就多了一个。符合，newin概率的那个排序往上走，再符合，就再选择一次，直到不符合了，fitin就往上走，这样p_fitvalue(fitin)这部分是不是变大了一点，变大了再看看这个值能被选择几次。其实每一次循环就是对于某一个个体来说，看能选择几次。比如假设fitin=9时，此时假设p_fitvalue(9)=0.5，被选择完了，下一次时fitin=10，我们假设离谱点，p_fitvalue(10)=0.7的话，那么在下一次去ms（newin）的话(假设现在已经选择了newin=30个个体了)，是不是一定为选择呀，然后第31个个体就是p_fitvalue(10)这个个体了，在32个，因为ms（newin）也是从小到大排序的，所以这次选择只是增加了一点点吧，假设为0.52（肯定大于0.5吧，不然也不会进来），0.52<0.7，又选择p_fitvalue(10)这个个体吧，这样一直下去，33,34，，，好，来到了第50个个体假设。此时ms（50）=0.71了，不行了，是不是不行了，fitin是不是加1变为11了，如果p_fitvalue(11)较p_fitvalue(10)变化不大，比如只变为了0.72，那么第11个个体的被选择的概率是不是只是0.02呀，这么小点，也只能是就选择那么2回就下一个了吧，而我们看看第10个个体，从0.5变到了0.7，它的概率是不是0.2呀，哇，0.2呀，好大呀，难怪能从第32个个体一直选择到第50个个体呀。好了到这里应该能明白了吧，再不明白就没招了，自己好好想吧。那么这么一系列循环后，是不是100个个体就被选择出来了呀，想想选择的结果是什么样子的呢？100个个体是不是会出现为：第2个第2个，第3个第3个第3个，第5个第5个第5个第5个第5个，第6个，第10个第10个第10个，等等等，总共100个个体吧（那些重复的数是在以前种群中的编号）。好多呀，应该讲明白了选择这块。

（6）怎么交叉

选择完后是不是就有新的相对优秀的100个个体了，那么对它们我们在来进行交叉变异才能形成新的个体，毕竟这100个个体还是先前100个个体里面挑好的挑出来的吧。

%-------------函数说明----------------    

%       输入变量：

%             pop:二进制的父代种群数

%             pc：交叉的概率

%       输出变量：

%             newpop：交叉后的种群数   

%---------------------------------------

function [newpop]=crossover(pop,pc)

[px,py]=size(pop);

newpop = ones(size(pop));

for i=1:2:px-1

    if(rand<pc)

        cpoint = round(rand*py);

        newpop(i,:) = [pop(i,1:cpoint),pop(i+1,cpoint+1:py)];

        newpop(i+1,:) = [pop(i+1,1:cpoint),pop(i,cpoint+1:py)];

    else

        newpop(i,:)=pop(i,:);

        newpop(i+1,:)=pop(i+1,:);

    end

end

这里涉及到了一个概率，就是交叉的概率pc，因为交叉不是一定发生，所以用一个概率来限制，部分还是要保持原来不变的。程序开始先生产一个与原种群同大小的全1矩阵newpop，然后循环for i=1:2:px-1，想想为什么间距是2呢？简单，每两个交叉吧，就像自然界伟大的规律一样，当然也有例外，不考虑吧。然后生产随机0~1数rand，在一比较就知道是不是该交叉，这里的pc = 0.6，那么想想每次rand的话，由于随机性，是不是它产生小于0.6的概率就是0.6呀。好了选择了交叉，就执行下面操作，其中还有个交叉点cpoint 计算，能明白吧，不解释，后面两句也能看懂。如果概率比0.6大的话（也就是0.4的概率rand出了0.6~1之间的数），就不选择交叉，直接赋值就行了吧。这样知道结束，有个问题，万一px为奇数怎么办？是不是相当于有一个光棍呀，是呀，咋办的呢？我也不知道，但是我们这里为100，还好是刚好1配1呀，很人性，大家编程序也要人性化一点呀，虽然那是数字，没有生命，但是体现了你的那个什么，道德修养，对世界的爱吧~_~。

（7）关于变异

     好了交叉完后就生成了新一代的种群了吧，那么在对这个种群再一次进行变异，实际情况下交叉变异应该是同时发生的吧，这里先后无所谓了。那么怎么变异呢？

%-------------函数说明----------------  

%       输入变量：

%                 pop  ：  二进制种群

%                 pm  ：  变异概率

%       输出变量：

%                newpop  ：  变异以后的种群

%---------------------------------------

function [newpop] = mutation(pop,pm)

[px,py] = size(pop);

newpop = ones(size(pop));

for i=1:px

    if(rand<pm)

        mpoint = round(rand*py);

        if mpoint<=0

            mpoint=1;

        end

        newpop(i,:) = pop(i,:);

        if newpop(i,mpoint)==0

            newpop(i,mpoint)=1;

        else newpop(i,mpoint)=0;

        end

    else

        newpop(i,:)=pop(i,:);

    end

end

个人感觉变异更简单，就是找个点（高级一点的找多个点），把它的值（0或1）取反不就得了嘛，程序就像上面的那样，讲到这里了感觉真不需要解释了，略过。

（8）最后一点了

至此，主体部分算是完了吧，剩下的就是显示、处理部分了，回到主程序去看吧，接下来就是更新种群吧，然后再变换十进制去，提取最好的y值和对应的x值，这部分还有个求最优适应度的函数：

%-------------函数说明---------------- 

%       输入变量：  

%                pop  :  种群

%                fitvalue  ：  种群适应度 

%

%       输出变量：

%                bestindividual  ： 最佳个体（二进制个体）

%                bestfit     ：   最佳适应度值

%---------------------------------------

function [bestindividual,bestfit]=best(pop,fitvalue)

[px,py]=size(pop);

bestindividual = pop(1,:);

bestfit = fitvalue(1);

for i=2:px

    if fitvalue(i)>bestfit

        bestindividual = pop(i,:);

        bestfit = fitvalue(i);

    end

end

感觉也没有什么号解释的，就是比较呗，最优值在bestfit里面，如果有比它大的，就更新，否则略过，不解释。

        最后是在matlab环境下画出来，if  mod(i,10)==0求余操作吧，也就是每迭代10代我画一次图，观察结果吧，后面的fprintf，打印出来，很简单，自己看吧。。。

至此，整个遗传算法的大致原理讲解完了，程序的各部分代码都在这里了，把每个小部分都编写一个m文件存起来再一起运行就ok了。所有的程序排列以及这部分实验结果在下一个里面再贴出来吧~_~

第一部分：通俗解释matlab之遗传算法求最优值(一)

第三部分：通俗解释matlab之遗传算法程序汇总与结果显示(三)