从零开始：如何用NSGA-II算法解决你的第一个多目标优化问题-平芜编程栈

从零开始：如何用NSGA-II算法解决你的第一个多目标优化问题

1. 多目标优化与NSGA-II算法基础

在工程设计和科学研究中，我们经常面临需要同时优化多个相互冲突目标的场景。比如汽车设计中需要平衡燃油经济性和动力性能，芯片设计需要权衡功耗和计算能力，这些都是典型的多目标优化问题(Multi-Objective Optimization Problems, MOPs)。

Pareto最优解是多目标优化中的核心概念。它指的是一组解，在这些解中任何一个目标的改进都会导致至少一个其他目标的恶化。这些解构成了所谓的Pareto前沿(Pareto Front)，为决策者提供了多种权衡选择。

NSGA-II(Non-dominated Sorting Genetic Algorithm II)是当前最流行的多目标优化算法之一，相比第一代NSGA算法，它有三项关键改进：

快速非支配排序：将时间复杂度从O(MN³)降低到O(MN²)，M为目标数，N为种群规模
精英保留策略：防止优秀个体在进化过程中丢失，加速收敛
拥挤度比较算子：无需指定共享参数即可保持种群多样性

% NSGA-II基本参数设置示例 pop_size = 200; % 种群规模 max_gen = 500; % 最大迭代次数 num_obj = 2; % 目标函数数量 num_var = 30; % 决策变量维度

2. NSGA-II算法实现详解

2.1 算法整体流程

NSGA-II的核心流程可以分为以下几个步骤：

初始化种群：随机生成初始解集
非支配排序：将种群分成不同Pareto等级
拥挤度计算：评估同一Pareto等级中个体的分布密度
选择操作：基于锦标赛选择机制挑选父代
遗传操作：通过交叉和变异产生子代
精英保留：合并父代和子代，选择新一代种群

function nsga_2_optimization % 初始化参数 pop = 200; gen = 500; M = 2; V = 30; min_range = zeros(1, V); max_range = ones(1,V); % 初始化种群并排序 chromosome = initialize_variables(pop, M, V, min_range, max_range); chromosome = non_domination_sort_mod(chromosome, M, V); % 主循环 for i = 1 : gen % 选择父代 parent_chromosome = tournament_selection(chromosome, round(pop/2), 2); % 生成子代 offspring_chromosome = genetic_operator(parent_chromosome, M, V, 20, 20, min_range, max_range); % 合并种群并选择新一代 intermediate_chromosome = [chromosome; offspring_chromosome]; chromosome = replace_chromosome(intermediate_chromosome, M, V, pop); end % 结果可视化 if M == 2 plot(chromosome(:,V+1), chromosome(:,V+2),'*'); xlabel('f1'); ylabel('f2'); title('Pareto Optimal Front'); end end

2.2 关键组件实现

快速非支配排序

非支配排序是NSGA-II的核心操作，它将种群划分为多个Pareto等级。第一等级包含不被任何其他个体支配的解，第二等级包含仅被第一等级支配的解，以此类推。

function f = non_domination_sort_mod(x, M, V) [N, ~] = size(x); front = 1; F(front).f = []; individual = []; % 计算支配关系 for i = 1 : N individual(i).n = 0; % 被支配计数 individual(i).p = []; % 支配的个体集合 for j = 1 : N % 比较目标函数值 dom_less = sum(x(i,V+1:V+M) < x(j,V+1:V+M)); dom_equal = sum(x(i,V+1:V+M) == x(j,V+1:V+M)); if dom_less == 0 && dom_equal ~= M individual(i).n = individual(i).n + 1; % i被j支配 elseif sum(x(i,V+1:V+M) > x(j,V+1:V+M)) == 0 && dom_equal ~= M individual(i).p = [individual(i).p j]; % i支配j end end % 第一前沿 if individual(i).n == 0 x(i,M+V+1) = 1; F(front).f = [F(front).f i]; end end % 分层排序 while ~isempty(F(front).f) Q = []; for i = 1 : length(F(front).f) if ~isempty(individual(F(front).f(i)).p) for j = 1 : length(individual(F(front).f(i)).p) individual(individual(F(front).f(i)).p(j)).n = ... individual(individual(F(front).f(i)).p(j)).n - 1; if individual(individual(F(front).f(i)).p(j)).n == 0 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 % 计算拥挤度 [~,index] = sort(x(:,M+V+1)); sorted_based_on_front = x(index,:); current_index = 0; for front = 1 : (length(F)-1) distance = 0; previous_index = current_index + 1; y = sorted_based_on_front(current_index+1:current_index+length(F(front).f),:); current_index = current_index + length(F(front).f); % 对每个目标计算拥挤度 for i = 1 : M [sorted_obj, obj_index] = sort(y(:,V+i)); f_max = sorted_obj(end); f_min = sorted_obj(1); y(obj_index(1),M+V+1+i) = Inf; y(obj_index(end),M+V+1+i) = Inf; for j = 2 : length(obj_index)-1 next_obj = sorted_obj(j+1); prev_obj = sorted_obj(j-1); if (f_max - f_min == 0) y(obj_index(j),M+V+1+i) = Inf; else y(obj_index(j),M+V+1+i) = (next_obj - prev_obj)/(f_max - f_min); end end end % 综合各目标拥挤度 distance = sum(y(:,M+V+2:M+V+1+M),2); y(:,M+V+2) = distance; sorted_based_on_front(previous_index:current_index,:) = y; end f = sorted_based_on_front; end

遗传操作实现

NSGA-II采用模拟二进制交叉(SBX)和多项式变异来生成新个体，这两种操作能有效保持种群的多样性。

function f = genetic_operator(parent_chromosome, M, V, mu, mum, l_limit, u_limit) [N,~] = size(parent_chromosome); child = []; for i = 1 : N/2 % 选择两个不同的父代 parent1 = parent_chromosome(randi(N),1:V); parent2 = parent_chromosome(randi(N),1:V); while isequal(parent1, parent2) parent2 = parent_chromosome(randi(N),1:V); end % SBX交叉 child1 = zeros(1,V); child2 = zeros(1,V); for j = 1 : V u = rand; if u <= 0.5 beta = (2*u)^(1/(mu+1)); else beta = (1/(2*(1-u)))^(1/(mu+1)); end child1(j) = 0.5*((1+beta)*parent1(j) + (1-beta)*parent2(j)); child2(j) = 0.5*((1-beta)*parent1(j) + (1+beta)*parent2(j)); % 边界处理 child1(j) = min(max(child1(j), l_limit(j)), u_limit(j)); child2(j) = min(max(child2(j), l_limit(j)), u_limit(j)); end % 多项式变异 for j = 1 : V if rand < 1/V % 变异概率 delta = (2*rand)^(1/(mum+1)) - 1; child1(j) = child1(j) + delta*(u_limit(j)-l_limit(j)); child1(j) = min(max(child1(j), l_limit(j)), u_limit(j)); end if rand < 1/V delta = 1 - (2*(1-rand))^(1/(mum+1)); child2(j) = child2(j) + delta*(u_limit(j)-l_limit(j)); child2(j) = min(max(child2(j), l_limit(j)), u_limit(j)); end end % 评估目标函数 child1(:,V+1:V+M) = evaluate_objective(child1, M, V); child2(:,V+1:V+M) = evaluate_objective(child2, M, V); child = [child; child1; child2]; end f = child; end

3. 实战案例：ZDT1测试问题

ZDT1是多目标优化领域常用的基准测试函数，包含两个相互冲突的目标：

f₁(x) = x₁
f₂(x) = g(x)[1 - √(x₁/g(x))]

其中g(x) = 1 + 9(∑xᵢ)/(n-1)，i=2,...,n

function f = evaluate_objective(x, M, V) % ZDT1测试函数 f = zeros(1,M); f(1) = x(1); g = 1 + 9*sum(x(2:end))/(V-1); f(2) = g * (1 - sqrt(x(1)/g)); end

3.1 参数设置与结果分析

对于ZDT1问题，我们使用以下参数配置：

参数	值	说明
种群大小	200	每代个体数量
最大代数	500	迭代次数
交叉概率	0.9	SBX交叉概率
变异概率	1/V	每个变量的变异概率
分布指数(η_c)	20	交叉分布指数
分布指数(η_m)	20	变异分布指数

运行NSGA-II后，我们可以得到近似Pareto前沿。理想情况下，ZDT1的Pareto前沿是凸的，对应x₁∈[0,1]且x₂=...=xₙ=0。

提示：在实际应用中，可以通过增加种群规模或迭代次数来获得更接近真实Pareto前沿的解集。但也要注意计算成本与解的质量之间的权衡。

4. 进阶技巧与常见问题

4.1 约束处理

NSGA-II可以通过罚函数法处理约束条件。基本思路是将约束违反程度转化为惩罚项加入目标函数：

function f = evaluate_constrained_objective(x, M, V) % 计算原始目标函数 f = evaluate_objective(x, M, V); % 约束条件示例：g(x) <= 0 constraint_violation = max(0, sum(x) - 5); % 假设约束为∑x <= 5 % 惩罚项 penalty = 1e6 * constraint_violation; f = f + penalty; end