第一次使用matlab建模,所以摸索了很长久这个软件,后来发现其实就可以把它当做以前使用的VC++或者是eclipse编辑器就可以,没有我想的那么困难。
遗传算法步骤
我首先是根据老师给的代码,读了每一个步骤的每一个算法,然后修改了代码中的函数部分得到了作业的结果,但是我对最后结果所呈现的图像并不是非常理解,读了一些博客,才有了自己的一点点理解。
遗传算法中的每一条染色体,对应着遗传算法的一种基因组编码也就是一种解决方案,我们使用适应性函数去衡量这个解决方案的优劣。从一种编码到它的适应度形成了一种映射,相当于在多元函数中求最优解。可以想象,多维曲面里面有数不清的“山峰”,每一个”小山“对应的是函数自变量的某个区间,"山峰"对应的就是在某个区间内的最优解。这个多维曲面里面海拔最高的"山峰"就是全局的最优解。遗传算法的任务就是尽量去爬到最高峰,而不是陷落在 一些小山峰里。
作业中的问题描述:
我使用matlab所呈现的图像:
遗传算法的实现过程
遗传算法实际是一种对潜在问题进行"数字化"编码的方案,去建立表现型和基因型的映射关系。
-
个体编码
本题中使用无符号二进制整数来表示
-
初始群体的产生
(1)首先是种群大小的选择,这里取popsize=100,即种群由100个个体组成,染色体长度取chromlength=10,也就是二进制编码长度为10。
(2)用随机函数rand()返回一个100行*10列的随机矩阵,矩阵中的每一个数都是0-1的随机数,再用round()函数把这些随机数四舍五入为0或1,这样就形成了100个长度为10的二进制字符串
function pop=initpop(popsize,chromlength)
pop = round(rand(popsize,chromlength));
-
适应度的计算
(1)把二进制转换成十进制
首先是建立表现型到基因型的映射关系,也就是说用二进制编码来表现图中的横坐标x,x也是所要解决问题中的函数自变量,题中要求x是在[5,20]的范围上。
function pop2 = binary2decimal(pop)
lb=5; %x的取值范围下界
ub=20; %x的取值范围上界
[px,py]=size(pop) %px是矩阵的行,py是矩阵的列数,这里px是种群个数,py为染色体长度
%[r,c]=size(A),当有两个输出参数时,size函数将矩阵的行数返回到第一个输出变量r,将矩阵的列数返回到第二个输出变量c。
for i = 1:py
pop1(:,i) = 2.^(py-i).*pop(:,i); %冒号代表对所有行都有这样的操作
%py-i相当于说第一列到最后一列的顺序倒置过来
end
%sum(.,2)对行求和,得到列向量
temp = sum(pop1,2);
pop2 = 5+(temp*(ub-lb)/1023)
%根据X的取值范围进行2进制到10进制的转换,(ub-lb)相当于一个线段,temp/1023相当于二进制在这个区间上的占比,得到的结果再乘以(ub-lb)得到其在这段线段上的位置,再加上下界,就转换成了这个区间上的十进制。
(2)带入函数得到适应度
function [objvalue] = cal_objvalue(pop)
x = binary2decimal(pop);%转化二进制数为x变量的变化域范围的数值
objvalue= x+10*sin(5*x)+7*cos(4*x);
-
选择运算(轮盘赌算法)
(1)先计算出群体中所有个体的适应度总和
(2)其次计算出每个个体的相对适应度的大小,把适应度总和想象成一个轮盘,相对适应度就是每个个体适应度在适应度总和中所占的百分比,即每个个体被遗传到下一代群体中的概率
(3)产生一个0-1之间的随机数,依据该随机数出现在上述哪一个概率区域内来确定各个个体被选中的次数
%如何选择新的个体
%输入变量:pop二进制种群,fitvalue:适应度值
%输出变量:newpop选择以后的二进制种群
function [newpop] = selection(pop,fitvalue)
%构造轮盘
[px,py] = size(pop);
totalfit = sum(fitvalue);%总的适应度值
p_fitvalue = fitvalue/totalfit;%每个适应度值在轮盘中的占比,也就是相对适应度的大小,即为每个个体被遗传到下一代群体中的概率
p_fitvalue = cumsum(p_fitvalue);%概率求和排序 cumsum()得到输入矩阵的每个元素对应的列向上(行向左)求和矩阵
ms = sort(rand(px,1));%从小到大排列 生成种群大小的0-1的随机数字,依据该随机数出现在上述哪一个概率区域内来确定各个个体被选
%中的次数
fitin = 1;
newin = 1;%行选择标量
%选出新的种群
while newin<=px %初始行向量为1,,从第一行开始遍历
if(ms(newin))<p_fitvalue(fitin)%当产生的随机数字<个体被遗传到下一代群体中的概率时,即为可遗传
newpop(newin,:)=pop(fitin,:);%遗传给新的种群
newin = newin+1;%遍历下一行
else
fitin=fitin+1;%如果不可遗传,就检查下一个更大的区域
end
end
-
交叉运算
循环遍历,范围控制在1到数第二行,步长为2。当随机数小于交叉概率时,生成一个随机的交叉点位置。先截取第i行交叉点前面的和第i+1行交叉点后面的编码进行连接;再截取第i+1行交叉点前面编码和第i行交叉点后面的编码连接。如果随机数大于交叉概率,不交叉,直接把pop中的种群数赋给交叉后的种群数
function [newpop] = crossover(pop,pc)
[px,py] = size(pop);
newpop = ones(size(pop));
for i = 1:2:px-1 %for循环的范围在1到px-1,步长为2,也就是从第一行到第三行到第五行一直到第px-1行
if(rand<pc)%当随机数小于交叉概率时才进行交叉
%rand没有参数的时候就是从均匀分布中生成一个0-1之间的随机数
cpoint = round(rand*py);%round是一个四舍五入的函数 用随机数*列数是在随机设置交叉点的位置
newpop(i,:) = [pop(i,1:cpoint),pop(i+1,cpoint+1:py)];%截取第i行交叉点前面的和第i+1行交叉点后面的进行连接
newpop(i+1,:) = [pop(i+1,1:cpoint),pop(i,cpoint+1:py)];%截取第i+1行交叉点前面的和第i行交叉点后面的进行连接
else %如果随机数大于交叉概率,也就是不进行交叉直接把pop中的种群数赋给交叉后的种群数
newpop(i,:) = pop(i,:);
newpop(i+1,:) = pop(i+1,:);
end
end
-
变异运算
从第一行遍历到最后一行,如果产生的随机数小于变异的概率就发生变异,随即产生变异点的位置,对i行变异点的编码取反。否则直接把种群赋给新的种群
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) == 1
newpop(i,mpoint) = 0;
%对第i行变异点取反
end
else newpop(i,:) = pop(i,:);
end
end
-
寻找最优解
首先初始化最佳个体和最佳适应度值bestindividual和bestfit默认为种族中的第一个个体,从第二行开始如果该行的适应度值大于最佳适应度值,就把该行作为最佳个体,该行适应值作为最佳适应值。
function [bestindividual, bestfit] = best(pop,fitvalue)
[px,py] = size(pop); %获取矩阵的行数和列数,size()将矩阵的行数给px,列数给py
bestindividual = pop(1,:);%把第一行的所有列弹出 赋给bestindividual,pop() 函数用于移除列表中的一个元素(默认最后一个元素),并且返回该元素的值。
bestfit = fitvalue(1); %这两行相当于在初始化最佳个体和最佳适应度值
for i = 2:px %for循环从第二行开始,
if fitvalue(i)>bestfit %如果改行的适应值大于最佳适应值,
bestindividual = pop(i,:);%把该行作为最佳个体,
bestfit = fitvalue(i);%该行适应值为最佳适应值
end
end
- mian.m
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,50) == 0
figure; %Matlab中的figure命令,能够创建一个用来显示图形输出的一个窗口对象
fplot('x+10*sin(5*x)+7*cos(4*x)',[5 20]);%自适应绘图,第一个参数为函数,第二个参数为坐标轴取值范围
hold on;% hold on作用是保持原图并接受此后绘制的新的曲线,叠加绘图; 多次叠绘: plot命令可以同时绘制多条曲线
plot(x1,y1,'*');%最优解用*标注
title(['迭代次数为n=' num2str(i)]);%num2str函数将数值转换为字符串
%plot(x1,y1,'*');
end
end
fprintf('The best X is --->>%5.2f\n',x2); %指定数据输出时格式为小数形式
fprintf('The best Y is --->>%5.2f\n',bestfit);