NSGA_2总结梳理附代码按行详细注解_nsga2

做实验需要解决多目标优化问题，之前也没用过Matlab，看代码也是学习Matlab语法的过程，所以很详细的注解了基本上每一行代码，下面代码亲测可以直接运行，如果有问题的地方欢迎指正。

下面代码可能有些长，主要是注释加的比较多,如果想要替换函数的话，直接在evaluate_objective里替换，在主函数里修改M和V即可

目录

一.NSGA-2算法简介

二.NSGA-2算法整体流程图

三.算法及各函数讲解

1.主函数：nsga_2_optimization

2.目标函数：evaluate_objective

3.初始化代码：initialize_variables

4.快速非支配排序和拥挤度计算代码：non_domination_sort_mod

5.锦标赛选择过程：tournament_selection

6.交叉 变异代码：genetic_operator

7.生成新的种群（精英策略):replace_chromosome

---

一.NSGA-2算法简介

SGA2主要是对NSGA算法的改进。NSGA是N. Srinivas 和 K. Deb在1995年发表的一篇名为《Multiobjective function optimization using nondominated sorting genetic algorithms》的论文中提出的该算法在快速找到Pareto前沿和保持种群多样性方面都有很好的效果，不过在这么多年的应用中也出现了如下的一些问题：

1。非支配排序的时间复杂的很大，为O（MN3）。其中M为目标函数的数量，N为种群规模。

2。不支持精英策略。精英策略在保持好的个体及加速向Pareto前沿收敛方面都有很好的表现。

3。需要自己指定共享参数。该参数将对种群的多样性产生很大的影响。

关于MOP的理解可参考另一篇博客：https://blog.csdn.net/qq_39552268/article/details/111656270

二.NSGA-2算法整体流程图

流程图各个步骤旁边的是函数名对应下面的各个算法函数名

这是最终跑出的结果:

三.算法及各函数讲解

1.主函数：nsga_2_optimization

function nsga_2_optimization
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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;
    %% parent_chromosome 竞标赛选择的适合繁殖的父代 M 目标函数数量 V维度（决策变量的个数）      mu = 20;%交叉和变异算法的分布指数    mum = 20; 
    %%min_range = zeros(1, V); %下界 生成1*30的个体向量 全为0 
    %%max_range = ones(1,V); %上界 生成1*30的个体向量 全为1  这个在这里用来约束解（Xi）的范围
    %%offspring_chromosome不一定生成了500个后代
    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);%对新的种群进行快速非支配排序
    %%精英选择 从子代和父代中选出Pop个
    chromosome = replace_chromosome(intermediate_chromosome, M, V, pop);%选择合并种群中前N个优先的个体组成新种群
   %%每计算100代清空下控制台
    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

2.目标函数：evaluate_objective

这里采用的是两目标函数:ZDT1  ZDT1是MOP中常用的测试函数

%%目标函数 ZDT1是MOP中常用的测试函数
%%这里有两个目标函数 f1 f2 这里就是求f1 f2值得公式
function f = evaluate_objective(x, M, V)%%计算每个个体的M个目标函数值
%%决策变量就是解的个数
%%f(i,V + 1: K) = evaluate_objective(f(i,:), M, V); % M是目标函数数量 V是决策变量个数 f(:,1)就是取f 矩阵的第1列。
f = [];
f(1) = x(1);
g = 1;
sum = 0;
for i = 2:V
    sum = sum + x(i);
end
sum = 9*(sum / (V-1));
g = g + sum;
f(2) = g * (1 - sqrt(x(1) / g));
end

3.初始化代码：initialize_variables

使用在指定范围内的随机值初始化群体。每条染色体由决策变量组成。此外，目标函数，等级和拥挤距离信息的值也被添加到染色体向量中，但是仅对具有决策变量的向量的元素进行操作以执行诸如交叉和变异的遗传操作。

function f = initialize_variables（N，M，V，min_tange，max_range）
N - 目标空间总体大小
M - 目标函数的数量
V - 决策变量的数量

约束条件：
min_range - 十进制值的向量，指示每个决策变量的最小值。
max_range - 决策变量的最大可能值的向量。

初始化决策变量基于可能的最大和最小值，V 是决策变量的个数，从最大最小值之间随机选出一个值作为每一个决策变量
对于简化计算处理染色体的数据和目标方程有着串联的关系，V+1到K 具有目标方程的值，一次目标评价函数带一个染色体，事实上，只有决策变量被传递给函数关于目标函数的数量处理并返回。
目标函数的值

function f = initialize_variables(N, M, V, min_range, max_range)%f是一个由种群个体组成的矩阵
%M目标函数数量
%V维度（决策变量的个数）
%N 种群数量
min = min_range; %下界 生成1*30的个体向量 全为0
max = max_range; %上界 生成1*30的个体向量 全为1
 
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个决策变量， %%-ppppppppppp
                                                    %这行代码为每个个体的所有决策变量在约束条件内随机取值
    end
    f(i,V + 1: K) = evaluate_objective(f(i,:), M, V); % M是目标函数数量 V是决策变量个数 f(:,1)就是取f 矩阵的第1列。
                                                    %为了简化计算将对应的目标函数值储存在染色体的V + 1 到 K的位置。
                                              %%k+1-V 存的是目标函数的值
end

4.快速非支配排序和拥挤度计算代码：non_domination_sort_mod

unction f = non_domination_sort_mod（x，M，V）
此函数根据非支配对当前popultion进行排序。
第一个前面的所有个体的等级为1，
第二个前面的个体被赋予等级2，
依此类推。在分配等级之后，计算每个前沿中的拥挤度。
N - 目标空间总体大小
M - 目标函数的数量 2
V - 决策变量的数量

%% 对初始种群开始排序 快速非支配排序
% 使用非支配排序对种群进行排序。该函数返回每个个体对应的排序值和拥挤距离，是一个两列的矩阵。  
% 并将排序值和拥挤距离添加到染色体矩阵中 
%M + V + 1 M + V + 2 存的是分层等级
%x 进来时是200行（种群个数）*32列（30列x1-x30 2列 f1 f2） 
%x出去时33列存的是当前个体的层级 1为最高
function f = non_domination_sort_mod(x, M, V)
%%chromosome = non_domination_sort_mod(chromosome, M, V); x/chromosome是解
%%其中包括了30个解和 目标函数的值 共32个 M 目标函数数量 V决策变量个数
[N, ~] = size(x);% N为矩阵x的行数，也是种群的数量
clear m
front = 1;%front记录了当前正在筛选那一层级的个体
F(front).f = []; %f是一个数组 用于存放当前front层级的个体
individual = [];%%存放个体i的支配个体的信息
 
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;%1代表最高等级 改个体i的层级为1
        F(front).f = [F(front).f i];%等级为1的非支配解集 f是个矩阵 在F(front).f 矩阵后面加上 i 赋值给F(front).f 
    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(F(front).f(i)).p(j) 代表front层个体所支配的一个个体
            	individual(individual(F(front).f(i)).p(j)).n = ...%...表示的是与下一行代码是相连的， 这里表示个体j的被支配个数减1  
                	individual(individual(F(front).f(i)).p(j)).n - 1; %因为层级最高为1(层级为1 即这时n为0)的在上面已经筛选完 这个循环里面的n最少为1 都被支配
        	   	   %代表去掉front层的个体后,front层个体所支配的一个个体的被支配的个数为0时，这时它是这时候的优解之一
                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
 %x(:,M + V + 1)就是取x矩阵的第M + V + 1列 temp在这里没用 为了就是得到index_of_fronts
 %这里想让这个矩阵按照第M+V+1也就层级进行排序，因为matlab没有直接排序的方法，所以采取这样的方法 x本身的顺序没有变
 %sorted_based_on_front存储了排序后的矩阵
[temp,index_of_fronts] = sort(x(:,M + V + 1));%对个体的代表排序等级的列向量进行升序排序 index_of_fronts表示排序后的值对应原来的未排序的矩阵在当前列的行下标
%x(:,M + V + 1) 就是取这个矩阵的M + V + 1列
for i = 1 : length(index_of_fronts)
    sorted_based_on_front(i,:) = x(index_of_fronts(i),:);%sorted_based_on_front中存放的是x矩阵按照排序等级升序排序后的矩阵 index_of_fronts(i)放到第i行
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; %% current_index记录的是每一层级的第一个个体在sorted_based_on_front的下标-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));%按照目标函数值排序 先按目标f1的值排序 再按f2排序 
        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存放按照目标函数值排序后的y矩阵 在这里的y已经对层级进行过排序
        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) %%如果这一层级里只有一个个体 则这一层的个体fi的距离为无穷
                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); %%这里是归一化处理 用最简单的标准化方法 让间距的值小一点 y的M + V + 1 + 1记录f1的拥挤度  M + V + 1 + 2记录
            end
         end
    end
    distance = [];
    distance(:,1) = zeros(length(F(front).f),1);%%初始化 这个列向量用于保存当前front层的每一个体的拥挤度
    for i = 1 : M
        distance(:,1) = distance(:,1) + y(:,M + V + 1 + i);  %%将y的f1的距离和f2的距离相加作为该个体的拥挤度
    end
    y(:,M + V + 2) = distance;
    y = y(:,1 : M + V + 2); %%应该是可省的
    z(previous_index:current_index,:) = y; %%在这里 previous_index记录该层级中第一个个体在矩阵中的行下标 current_index记录该层级中最后一个个体在矩阵中的行下标
end
f = z();%得到的是已经包含等级和拥挤度的种群矩阵 并且已经按等级和拥挤距离排序

5.锦标赛选择过程：tournament_selection

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

该算法是为交配池选择个体的选择策略。选择基于锦标赛选择。参数染色体是在进行比赛选择后，从当前的世代种群中选择个体形成一个size pool_size的交配池，比赛的size变为tour_size。通过改变比赛规模，选择压力可以调整。但是对于NSGA-II, tour_size被固定为2，但是用户可以自由地尝试不同的锦标赛大小。此外，人们还注意到，一场锦标赛规模超过5是没有任何意义的。
锦标赛选择过程
在一个锦标赛的选择过程中，n个人被随机选择，其中n等于tour_size。从这些个体中只选择一个，并添加到交配池中，其中交配池的大小为pool_size。选择是基于两个标准执行的。首先也是最重要的是解决方案所处的位置。选择级别较低的个人。其次，如果两个个体的秩相同，则比较拥挤距离。选择拥挤距离较大的个体。

function f = tournament_selection(chromosome, pool_size, tour_size)
%chromosome记录了整个种群
%%tour_size=2;%竞标赛  参赛选手个数
%%pool_size round(pop/2); 种群的一半 交配池大小
[pop, variables] = size(chromosome);%获得种群的个体数量和决策变量数量
rank = variables - 1;%个体向量中排序值所在位置 层级位置
distance = variables;%个体向量中拥挤度所在位置
%竞标赛选择法，每次随机选择两个个体，优先选择排序等级高的个体，如果排序等级一样，优选选择拥挤度大的个体
for i = 1 : pool_size
    %%这里随机从交配池选择了两个个体 存到candidate
    for j = 1 : tour_size
        candidate(j) = round(pop*rand(1));%随机选择参赛个体 
        if candidate(j) == 0 %随机到0时处理下 因为参赛个体下标从1开始
            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返回该参赛者的索引 
    %%c_obj_rank是一个一维数组，这里返回的就是这个数组里值等于min(c_obj_rank)的数的下标，可能返回一个数组
    if length(min_candidate) ~= 1%如果两个参赛者的排序等级相等 则继续比较拥挤度 优先选择拥挤度大的个体 
        %！=1时意味着两个个体的层级相等 这时min_candidate里这两个个体有相同层级时 拥挤距离大的  觉得代码写的比较冗余
        max_candidate = ...
        find(c_obj_distance(min_candidate) == max(c_obj_distance(min_candidate)));%选择出距离大的一个
    % find(c_obj_distance(min_candidate)
    % ==max(c_obj_distance(min_candidate))) 这里两边都把(min_candidate)去了应该没影响
    % 没有测试
        if length(max_candidate) ~= 1 %若还等于1 则说明两个层级和拥挤度都一样 互相不能支配
            max_candidate = max_candidate(1);
        end
        %%总之 min_candidate记录的是 比较层级时 候选个体中层级小的
        %%  max_candidate记录的是 层级相等的时候 候选个体中拥挤度
        f(i,:) = chromosome(candidate(min_candidate(max_candidate)),:);
    else
        f(i,:) = chromosome(candidate(min_candidate(1)),:);
    end
    %%算法到 这里共筛选出了 总群大小/2的非劣个体
end

6.交叉 变异代码：genetic_operator

交叉算法选择的是模拟二进制交叉，变异算法选择的是多项式变异

模拟二进制交叉（SBX）：

多项式变异:

function f  = genetic_operator(parent_chromosome, M, V, mu, mum, l_limit, u_limit)
 
%% 交叉变异代码
%% 交叉算法选择的是模拟二进制交叉，变异算法选择的是多项式变异
%% parent_chromosome 竞标赛选择的适合繁殖的父代 M 目标函数数量 V维度（决策变量的个数）     
%% mu = 20;%交叉和变异算法的分布指数    mum = 20; 
%% min_range = zeros(1, V); %下界 生成1*30的个体向量 全为0 
%% max_range = ones(1,V); %上界 生成1*30的个体向量 全为1 
[N,m] = size(parent_chromosome);%N是交配池中的个体数量 原种群数量的一半 这里是250
 %%从原种群数量的一半的较优父代中进行交叉变异生成后代个体,下一轮选择之前,子代个体最大是原种群大小
clear m
p = 1;
was_crossover = 0;%是否交叉标志位
was_mutation = 0;%是否变异标志位
 
%%这里设置90%的概率交叉 10%概率变异 
for i = 1 : N%这里虽然循环N次，但是每次循环都会有概率产生2个或者1个子代，所以最终产生的子代个体数量大约是2N个
    if rand(1) < 0.9%交叉概率0.9 也就是有几率不交叉，只变异 可调
        child_1 = [];
        child_2 = [];
        %%初始化子代种群
        
        %%利用SBX进行交叉
        %随机选取两个父代个体
        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_1 = parent_chromosome(parent_1,:);
        parent_2 = parent_chromosome(parent_2,:);
        %%计算贝特 mu就是公式中的恩特 是人为设定的
        %%这里循环是让每个个体的x都交叉
        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));
            %%u_limit就是设置的x经过归一化的最大取值 l_limit最小取值 让交叉过后的后代x在规定的范围内
            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
        %% 将x带入公式求f1 f2
        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
 
    %%child保存了所有变异或交叉过后的子代个体 p作为 child的索引
    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;

7.生成新的种群（精英策略):replace_chromosome

1.首先将父代种群C i C_iCi​和子代种群D i D_iDi​合成种群R i R_iRi​。
2.根据以下规则从种群R i R_iRi​生成新的父代种群C i + 1 C_{i+1}Ci+1​：
      ①根据Pareto等级从低到高的顺序，将整层种群放入父代种群C i + 1 C_{i+1}Ci+1​，直到某一层该层个体不能全部放入父代种群C i + 1 C_{i+1}Ci+1​；
      ②将该层个体根据拥挤度从大到小排列，依次放入父代种群C i + 1 C_{i+1}Ci+1​中，直到父代种群C i + 1 C_{i+1}Ci+1​填满。

function f  = replace_chromosome(intermediate_chromosome, M, V,pop)%精英选择策略
%%replace_chromosome(intermediate_chromosome, M, V, pop)
%%intermediate_chromosom 子代和父代两代种群放在一起进行非支配排序后的矩阵
 %%生成新的种群（精英策略）
[N, m] = size(intermediate_chromosome);
%%获得按层级进行排序的种群索引
[temp,index] = sort(intermediate_chromosome(:,M + V + 1)); %%Matlab中没有办法直接让矩阵按照某一列直接排序所以只养血
 
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));%%这样写避免了遍历寻找当前层级为i的个体 直接得到当前层级个体最大的索引，因为已经排过序
    %%所以根据current_index就可以得到当前层级的所有个体
    %%current_inde的位置就代表了当前已经筛选出了多少个个体了 因为前面排过序
    %当该层级的个体大于所需的种群大小时
    if current_index > pop
        %%总共需要筛选出pop个个体,remaining保存还需筛选出多少个体新一代种群才达到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');
        %%获得该层级中拥挤距离大的remaining个 个体
        for j = 1 : remaining
            f(previous_index + j,:) = temp_pop(temp_sort_index(j),:);
        end
        return;
    elseif current_index < pop
        %%小于Pop 说明新种群中个体数量<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

参考:

      1>NSGA2 算法MATLAB完整代码 中文注释详解_joe的博客-CSDN博客_nsga2算法matlab

      2>NSGA_2 Matlab 算法详解完整代码 中文注释详解_机器学习初学者必看，关注我，一起了解机器学习-CSDN博客_nsga2算法matlab代码