NSGA-II多目标优化:非支配排序与拥挤距离的MATLAB实现解析 📅 发布时间:2026/9/17 1:56:45 👁 浏览次数: 简介这份MATLAB资源实现了经典的非支配排序遗传算法NSGA-II用于求解多目标优化问题。面向需要同时权衡多个冲突目标的研究人员、工程师及算法学习者可在工程设计、投资决策、资源分配等场景下直接使用或二次开发。压缩包共13个文件核心为12个.m脚本包含快速非支配排序、拥挤距离计算、锦标赛选择、交叉与变异等完整模块并附有适应度函数示例、种群初始化与帕累托前沿绘图工具另有1个来源链接文件整体仅6KB轻量且结构清晰。目前已有989人学习下载。使用者只需定义自己的目标函数并设置种群规模、进化代数等参数即可运行主程序观察解集迭代收敛过程与前沿分布特别适合作为理解NSGA-II机制、开展对比实验或搭建自研多目标优化框架的起步模板。1. 一个被权重害惨的多目标问题NSGA-II 是怎么绕过去的工程优化里最常见的做法是把几个互相打架的指标用权重压成一个目标函数。权重一旦拍错优化器只会沿着某条签名路径冲到底最后给你的往往是一个平庸折中甚至完全不能用的方案。NSGA-II 的高明之处在于不打权重而是直接在多个目标之间建立“谁支配谁”的关系并用非支配排序把解集一层层分成不同等级的帕累托前沿。MATLAB 这份 NSGA-II 实现把整个流程拆成 NonDominatedSorting.m、CalcCrowdingDistance.m、Crossover.m、Mutate.m 等独立脚本改目标函数不需要动排序逻辑适合要快速验证多目标算法的人。你不需要先啃完一本多目标教材只要弄懂支配、拥挤距离和精英保留这三个关键设计收敛性和多样性就能同时抓在手里。2. 非支配排序与拥挤距离NSGA-II 的排序机制拆解这一章从这份包里最核心的几个 m 文件入手先讲支配判断怎么写再讲层数怎么切分最后讲拥挤距离如何保多样性。2.1 支配判断NSGA-II Dominates.m 是整个算法的地基在多目标优化里解与解之间的关系不是简单的“大”或者“小”。假设我们做最小化目标函数值向量为Cost那么解 a 支配解 b 的条件是a 在所有目标上的表现都不比 b 差并且至少在一个目标上严格优于 b。如果两个解互不支配它们就有可能是帕累托前沿上的候选点。这份包里的NSGA-II Dominates.m实现的就是这一句逻辑核心代码大概长这样function d Dominates(x, y) % 目标默认按最小化处理 % allx 在所有目标上都不比 y 差 % anyx 至少在某个目标上严格优于 y d all(x.Cost y.Cost) any(x.Cost y.Cost); end这里x.Cost和y.Cost必须是同尺寸的向量顺序必须与你在Cost.m里定义的目标函数一致。需要注意如果某个目标要求最大化不要在Dominates.m里改判断方向更稳妥的做法是在Cost.m里对那个目标取负号。这样整个排序模块始终保持“最小化”语义后续所有代码都不需要分裂成两套逻辑。这段代码看起来短但它决定了非支配排序的正确性。一个常见的错误是把all写成any或者把严格小于漏掉结果就是算法稳定地退化成普通遗传算法前沿面上会出现大量明显被支配的点。2.2 快速非支配排序NonDominatedSorting.m 如何切出第一层解集有了支配判断接下来就是把整个种群划分成不同非支配层。朴素的 NSGA 做法是不断循环“找出当前剩余个体中的所有非支配解标记一层后从集合里删掉”复杂度较高。NSGA-II 的快速非支配排序通过维护两个数据结构来加速一个是dominanceCount记录有多少个个体支配当前个体另一个是dominatedSet记录当前个体支配了哪些个体。伪代码简化后如下function [pop, F] NonDominatedSorting(pop) n numel(pop); dominanceCount zeros(n, 1); dominatedSet cell(n, 1); F {}; for i 1:n for j i1:n if Dominates(pop(i), pop(j)) dominatedSet{i}(end1) j; dominanceCount(j) dominanceCount(j) 1; elseif Dominates(pop(j), pop(i)) dominatedSet{j}(end1) i; dominanceCount(i) dominanceCount(i) 1; end end if dominanceCount(i) 0 pop(i).Rank 1; F{1}(end1) i; end end k 1; while ~isempty(F{k}) Q []; for i F{k} for j dominatedSet{i} dominanceCount(j) dominanceCount(j) - 1; if dominanceCount(j) 0 pop(j).Rank k 1; Q(end1) j; end end end k k 1; F{k} Q; end end这段代码中第一层是所有dominanceCount为 0 的个体也就是当前种群中任何解都支配不了的那些候选点。剥离完第一层后维护的计数会被更新再找出新的“准第一层”成为第二层。这样每个个体最终都会得到一个Rank字段Rank 1代表最接近帕累托前沿的个体Rank越大越靠后。这里的复杂度看起来仍然是嵌套循环但因为每层剥离时只访问被支配集合中的个体整体比较次数要远小于每次全量比较的朴素实现。种群规模在 500 以下时两者差别不大一旦超过 1000快速排序的优势会非常明显。2.3 拥挤距离与 SortPopulation.m让解集不抱团非支配排序保证了“质量”但如果有几十个个体都落在Rank 1同一层内下一轮选择时该优先保留谁NSGA-II 的答案是拥挤距离。拥挤距离计算同一非支配层里每个个体周围目标空间中被邻居包围的密度距离越大代表周围越“空旷”也就越值得保留。CalcCrowdingDistance.m的常见实现如下function pop CalcCrowdingDistance(pop, F) nObj numel(pop(1).Cost); for i 1:numel(F) members F{i}; n numel(members); if n 3 continue; end dist zeros(1, n); for j 1:nObj [~, idx] sort([pop(members).Cost(j)], ascend); dist(idx(1)) inf; dist(idx(end)) inf; range pop(members(idx(end))).Cost(j) - pop(members(idx(1))).Cost(j); if range 0 continue; end for k 2:n-1 dist(idx(k)) dist(idx(k)) ... (pop(members(idx(k1))).Cost(j) - pop(members(idx(k-1))).Cost(j)) / range; end end for k 1:n pop(members(k)).CrowdingDistance dist(k); end end end这段逻辑里有两个细节值得注意。第一每个目标维度上排序后最边界的两个个体距离被直接赋成了inf这是刻意保留极端解的技巧否则帕累托前沿的两端会被中间个体挤掉。第二各目标的距离贡献做了归一化避免量纲差异过大的目标主导距离计算。排序动作最终交给SortPopulation.m它会按照“Rank 升序、CrowdingDistance 降序”把种群排成一个线性序列。之后BinaryTournamentSelection.m只需随机抓两个个体挑排在前面的那个进交配池选择压力就自然形成了。下面是朴素 NSGA 与 NSGA-II 排序机制的一个粗略对比维度朴素 NSGANSGA-II 快速非支配排序分层方式每轮全量扫描后剔除维护计数 逐层剥离同层多样性共享函数计算开销大拥挤距离 边界 inf时间开销每轮重复比较多比较总次数明显下降参数敏感性需要调共享半径无需额外距离参数SortPopulation.m之后被nsga2.m主循环复用合并父子种群后截断到原来的种群规模这一步就是 NSGA-II 的精英保留机制。3. 从 CreateEmptyIndividuals.m 到 Cost.m把自定义目标函数挂进 NSGA-II这一章处理最关键的一步如何把自己要优化的目标函数填进这套算法中以及初始化、取值和显示三个文件之间是怎么配合的。3.1 CreateEmptyIndividuals.m 里的个体结构不是数组而是结构体很多第一次接触这份代码的人会问为什么不能用一个大矩阵装着所有设计变量原因在于多目标优化中每个个体除了位置还要记录目标值、非支配层级、拥挤距离有时还需要记录约束违反度。把这些字段打包成一个结构体数组是 MATLAB 里最直观的做法。CreateEmptyIndividuals.m的职责就是预先建好这样一个空结构体数组function pop CreateEmptyIndividuals(n) % 建立一个包含 n 个结构体元素的空种群 empty_individual.Position []; empty_individual.Cost []; empty_individual.Rank []; empty_individual.CrowdingDistance []; % repmat 复制 n 份得到 n 行 1 列的种群 pop repmat(empty_individual, n, 1); end这段代码里Position是设计变量向量Cost是目标函数值向量Rank和CrowdingDistance是排序阶段填写的字段。结构体数组的好处是访问单个个体时可以直接写pop(i).Cost在多目标模块间传递数据时语义非常清晰。代价是性能不如纯数值矩阵但 MATLAB 里结构体数组对于几百个个体、几十维目标的问题完全够用。3.2 Cost.m 与 GetCosts.m多目标输出怎么组织Cost.m是整套代码里唯一需要你动手改的核心文件。它接收一个设计变量向量返回一个目标函数值向量。以一个双目标优化问题为例function z Cost(x) % 示例两个互相冲突的目标 % x 是 Position 向量比如 x [x1, x2] z(1) x(1)^2 (x(2) - 1)^2; z(2) (x(1) - 1)^2 x(2)^2; end这里z(1)、z(2)构造的是一个 1×2 行向量。Dominates.m里的all(x.Cost y.Cost)在 MATLAB 中按元素比较时并不要求行向量与列向量只要两个Cost尺寸一致即可。不过为了统一建议始终让Cost.m返回行向量因为后续GetCosts.m会按“每个个体占一行”的方式重组数据如果混用行列向量提取结果会非常混乱。GetCosts.m正是负责把种群中所有人的Cost字段批量提出来形成一个矩阵供绘图和计算指标使用function costs GetCosts(pop) % 从结构体数组中提取所有个体的目标值 % numel(pop(1).Cost) 是目标函数的个数 costs reshape([pop.Cost], numel(pop(1).Cost), []); end这里[pop.Cost]会把所有个体的目标值横向拼在一起画图时常见的问题是发现costs的维度反了。加上reshape并转置后得到的是nPop × nObj的矩阵第一列是第一目标第二列是第二目标。这个约定直接服务于PlotFronts.m也让你自己写统计时不用反复猜数据排列。3.3 PlotFronts.m 与个体字段的使用约定PlotFronts.m在代码包里更多是结果显示模块。它接收GetCosts.m提取出的矩阵再根据个体的Rank字段给散点上色。典型的调用方式是在nsga2.m跑完后把所有代的第一个前沿面画在一张图上观察解从前代到后代如何往前沿靠拢。下面这张表把个体字段与使用场景对应起来方便你在改代码时快速定位字段类型写入时机消费方Position行向量CreateEmptyIndividuals / CrossoverMutate、CostCost行向量Cost.mDominates、GetCostsRank标量NonDominatedSortingSortPopulation、PlotFrontsCrowdingDistance标量CalcCrowdingDistanceSortPopulation、锦标赛选择需要提醒的是如果你要加一个约束条件不要直接在Cost.m里返回一个很大的惩罚数。挤压到一个目标上会让 Pareto 排序失真更常见的做法是给个体增加一个Violation字段在排序时先比较约束违反度再比较Rank。这份包原本没有这个字段你可以按上面结构体定义的方式手动加不影响核心排序逻辑。4. 交叉、变异与代际循环nsga2.m 里如何实现和调整遗传算子NSGA-II 的核心框架说到底是遗传算法。交叉负责在已知区域附近探索变异负责跳出局部热点而代际循环负责把新解重新纳入排序体系。这一章把三个部分串起来讲清楚。4.1 Crossover.m 与 Mutate.m模拟二进制交叉与多项式变异多目标优化里最常用的连续变量交叉算子不是单点交叉而是模拟二进制交叉SBX。它能保证两个父代产生两个子代的均值与父代一致适合实数编码。一份典型实现如下function [y1, y2] Crossover(x1, x2, etaC) % x1, x2 是父代设计变量向量etaC 是交叉分布指数 r rand(size(x1)); mask r 0.5; beta (2*r).^(1/(etaC1)) .* mask ... (1./(2*(1-r))).^(1/(etaC1)) .* (1-mask); y1 0.5 * ((1beta).*x1 (1-beta).*x2); y2 0.5 * ((1-beta).*x1 (1beta).*x2); end这里beta是决定子代偏离父代程度的系数。etaC越大beta分布越集中子代离父代越近etaC越小子代搜索范围越大。一般工业问题的etaC取 10 到 20 之间太大会导致算法在一代之内几乎无法产生新的探索方向。变异算子则常用多项式变异核心是在当前变量值上叠加一个小扰动function y Mutate(x, pm, etaM, VarMin, VarMax) % pm 是每个基因的变异概率etaM 是变异分布指数 y x; for i 1:numel(y) if rand pm u rand; if u 0.5 delta (2*u)^(1/(etaM1)) - 1; else delta 1 - (2*(1-u))^(1/(etaM1)); end y(i) x(i) delta * (VarMax(i) - VarMin(i)); y(i) max(min(y(i), VarMax(i)), VarMin(i)); end end end变异代码最后一行特别重要边界裁剪必须在扰动计算完成之后做。如果先裁剪再计算扰动方向位于边界附近的变量会失去向边界外试探的机会。pm的取值和变量维度有关通常取1 / numel(x)到 0.1 之间过大容易把已经收敛的好解打散。4.2 nsga2.m 的主循环精英保留与选择压力nsga2.m是整套代码的入口。它的主循环逻辑比一般遗传算法多两步第一步是合并父代和子代第二步是重新排序并截断。伪代码如下for gen 1:MaxGen % 1. 用二进制锦标赛从 pop 中选出父代 parent BinaryTournamentSelection(pop); % 2. 对父代两两执行 Crossover 和 Mutate生成 offspring offspring Crossover(parent(1:2:end), parent(2:2:end), etaC); offspring Mutate(offspring, pm, etaM, VarMin, VarMax); % 3. 父子合并 merged [pop; offspring]; % 4. 重新排序并计算拥挤距离 [merged, F] NonDominatedSorting(merged); merged CalcCrowdingDistance(merged, F); [merged, ~] SortPopulation(merged); % 5. 只保留前 nPop 个个体形成下一代 pop merged(1:nPop); end这个循环里最有价值的点是父子合并后的截断机制。普通遗传算法是父子竞争只有子代可能替换父代NSGA-II 则让父代和子代站到同一条线上重新竞争所以上一代Rank 1的个体不可能被这一代运气差就丢掉。这就是精英保留的显式实现也是 NSGA-II 收敛性比早期 NSGA 稳定的关键原因。BinaryTournamentSelection.m的逻辑是随机挑两个个体比较Rank和CrowdingDistance选择排序更靠前的。这一步不需要额外设置选择概率选择压力完全由排序结果决定代码维护也更简单。4.3 参数与收敛性种群、代数、交叉概率与变异概率怎么调这份包本身没有复杂的配置文件参数全部写在nsga2.m头部。下面是实际调参时常用的参考范围参数常见范围对收敛性的影响nPop50 ~ 200太小时前沿断裂太大时每代排序变慢MaxGen100 ~ 500决定前沿能推进多远pc0.8 ~ 1.0控制子代生成的频率pm1/nVar ~ 0.2太大破坏收敛太小陷入局部etaC10 ~ 20控制交叉步长etaM15 ~ 30控制变异步长一个常用的调参启动点是nPop 100MaxGen 200pc 0.9pm 0.1etaC 15etaM 20。如果发现最终前沿上解分布有明显空洞优先把 nPop 调大而不是加 MaxGen。如果发现前沿面长期不移动检查一下是不是pc太小或者变异概率过高导致好基因被频繁打散。对于目标函数计算量很大的工程问题你可以把评估函数做成批量向量化但注意不要破坏Cost.m返回行向量的约定。否则GetCosts.m和支配比较都会出现隐性错误。更简单的做法是在循环外预计算一部分与种群无关的数据而不是急着改主循环并行逻辑。5. 用 PlotFronts.m 和重复实验验证 NSGA-II 的收敛性最后一章不谈原理只讲两个能直接用起来的验证技巧一个是画图一个是用多次独立运行代替单次观察。5.1 动态刷新前沿边跑边看是否卡住很多人在nsga2.m跑完之后才调PlotFronts.m这时如果前沿很差你根本不知道是算法的问题还是目标函数定义的问题。更有效的做法是在代际循环内部每隔几代刷新一次图形if mod(gen, 10) 0 costs GetCosts(pop); plot(costs(:,1), costs(:,2), o); drawnow; end这里mod(gen, 10)的意思是每十代刷新一次既不影响主循环也能看到前沿从初始散点逐步收敛到窄带的过程。判断收敛性的一个直觉标准是最后几十代点只在一个很小的带状区域内移动而不是还在大范围跳动。如果 200 代后前沿依然像一团雾多半是pm过大或精英保留被截断逻辑破坏了。5.2 用多轮独立运行画箱线图单次运行只能证明这次问题收敛了不能说明算法稳定。工程上我会写一个外层脚本重复运行多次并记录每轮的最优目标值变化然后用箱线图检验波动范围nRuns 10; history zeros(nRuns, 2); for r 1:nRuns pop RunNSGA2(); % 把主循环封装成一个函数 costs GetCosts(pop); history(r, :) min(costs, [], 1); end boxplot(history, {目标1, 目标2});这里的min(costs, [], 1)是取每个目标在最终代中的最小值用来观察算法是否有偶尔失手的情况。如果箱体宽度明显大于前沿本身宽度说明你选用的nPop或MaxGen对这个问题偏小需要增加代数或增大交叉分布指数。如果你想让多轮运行结果更可复现可以在每轮运行前固定随机数种子这样对比参数变化时排除随机因素。本文还有配套的精品资源点击获取