【免费源码】基于MATLAB的精英非支配排序多目标遗传算法NSGAⅡ,非工具箱

作品简介


## 1 引言


在文献[1]中,作者提出了一种基于非支配排序的多目标进化算法(MOEA),称为非支配排序遗传算法 II (NSGA-II),它缓解了进化算法以下三个困难:


1. 时间复杂度为$O(MN^3)$,其中$M$为求解目标数,$N$为种群数目

2. 非精英主义方法

3. 需要指定共享参数


具体来说,提出了一种计算复杂度为$O(MN^2)$的快速非支配排序方法。此外,还提出了一个选择算子,它通过结合父母和后代种群并选择最佳(关于适应度和传播)解决方案来创建交配池。对困难测试问题的仿真结果表明,与 Pareto 归档进化策略和强度 Pareto EA 相比,在大多数问题中,所提出的 NSGA-II 能够在真正的 Pareto 最优前沿附近找到更好的解分布和更好的收敛性——另外两个精英 MOEAs 特别关注创造一个多样化的帕累托最优前沿。


原则上,一个问题中存在多个目标会产生一组最优解(主要是称为帕累托最优解),而不是单个最优解。 在没有任何进一步信息的情况下,不能说这些帕累托最优解中的一个比另一个更好。 这要求用户找到尽可能多的帕累托最优解。 经典优化方法(包括多准则决策方法)建议通过一次强调一个特定的帕累托最优解将多目标优化问题转换为单目标优化问题。 当这种方法用于寻找多个解决方案时,必须多次应用,希望在每次模拟运行时找到不同的解决方案。


在过去几十年中,已经提出了许多多目标进化算法 (MOEA)。由于进化算法 (EA) 与大量解决方案一起工作,因此可以扩展一个简单的 EA 以维护一组不同的解决方案。 强调向真正的帕累托最优区域移动,EA 可用于在一次模拟运行中找到多个帕累托最优解。


 N. Srinivas,Kalyanmoy Deb提出的非支配排序遗传算法 (NSGA) 是最早的此类 EA 之一。 多年来,对 NSGA 方法的主要批评如下:


1. 非支配排序计算复杂度高:NSGA非支配排序算法的计算复杂度为$O(MN^3)$(其中$M$为目标数,$N$为种群规模)。 这使得NSGA 对于大种群规模的计算成本很高。 之所以会出现这种大的复杂性,是因为每一代的非支配排序过程都涉及到复杂性。

2. 缺乏精英主义:最近的研究结果表明,精英主义可以显着加快 GA 的性能,这也有助于防止一旦找到好的解决方案就丢失。

3. 需要指定共享参数$\delta_{share}$:确保种群多样性以得到多种等效解决方案的传统机制主要依赖于共享的概念。 共享的主要问题是它需要指定共享参数($\delta_{share}$)。


本文的NSGA-Ⅱ解决了所有这些问题。 从对许多困难测试问题的模拟结果来看,我们发现 NSGA-II 在寻找多样化集合方面优于其他两个当代 MOEA:帕累托存档进化策略 (PAES)和强度帕累托 EA (SPEA)的解决方案,并在真正的帕累托最优集附近收敛。


## 2 算法介绍


### 2.1 快速非支配排序方法


为了清楚起见,作者首先描述了一个简单而缓慢的过程,将搜索代理分类为不同的非支配级别。 此后,作者描述了一种快速方法。


在一种朴素的方法中,为了识别大小为 $N$ 的种群中第一个非支配前沿的解,可以将每个解与种群中的所有其他解进行比较,以确定它是否被支配。这需要对每个解决方案进行 $O(MN)$ 比较,其中 $M$ 是目标的数量。当这个过程继续寻找种群中第一非支配层的所有成员时,总复杂度为$O(MN^2)$。在这个阶段,第一非支配前沿的所有个体都被找到。为了找到下一个非支配前沿的个体,将第一个前沿的解暂时忽略,重复上述过程。在最坏的情况下,找到第二前沿的任务也需要 $O(MN^2)$ 计算,特别是当 $O(N)$ 数量的解决方案属于第二和更高的非支配水平时。这个论点适用于寻找第三和更高层次的非支配。因此,最坏的情况是当有 $N$ 个前沿并且每个前沿仅存在一个解决方案时。这需要总体 $O(MN^3)$ 计算。请注意,此过程需要 $O(N)$ 存储。在下面的段落和方程式中,作者描述了一种需要 $O(MN^2)$ 计算的快速非支配排序方法。


首先,对于每个解决方案,需要计算两个实体:


- 支配计数$n_p$,支配解决方案 $p$ 的数量

- $S_p$,解决方案 $p$ 支配的一组解决方案。 


这需要 $O(mn^2)$ 比较。


第一个非支配前沿中的所有解决方案都将其支配计数设为零。现在,对于 $n_p=0$ 的每个解 $p$,我们访问其集合 $S_p$ 的每个成员($q$),并将其支配计数减一。 这样做时,如果任何成员 $q$ 的支配计数变为零,我们将其放在单独的列表 $Q$ 中。这些成员属于第二非支配前沿。现在,对 $Q$ 的每个成员继续上述过程,并确定第三条前沿。 这个过程一直持续到所有前沿都被确定为止。


对于第二或更高非支配级别的每个解决方案,支配计数 $n_p$ 最多可以是 $N-1$。因此,每个解决方案 $p$ 在其支配计数变为零之前将被访问最多 $N-1$ 次。此时,该解决方案被分配一个非支配级别,并且永远不会再被访问。由于最多有 $N-1$ 个这样的解,总复杂度为 $O(N^2)$。因此,整个过程的复杂度为 $O(MN^2)$。计算这个复杂度的另一种方法是实现第一个内循环的主体(其中$p\in F_i$)精确执行 $N$ 次,因为每个个体最多可以是一个边沿的成员,并且第二个内部循环(对于$S_p$中的每个 $q$)可以对每个个体执行最多 $(N-1)$ 次[每个个体最多支配 $(N-1)$ 个个体,并且每个支配检查最多需要 $M$ 次比较] 导致总体 $O(MN^2)$ 计算。需要注意的是,虽然时间复杂度降低到了$O(MN^2)$,但是存储需求却增加到了$O(N^2)$。


快速非支配排序伪代码如下图所示:

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-fpkNKk6P-1672494971903)(images/NSGA-Ⅱ-2.png)]


### 2.2 多样性保护


我们之前提到,随着收敛到帕累托最优集,还希望 EA 在获得的解决方案集中保持良好的解决方案传播。 最初的 NSGA 使用众所周知的共享函数方法,该方法已被发现通过适当设置其相关参数来维持种群的可持续多样性。 共享函数方法涉及共享参数 $\delta_{share}$,它设置问题中所需的共享程度。 此参数与选择用于计算两个总体成员之间的邻近度度量的距离度量有关。 参数 $\delta_{share}$ 表示任何两个解决方案共享彼此适应度的距离度量的最大值。 此参数通常由用户设置,尽管存在一些指导原则。 


这种共享函数方法以下两个困难:


- 共享函数方法在维持解的传播方面的性能在很大程度上取决于所选的 $\delta_{share}$ 值。

- 由于每个解决方案都必须与总体中的所有其他解决方案进行比较,因此共享函数方法的总体复杂度为$O(n^2)$。


在提出的 NSGA-II 中,我们将共享函数方法替换为拥挤比较方法,在一定程度上消除了上述两个困难。 新方法不需要任何用户定义的参数来维持种群成员之间的多样性。 此外,建议的方法具有更好的计算复杂度。 为了描述这种方法,我们首先定义一个密度估计度量,然后提出拥挤比较算子。


---

- 1. 密度估计:为了估计总体中特定解决方案周围的解决方案密度,我们计算沿每个目标的该点两侧的两个点的平均距离。 这个量$i_{distance}$用作对使用最近邻居作为顶点形成的长方体周长的估计(称为拥挤距离)。 在图 1 中,第 $i$ 个解在其前面的拥挤距离(用实心圆圈标记)是长方体的平均边长(用虚线框表示)。

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-KQi5lnvJ-1672494971924)(images/NSGA-Ⅱ-1.png)]

拥挤距离计算需要根据每个目标函数值按数量级升序对总体进行排序。 此后,对于每个目标函数,边界解(具有最小和最大函数值的解)被分配一个无限距离值。 所有其他中间解都分配了一个距离值,该距离值等于两个相邻解的函数值的绝对归一化差。 使用其他目标函数继续该计算。 总体拥挤距离值计算为每个目标对应的各个距离值的总和。 在计算拥挤距离之前,对每个目标函数进行归一化。 页面底部显示的算法概述了非支配集 $I$ 中所有解的拥挤距离计算过程。

这里 $I[i]\cdot m$ 是指集合 $I$ 中第 $i$ 个个体的第 $m$ 个目标函数值,参数 $f_m^{max}$ 和 $f_m^{min}$ 是第 $m$ 个目标函数的最大值和最小值。 该过程的复杂性由排序算法控制。 由于涉及到最多 $n$ 个解的 $m$ 个独立排序(当所有种群成员都在一个前面 $I$ 时),所以上述算法具有 $O(mn\log{n})$ 计算复杂度。

在为集合 $I$ 中的所有人口成员分配了距离度量之后,我们可以比较两个解决方案与其他解决方案的接近程度。 在某种意义上,这个距离度量值较小的解决方案被其他解决方案更拥挤。 这正是我们在提议的拥挤比较算子中比较的,如下所述。 尽管图 1 说明了两个目标的拥挤距离计算,但该过程也适用于两个以上的目标。

---

- 2. 拥挤比较算子:拥挤比较算子 ($\prec_n$) 将算法各个阶段的选择过程引导到均匀展开的帕累托最优前沿。 假设代理中的每个个体都有两个属性:

 - 非支配等级($i_{rank}$)

 - 拥挤距离($i_{distance}$)

我们现在将偏序$\prec_n$定义为:

$$

i\prec_n j\quad if(i_{rank}

or\quad((i_{rank}=j_{rank})and(i_{distance}>j_{distance}))

$$

也就是说,在具有不同非支配等级的两个解决方案之间,我们更喜欢具有较低(更好)等级的解决方案。 否则,如果两个解决方案属于同一前沿,那么我们更喜欢位于不太拥挤区域的解决方案。 


拥挤距离计算伪代码如下:

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-Fb6VA7fz-1672494971926)(images/NSGA-Ⅱ-3.png)]


有了这三个新的创新——一个快速的非支配排序过程、一个快速的拥挤距离估计过程和一个简单的拥挤比较算子,我们现在可以描述 NSGA-II 算法了。


### 2.3 算法主体


最初,创建一个随机父群体 $P_0$。 代理是根据非支配性排序的。 每个解决方案都被分配一个与其非支配级别相等的适应度(或等级)(1 是最佳级别,2 是次佳级别,依此类推)。 因此,假设适应度最小化。 首先,使用通常的二元锦标赛选择、重组和变异算子来创建大小为 $n$ 的后代种群 $Q_0$。由于通过将当前种群与先前找到的最佳非支配解进行比较来引入精英主义,因此在初始生成之后过程会有所不同。 我们首先描述所提出算法的第 $t$ 代,伪代码如下图所示:

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-jQLBXkg2-1672494971928)(images/NSGA-Ⅱ-4.png)]


一步一步的过程表明 NSGA-II 算法简单明了。首先,形成一个组合种群 $R_t=P_t∪Q_t$。种群 $R_t$ 的大小为 $2N$。然后,根据非支配对种群 $R_t$ 进行排序。由于所有以前和当前的人口成员都包含在 $R_t$ 中,因此确保了精英主义。现在,属于最佳非支配集 $F_1$ 的解在组合种群中是最佳解,并且必须比组合种群中的任何其他解更加强调。如果 $F_1$ 的大小小于 $N$,我们肯定会选择集合 $F_1$ 的所有成员作为新种群 $P_{t+1}$。种群 $P_{t+1}$ 的其余成员按其排名顺序从后续的非支配前沿中选择。因此,接下来选择集合 $F_2$ 中的解,然后选择集合 $F_3$ 中的解,依此类推。继续此过程,直到无法容纳更多集合。假设集合 $F_1$ 是最后一个非支配集合,超出该集合不能容纳其他集合。一般来说,从 $F_1$ 到 $F_l$ 的所有集合中的解决方案的数量都会大于人口规模。为了精确选择 $N$ 个人口成员,我们使用拥挤比较算子 $\prec_n$ 对最后一个前 $F_l$ 的解决方案进行降序排序并选择填补所有人口空缺所需的最佳解决方案。 NSGA-II 程序也如图 2 所示。大小为 $N$ 的新种群 $P_{t+1}$ 现在用于选择、交叉和突变以创建大小为 $N$ 的新种群 $Q_{t+1}$。重要的是要注意,我们使用二元锦标赛选择算子,但选择标准现在基于拥挤比较算子 $\prec_n$。由于该算子需要总体中每个解决方案的秩和拥挤距离,我们在形成总体 $P_{t+1}$ 时计算这些数量,如上述算法所示。


[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-nvbjlsdP-1672494971943)(images/NSGA-Ⅱ-5.png)]


考虑整个算法的一次迭代的复杂性。 基本操作及其最坏情况的复杂性如下:


- 非支配排序复杂度为:$O(M(2N)^2))$

- 拥挤距离分配复杂度:$O(M(2N)\log(2N))$

- $\prec_n$排序复杂度为:$O(2N\log(2N))$


该算法的总体复杂度为$O(MN^2)$,由算法的非支配排序部分控制。 如果仔细执行,大小为 $2N$ 的完整种群不需要根据非支配排序。 一旦排序过程找到足够数量的前沿以在 $P_{t+1}$ 中有 $N$ 个成员,就没有理由继续排序过程。


通过使用拥挤比较程序引入非支配解决方案之间的多样性,该程序用于锦标赛选择和人口减少阶段。 由于解决方案与其拥挤距离(邻域中解决方案的密度度量)竞争,因此不需要额外的利基参数(例如 NSGA 中需要的 $\delta_{share}$)。 虽然拥挤距离是在目标函数空间中计算的,但如果需要的话,它也可以在参数空间中实现。 然而,在本研究中进行的所有模拟中,我们都使用了目标函数空间细分。


## 3 算法实现


```bash

%% 主函数代码

clear

pop = 200; %种群数量

gen = 500; %迭代次数

M = 2; %目标函数数量

V = 30; %维度(决策变量的个数)

min_range = zeros(1, V); %下界 生成1*30的个体向量 全为0

max_range = ones(1,V); %上界 生成1*30的个体向量 全为1

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

chromosome = initialize_variables(pop, M, V, min_range, max_range);%初始化种群

chromosome = non_domination_sort_mod(chromosome, M, V);%对初始化种群进行非支配快速排序和拥挤度计算



for i = 1 : gen

  pool = round(pop/2);%round() 四舍五入取整 交配池大小

  tour = 2;%竞标赛 参赛选手个数

  parent_chromosome = tournament_selection(chromosome, pool, tour);%竞标赛选择适合繁殖的父代

  mu = 20;%交叉和变异算法的分布指数

  mum = 20;

  offspring_chromosome = genetic_operator(parent_chromosome,M, V, mu, mum, min_range, max_range);%进行交叉变异产生子代 该代码中使用模拟二进制交叉和多项式变异 采用实数编码

  [main_pop,~] = size(chromosome);%父代种群的大小

  [offspring_pop,~] = size(offspring_chromosome);%子代种群的大小

   

  clear temp

  intermediate_chromosome(1:main_pop,:) = chromosome;

  intermediate_chromosome(main_pop + 1 : main_pop + offspring_pop,1 : M+V) = offspring_chromosome;%合并父代种群和子代种群

  intermediate_chromosome = non_domination_sort_mod(intermediate_chromosome, M, V);%对新的种群进行快速非支配排序

  chromosome = replace_chromosome(intermediate_chromosome, M, V, pop);%选择合并种群中前N个优先的个体组成新种群

  if ~mod(i,100)

    clc;

    fprintf('%d generations completed\n',i);

  end

end



if M == 2

  plot(chromosome(:,V + 1),chromosome(:,V + 2),'*');

  xlabel('f_1'); ylabel('f_2');

  title('Pareto Optimal Front');

elseif M == 3

  plot3(chromosome(:,V + 1),chromosome(:,V + 2),chromosome(:,V + 3),'*');

  xlabel('f_1'); ylabel('f_2'); zlabel('f_3');

  title('Pareto Optimal Surface');

end


%% 1 初始化代码


function f = initialize_variables(N, M, V, min_range, max_range)%f是一个由种群个体组成的矩阵

min = min_range;

max = max_range;

K = M + V;%%K是数组的总元素个数。为了便于计算,决策变量和目标函数串在一起形成一个数组。  

%对于交叉和变异,利用目标变量对决策变量进行选择

for i = 1 : N

  for j = 1 : V

    f(i,j) = min(j) + (max(j) - min(j))*rand(1);%f(i j)表示的是种群中第i个个体中的第j个决策变量,

                          %这行代码为每个个体的所有决策变量在约束条件内随机取值

  end

  f(i,V + 1: K) = evaluate_objective(f(i,:), M, V); % M是目标函数数量 V是决策变量个数

                          %为了简化计算将对应的目标函数值储存在染色体的V + 1 到 K的位置。

end

end

%% 2 快速非支配排序和拥挤度计算代码


%% 对初始种群开始排序 快速非支配排序

% 使用非支配排序对种群进行排序。该函数返回每个个体对应的排序值和拥挤距离,是一个两列的矩阵。  

% 并将排序值和拥挤距离添加到染色体矩阵中 

function f = non_domination_sort_mod(x, M, V)

[N, ~] = size(x);% N为矩阵x的行数,也是种群的数量

clear m

front = 1;

F(front).f = [];

individual = [];

 

for i = 1 : N

  individual(i).n = 0;%n是个体i被支配的个体数量

  individual(i).p = [];%p是被个体i支配的个体集合

  for j = 1 : N

    dom_less = 0;

    dom_equal = 0;

    dom_more = 0;

    for k = 1 : M    %判断个体i和个体j的支配关系

      if (x(i,V + k) < x(j,V + k))  

        dom_less = dom_less + 1;

      elseif (x(i,V + k) == x(j,V + k))

        dom_equal = dom_equal + 1;

      else

        dom_more = dom_more + 1;

      end

    end

    if dom_less == 0 && dom_equal ~= M % 说明i受j支配,相应的n加1

      individual(i).n = individual(i).n + 1;

    elseif dom_more == 0 && dom_equal ~= M % 说明i支配j,把j加入i的支配合集中

      individual(i).p = [individual(i).p j];

    end

  end  

  if individual(i).n == 0 %个体i非支配等级排序最高,属于当前最优解集,相应的染色体中携带代表排序数的信息

    x(i,M + V + 1) = 1;

    F(front).f = [F(front).f i];%等级为1的非支配解集

  end

end

%上面的代码是为了找出等级最高的非支配解集

%下面的代码是为了给其他个体进行分级

while ~isempty(F(front).f)

  Q = []; %存放下一个front集合

  for i = 1 : length(F(front).f)%循环当前支配解集中的个体

    if ~isempty(individual(F(front).f(i)).p)%个体i有自己所支配的解集

     for j = 1 : length(individual(F(front).f(i)).p)%循环个体i所支配解集中的个体

       individual(individual(F(front).f(i)).p(j)).n = ...%...表示的是与下一行代码是相连的, 这里表示个体j的被支配个数减1

         individual(individual(F(front).f(i)).p(j)).n - 1;

        if individual(individual(F(front).f(i)).p(j)).n == 0% 如果q是非支配解集,则放入集合Q中

         x(individual(F(front).f(i)).p(j),M + V + 1) = ...%个体染色体中加入分级信息

            front + 1;

          Q = [Q individual(F(front).f(i)).p(j)];

        end

      end

    end

  end

  front = front + 1;

  F(front).f = Q;

end

 

[temp,index_of_fronts] = sort(x(:,M + V + 1));%对个体的代表排序等级的列向量进行升序排序 index_of_fronts表示排序后的值对应原来的索引

for i = 1 : length(index_of_fronts)

  sorted_based_on_front(i,:) = x(index_of_fronts(i),:);%sorted_based_on_front中存放的是x矩阵按照排序等级升序排序后的矩阵

end

current_index = 0;

 

%% Crowding distance 计算每个个体的拥挤度

 

for front = 1 : (length(F) - 1)%这里减1是因为代码55行这里,F的最后一个元素为空,这样才能跳出循环。所以一共有length-1个排序等级

  distance = 0;

  y = [];

  previous_index = current_index + 1;

  for i = 1 : length(F(front).f)

    y(i,:) = sorted_based_on_front(current_index + i,:);%y中存放的是排序等级为front的集合矩阵

  end

  current_index = current_index + i;%current_index =i

  sorted_based_on_objective = [];%存放基于拥挤距离排序的矩阵

  for i = 1 : M

    [sorted_based_on_objective, index_of_objectives] = ...

      sort(y(:,V + i));%按照目标函数值排序

    sorted_based_on_objective = [];

    for j = 1 : length(index_of_objectives)

      sorted_based_on_objective(j,:) = y(index_of_objectives(j),:);% sorted_based_on_objective存放按照目标函数值排序后的x矩阵

    end

    f_max = ...

      sorted_based_on_objective(length(index_of_objectives), V + i);%fmax为目标函数最大值 fmin为目标函数最小值

    f_min = sorted_based_on_objective(1, V + i);

    y(index_of_objectives(length(index_of_objectives)),M + V + 1 + i)...%对排序后的第一个个体和最后一个个体的距离设为无穷大

      = Inf;

    y(index_of_objectives(1),M + V + 1 + i) = Inf;

     for j = 2 : length(index_of_objectives) - 1%循环集合中除了第一个和最后一个的个体

      next_obj = sorted_based_on_objective(j + 1,V + i);

      previous_obj = sorted_based_on_objective(j - 1,V + i);

      if (f_max - f_min == 0)

        y(index_of_objectives(j),M + V + 1 + i) = Inf;

      else

        y(index_of_objectives(j),M + V + 1 + i) = ...

           (next_obj - previous_obj)/(f_max - f_min);

      end

     end

  end

  distance = [];

  distance(:,1) = zeros(length(F(front).f),1);

  for i = 1 : M

    distance(:,1) = distance(:,1) + y(:,M + V + 1 + i);

  end

  y(:,M + V + 2) = distance;

  y = y(:,1 : M + V + 2);

  z(previous_index:current_index,:) = y;

end

f = z();%得到的是已经包含等级和拥挤度的种群矩阵 并且已经按等级排序排序

end

%% 3 竞标赛选择代码


function f = tournament_selection(chromosome, pool_size, tour_size)

[pop, variables] = size(chromosome);%获得种群的个体数量和决策变量数量

rank = variables - 1;%个体向量中排序值所在位置

distance = variables;%个体向量中拥挤度所在位置

%竞标赛选择法,每次随机选择两个个体,优先选择排序等级高的个体,如果排序等级一样,优选选择拥挤度大的个体

for i = 1 : pool_size

  for j = 1 : tour_size

    candidate(j) = round(pop*rand(1));%随机选择参赛个体

    if candidate(j) == 0

      candidate(j) = 1;

    end

    if j > 1

      while ~isempty(find(candidate(1 : j - 1) == candidate(j)))%防止两个参赛个体是同一个

        candidate(j) = round(pop*rand(1));

        if candidate(j) == 0

          candidate(j) = 1;

        end

      end

    end

  end

  for j = 1 : tour_size% 记录每个参赛者的排序等级 拥挤度

    c_obj_rank(j) = chromosome(candidate(j),rank);

    c_obj_distance(j) = chromosome(candidate(j),distance);

  end

  min_candidate = ...

    find(c_obj_rank == min(c_obj_rank));%选择排序等级较小的参赛者,find返回该参赛者的索引

  if length(min_candidate) ~= 1%如果两个参赛者的排序等级相等 则继续比较拥挤度 优先选择拥挤度大的个体

    max_candidate = ...

    find(c_obj_distance(min_candidate) == max(c_obj_distance(min_candidate)));

    if length(max_candidate) ~= 1

      max_candidate = max_candidate(1);

    end

    f(i,:) = chromosome(candidate(min_candidate(max_candidate)),:);

  else

    f(i,:) = chromosome(candidate(min_candidate(1)),:);

  end

end

end

%% 4、5 交叉 变异代码


function f = genetic_operator(parent_chromosome, M, V, mu, mum, l_limit, u_limit)

[N,m] = size(parent_chromosome);%N是交配池中的个体数量

 

clear m

p = 1;

was_crossover = 0;%是否交叉标志位

was_mutation = 0;%是否变异标志位

 

for i = 1 : N%这里虽然循环N次,但是每次循环都会有概率产生2个或者1个子代,所以最终产生的子代个体数量大约是2N个

  if rand(1) < 0.9%交叉概率0.9

    child_1 = [];

    child_2 = [];

    parent_1 = round(N*rand(1));

    if parent_1 < 1

      parent_1 = 1;

    end

    parent_2 = round(N*rand(1));

    if parent_2 < 1

      parent_2 = 1;

    end

    while isequal(parent_chromosome(parent_1,:),parent_chromosome(parent_2,:))

      parent_2 = round(N*rand(1));

      if parent_2 < 1

        parent_2 = 1;

      end

    end

    parent_1 = parent_chromosome(parent_1,:);

    parent_2 = parent_chromosome(parent_2,:);

    for j = 1 : V

      u(j) = rand(1);

      if u(j) <= 0.5

        bq(j) = (2*u(j))^(1/(mu+1));

      else

        bq(j) = (1/(2*(1 - u(j))))^(1/(mu+1));

      end

      child_1(j) = ...

        0.5*(((1 + bq(j))*parent_1(j)) + (1 - bq(j))*parent_2(j));

      child_2(j) = ...

        0.5*(((1 - bq(j))*parent_1(j)) + (1 + bq(j))*parent_2(j));

      if child_1(j) > u_limit(j)

        child_1(j) = u_limit(j);

      elseif child_1(j) < l_limit(j)

        child_1(j) = l_limit(j);

      end

      if child_2(j) > u_limit(j)

        child_2(j) = u_limit(j);

      elseif child_2(j) < l_limit(j)

        child_2(j) = l_limit(j);

      end

    end

    child_1(:,V + 1: M + V) = evaluate_objective(child_1, M, V);

    child_2(:,V + 1: M + V) = evaluate_objective(child_2, M, V);

    was_crossover = 1;

    was_mutation = 0;

  else%if >0.9

    parent_3 = round(N*rand(1));

    if parent_3 < 1

      parent_3 = 1;

    end

    child_3 = parent_chromosome(parent_3,:);

    for j = 1 : V

      r(j) = rand(1);

      if r(j) < 0.5

        delta(j) = (2*r(j))^(1/(mum+1)) - 1;

      else

        delta(j) = 1 - (2*(1 - r(j)))^(1/(mum+1));

      end

      child_3(j) = child_3(j) + delta(j);

      if child_3(j) > u_limit(j) % 条件约束

        child_3(j) = u_limit(j);

      elseif child_3(j) < l_limit(j)

        child_3(j) = l_limit(j);

      end

    end 

    child_3(:,V + 1: M + V) = evaluate_objective(child_3, M, V);

    was_mutation = 1;

    was_crossover = 0;

  end% if <0.9

  if was_crossover

    child(p,:) = child_1;

    child(p+1,:) = child_2;

    was_cossover = 0;

    p = p + 2;

  elseif was_mutation

    child(p,:) = child_3(1,1 : M + V);

    was_mutation = 0;

    p = p + 1;

  end

end

 

f = child;

end

%交叉算法选择的是模拟二进制交叉,变异算法选择的是多项式变异, 算法的具体过程大家可以在网上查阅一下


%% 8 生成新的种群(精英策略)


function f = replace_chromosome(intermediate_chromosome, M, V,pop)%精英选择策略

 

[N, m] = size(intermediate_chromosome);

[temp,index] = sort(intermediate_chromosome(:,M + V + 1));

 

clear temp m

for i = 1 : N

  sorted_chromosome(i,:) = intermediate_chromosome(index(i),:);

end

 

max_rank = max(intermediate_chromosome(:,M + V + 1));

 

previous_index = 0;

for i = 1 : max_rank

  current_index = max(find(sorted_chromosome(:,M + V + 1) == i));

  if current_index > pop

    remaining = pop - previous_index;

    temp_pop = ...

      sorted_chromosome(previous_index + 1 : current_index, :);

    [temp_sort,temp_sort_index] = ...

      sort(temp_pop(:, M + V + 2),'descend');

    for j = 1 : remaining

      f(previous_index + j,:) = temp_pop(temp_sort_index(j),:);

    end

    return;

  elseif current_index < pop

    f(previous_index + 1 : current_index, :) = ...

      sorted_chromosome(previous_index + 1 : current_index, :);

  else

    f(previous_index + 1 : current_index, :) = ...

      sorted_chromosome(previous_index + 1 : current_index, :);

    return;

  end

  previous_index = current_index;

end

end

 

```



## 4 参考文献


[1] Deb K , Pratap A , Agarwal S , et al. A fast and elitist multiobjective genetic algorithm: NSGA-II[J]. IEEE Transactions on Evolutionary Computation, 2002, 6(2):182-197.

## 5 代码免费获取


创作时间:2022-12-31 22:02:27