IT小白

这个屌丝很懒，什么也没留下！

热门标签

Matlab无基础快速上手1（遗传算法框架）

作者：IT小白 | 2024-04-28 23:13:19

踩

本文用经典遗传算法框架模板，对matlab新手友好，快速上手看懂matlab代码，快速应用实践，源代码在文末给出。

基本原理：

遗传算法（Genetic Algorithm，GA）是一种受生物学启发的优化算法，它模拟了生物进化中的基因遗传和自然选择过程，通过模拟“进化”过程来搜索问题的最优解或近似最优解。

遗传算法的基本思想是通过对候选解（称为个体或染色体）的不断进化和变异来搜索最优解。在遗传算法中，每个个体都代表了问题的一个潜在解，并且具有一个与其适应度相关联的数值，用于评价个体的优劣。

遗传算法的主要步骤包括：

初始化：随机生成一组初始个体（种群），并为每个个体计算适应度。
选择：根据个体的适应度值，使用选择算子（如轮盘赌算法、竞争选择等）选择一部分个体作为父代。
交叉：从父代中选取一对个体，通过交叉操作产生新的个体（子代），以增加种群的多样性。
变异：对新生成的个体进行变异操作，以引入新的基因组合，增加种群的多样性。
替换：使用替换策略（如保留最优个体、精英保留等）更新种群，产生下一代。
重复迭代：重复进行选择、交叉、变异和替换步骤，直到满足停止条件（如达到最大迭代次数、找到满意解等）。

通过不断重复这些步骤，遗传算法能够在解空间中搜索出较好的解，逐步逼近最优解。由于其并行性和全局搜索能力，遗传算法被广泛应用于各种领域的优化问题，如工程优化、机器学习、数据挖掘等。

参数与子函数定义：

此代码的待解决问题为function函数，也就是多元函数平方和求极小值问题，看懂代码后，可将function应用为自己的待解决问题。


function MySGA
NP = 50;                     %种群规模
Max_N = 1000;                %迭代次数
Pc = ones(1,NP)*0.8;         %ones表示创建一个全为1的数组，参数为数组大小
Pm = ones(1,NP)*0.1;         %变异概率
D = 500;
MinX = -100;
MaxX = 100;
Error = 1.0e-10;
X = MinX + (MaxX - MinX)*rand(NP,D);
F=fun(X);
[BestF,BestLow] = min(F);
 
function F = fun(X)
for i = 1:1:length(X(:,1))
    F(i) = sum(X(i,:).^2)
end

代码讲解：

ones函数：

这里出现了ones函数

用法为：A = ones(sz)
sz为参数，可以是一个表示数组大小的向量，可以是一个数字标量、一个向量或一个矩阵
这个函数表示将一个空间内的元素都置为1
ones（3）创建一个 3x3 的矩阵，所有元素都为 1。只有一个参数默认还是2维x2维的矩阵。尽量使用ones（3，3）来创建，增加代码可读性
ones（2,4）创建一个 2x4 的矩阵，所有元素都为 1。
ones（[2,3,4]）创建一个 2x3x4 的数组，所有元素都为 1。

标量，向量，矩阵

这里提到了标量，向量，矩阵，这里继续学习一下：

数字标量：就是一个普通的数
向量：向量是一种特殊的数组，它只有一个维度，有水平和垂直两种向量，也可以说行向量或者列向量，例如：行向量：v = [1, 2, 3, 4];列向量：v = [1; 2; 3; 4];
矩阵：矩阵是一个由数字按行和列排列成的二维数组。
创建一个矩阵：A = [1, 2, 3; 4, 5, 6; 7, 8, 9];

这里的Pc = ones(1,NP)*0.8就表示创建一个包含1行，NP列的数组，全为0.8，表示变异概率。

这里出现了Error = 1.0e-10;，这种表示方法为科学计数法表示 $1 \times 10^{-10}$

接下来是X = MinX + (MaxX - MinX)*rand(NP,D);这里先看rand函数

rand函数：

定义：rand函数用于生成指定大小的服从均匀分布的随机数矩阵或随机数向量
参数：两个参数表示生成mxn维的矩阵，矩阵元素为0-1之间的随机数，若有一个参数为n则表示生成一个包含n个0-1之间随机数的向量，若无参数则生成一个0-1之间的随机数
返回：例如r = rand(3, 3)，这会生成一个3x3的矩阵返回给r，且矩阵的元素为（0，1）之间的随机值

再看这个代码， (MaxX - MinX)*rand(NP,D)表示生成为（0，MaxX - MinX）直接的随机值的NP行，D列的矩阵，前面再加上MinX就表示生成范围为（MinX，MaxX）之间的NP行，D列的矩阵。

F=fun(X)

子函数代码：


function F = fun(X)
for i = 1:1:length(X(:,1))
    F(i) = sum(X(i,:).^2)
end

for函数：


for index = start:increment:end 
% 在这里执行循环体代码 
end

start是开始值，increment是步长，end是结束值，步长为1的话可以省略increment

length函数：

length函数用于返回数组的长度或矩阵的最大维度

例如：


vec = [1, 2, 3, 4, 5];
len = length(vec); % 返回 5，因为 vec 含有 5 个元素

len为5，因为有5个元素


mat = [1, 2, 3; 4, 5, 6];
len = length(mat); % 返回 3，因为矩阵的列数为 3

len为3，因为矩阵有3列

X(:,1)用法：

这表示选中X矩阵的所有行的第一个元素，：表示所有，也就是提取出X的第一列的内容

分析此函数功能

for i = 1:1:length(X(:,1))表示根据X的第一列的长度进行遍历，也就是循环次数为X的行数。

F(i) = sum(X(i,:).^2)这个代码中，X(i,:)表示选中第i行的所有列，也就是提取出第i行元素，.^2表示进行乘方运算，那么这个代码的运行结果就是将第i行的所有元素平方和相加赋值给F的第i个元素。

min函数：

定义：min函数用于找到数组或矩阵中的最小值及其对应的索引（若有多个列则按列索引）
格式：[M, I] = min(A)
参数：参数A为待索引矩阵或向量
返回值：返回值为一个一维矩阵，一行或一列矩阵返回该行或该列的最小值和所在行/列位置，多行矩阵的话，第一个元素为A第一列的最小值，第二个元素为A第一列最小值的索引位置

主函数实现：

选择：


for gen = 1:1:Max_N %循环迭代次数
    time(gen) = gen;
    tmpF = cumsum(F/sum(F);     %归一化
    disp(tmpF);
    CopyX = X;
    for i = 1:1:NP %通过适应度选择个体
       rnd = rand;%生成0-1之间的数
       for flag = 1:1:NP
          if(rnd  < tmpF(flag))
              break;
          end
          X(i,:) = CopyX(flag,:);
       end    
    end
    
end

如上述代码所示，首先第一个for循环是循环迭代次数，下面的代码为每一次迭代所需要的操作，其中，选择操作采用了轮盘赌选择法。

time数组记录当前迭代次数，tmpF是制作好的轮盘。

分析tmpF的实现：

首先轮盘赌选择法是基于适应度进行选择的，适应度与选中概率成正比，适应度越高那么被选中的概率越大。轮盘的制作分为两步：

1.根据适应度计算概率

设种群规模为 $N$ ， $F_{i}$ 表示第i个个体的适应值，那么这个个体被选中的概率为： $P_{i}=\frac{F_{i}}{\sum_{j=1}^{N}F_{j}}$

代码中F/sum(F)实现此功能。

2.归一化制作轮盘

要实现概率越高选中越大就要将这个所有个体的概率归一化，那么使用cumsum函数计算每个个体与它前所有元素的累加和，这样就制作好了一个轮盘，实现所有个体被选中的累计概率为1。如下图所示：因此第一个for循环代码表示遍历整个种群对每个个体进行更新，rnd=rand表示模拟轮盘转动。第二个for循环进行个体选择。

上述提到的是经典遗传算法问题的轮盘制作过程，本模板需要计算的是多元函数求极值问题，具体为多个未知数的平方和最小问题，因此这里用F表示适应度的话，F越小、适应度越小被选中概率越高，因此轮盘制作代码需要改为如下：

tmpF = cumsum((1./F)/sum(1./F));

交叉：


    Pc_rand=rand(1,NP);
    for i=1:2:(NP-1)
        if Pc(i) > Pc_rand(i)
            alfa=rand;
            CpX = X(i,:);
            X(i,:)=alfa*X(i+1,:)+(1-alfa)*X(i,:);
            X(i+1,:)=alfa*CpX+(1-alfa)*X(i+1,:);
        end
    end

Pc_rand为1行NP列，也就是一行NP个元素，每个元素为0-1随机值的矩阵。使用它判断每两个个体是否进行交叉操作，Pc为1行NP列，元素全为0.8的矩阵，代表进行交叉操作的概率。

开始循环，循环次数为种群的一半，代表一次是否选中两个个体进行交叉操作，循环里进行判断是否进行交叉操作，实现了交叉操作概率为0.8。

交叉操作：生成一个0-1之间的随机数alfa = rand；将第i个个体乘以（1-alfa）加上第i+1个个体乘以alfa更新为第i个个体，第i+1个个体同理。

变异：


    Pm_rand=rand(NP,D);
    for i=1:1:NP
        for j=1:1:D
            if Pm(i) > Pm_rand(i,j)
                X(i,j)=MinX+(MaxX-MinX)*rand;
            end
        end
    end

这个代码就很简洁了，遍历每个个体的每个参数，当0.01比这个值大时（0.01就是变异概率），就随机改变这个参数的值，也就是变异了。

后续处理：


    X(NP,:)=bestX;%保留最佳个体防止淘汰
    F=fun(X); %评估当前种群
    [BestF,BestLow] = min(F);
    bestX=X(BestLow,:);
    result(gen)=BestF;  %记录结果

X(NP,:)=bestX就是将上一次循环后的最优个体保留，防止种群的最后一个个体位置，确保每次迭代的最佳个体被保留在种群中，以防止其在后续的进化过程中被淘汰掉。

整体代码：


function MySGA
NP = 100;                     %种群规模
Max_N = 10000;                %迭代次数
Pc = ones(1,NP)*0.8;         %ones表示创建一个全为1的数组，参数为数组大小
Pm = ones(1,NP)*0.01;         %变异概率
flagc=[0,Max_N]; 
D = 50;
MinX = -100;
MaxX = 100;
Error = 1.0e-10;
X = MinX + (MaxX - MinX)*rand(NP,D);
F=fun(X);
[BestF,BestLow] = min(F);
bestX=X(BestLow,:);
%disp(F);
for gen = 1:1:Max_N %循环迭代次数
    time(gen) = gen;
    tmpF = cumsum((1./F)/sum(1./F));     %归一化
    %disp(tmpF);
    CopyX = X;
    for i = 1:1:NP %通过适应度选择个体
       rnd = rand;%生成0-1之间的数
       for flag = 1:1:NP
          if(rnd  < tmpF(flag))
              break;
          end
          X(i,:) = CopyX(flag,:);
       end    
    end
    %种群交叉
    Pc_rand=rand(1,NP);
    for i=1:2:(NP-1)
        if Pc(i) > Pc_rand(i)
            alfa=rand;
            CpX = X(i,:);
            X(i,:)=alfa*X(i+1,:)+(1-alfa)*X(i,:);
            X(i+1,:)=alfa*CpX+(1-alfa)*X(i+1,:);
        end
    end
    %种群变异
    Pm_rand=rand(NP,D);
    for i=1:1:NP
        for j=1:1:D
            if Pm(i) > Pm_rand(i,j)
                X(i,j)=MinX+(MaxX-MinX)*rand;
            end
        end
    end
    X(NP,:)=bestX;%保留最佳个体防止淘汰
    F=fun(X);
    %------------------------------ 种群评估 -------------------------------
    [BestF,BestLow] = min(F);
    bestX=X(BestLow,:);
    %--------------------------- 记录结果 ----------------------------------
    result(gen)=BestF;
    if (BestF<Error) & (flagc(1)==0)
        flagc(1)=1;flagc(2)=gen;
    end
    if mod(gen,1000)==0
        disp(['代数：',num2str(gen),'----最优：',num2str(BestF)]);
        %disp(X);
    end
end
plot(time,result)
function F = fun(X)
for i = 1:1:length(X(:,1))
    F(i) = sum(X(i,:).^2);
end

此代码在交叉操作中还可以采用别的方式，可以读懂代码后，继续尝试增加多种方式。

声明：本文内容由网友自发贡献，不代表【wpsshop博客】立场，版权归原作者所有，本站不承担相应法律责任。如您发现有侵权的内容，请联系我们。转载请注明出处：https://www.wpsshop.cn/w/IT小白/article/detail/504765