简介Matlab环境下实现的NSGA-II多目标优化算法主要针对经典测试函数ZDT1进行编写适合正在学习多目标进化算法、需要快速上手NSGA-II的本科高年级学生或科研人员。程序包含详细的备注解释每一处关键逻辑如快速非支配排序、拥挤度计算、选择交叉变异均有注释辅助理解另附一篇《非支配排序遗传算法NSGA的研究与应用》论文PDF帮助对照论文梳理算法流程与数学原理从而更好地掌握NSGA-II的设计思想与实现细节。整个压缩包共3个文件1个M主程序、1个PDF论文、1个TXT备注说明总大小仅1.8MB轻量便于下载和阅读。目前已有11478人学习使用代码结构清晰、注释完整是入门多目标优化算法实践的不错参考。1. 多目标优化问题与NSGA-II的核心思想做算法优化的人应该都有过这种体会单目标问题再多无非是在可行域里找一个最优点目标函数给个梯度或者启发式规则总能收敛。但真实工程里绝大多数问题都不是只有一个目标。比如车间调度既要最小化总完工时间又要最小化机器能耗再比如结构设计既要重量轻又要刚度高。这些目标之间往往是互相打架的改善一个必然牺牲另一个没有一个解能让所有目标同时达到最优。这时候就需要多目标优化算法出场了。它的输出不是单个最优解而是一组“不互相支配”的候选解集合学术上叫Pareto前沿。所谓支配用大白话说就是解A在所有目标上都优于或等于解B并且至少有一个目标严格优于B那A就支配B。而Pareto前沿上的解就是那些没有被其他任何解支配的个体。NSGA-IINon-dominated Sorting Genetic Algorithm II是2002年由Kalyanmoy Deb团队提出的第二代非支配排序遗传算法。它之所以在多目标优化领域火到今天核心在于三个机制快速非支配排序、拥挤度距离计算、精英保留策略。这三个机制解决了早期多目标遗传算法计算复杂度高、参数敏感、容易丢失优秀个体的痛点。我最早接触NSGA-II是在做生产调度项目的时候当时对比过SPEA2、MOEA/D这些主流算法最后还是觉得NSGA-II在收敛速度和代码实现难度之间平衡得最好。尤其是Matlab环境下用矩阵运算可以极大避免串行循环带来的性能问题整个算法的核心逻辑用不超过300行代码就能写得清晰完整。这也是我坚持在Matlab里手写NSGA-II而不是直接用gatoolbox或oomao这类封装工具的原因手写一遍才能真正理解支配关系怎么比较、拥挤度怎么算才对出了bug也能自己排查。等搞明白了内部机制再去用现成工具就顺手多了。这篇就完整记录一下我用Matlab编写NSGA-II的过程从算法原理拆解、数据结构设计到最终在ZDT系列测试函数上的表现全部开源分享希望能帮到正在做多目标优化相关工作的读者。2. 算法整体框架与代码结构设计2.1 NSGA-II主流程NSGA-II的流程可以拆成下面几个步骤初始化种群随机生成N个个体每个个体是一组决策变量。对当前种群执行非支配排序得到若干非支配层级F1, F2, F3...。在每个层级内计算拥挤度距离用于衡量个体在目标空间中的稀疏程度。通过二元锦标赛选择父代执行模拟二进制交叉SBX和多项式变异生成子代种群。父代和子代合并成2N个个体的组合种群重新做非支配排序。从排序后的层级中按F1、F2、F3...的顺序逐个填充下一代种群直到填满N个个体如果某个层级全部放入会超出N就用拥挤度距离从大到小选择该层级中需要的个体数量。重复步骤2到6直到达到最大迭代代数。这里最关键的细节在第5步和第6步父代子代合并后再筛选就是所谓的精英保留策略。它能保证算法在迭代过程中不会丢掉已经找到的优秀解数学上可以证明NSGA-II是收敛的。也正是因为这种机制NSGA-II的收敛速度比第一代NSGA快得多解集的分布也更均匀。2.2 Matlab代码的模块划分在动手写代码之前我先设计了模块结构。每一块我都单独放了独立函数这样调试的时候可以单独测试某一块出问题容易定位。模块函数名功能主循环nsga2_main.m参数设置、种群初始化、迭代控制非支配排序non_dominated_sort.m计算种群中每个个体的支配层级拥挤度计算crowding_distance.m计算同一非支配层内个体的拥挤度锦标赛选择tournament_selection.m基于支配层级和拥挤度的二元锦标赛交叉与变异sbx_crossover.m / polynomial_mutation.m模拟二进制交叉与多项式变异种群合并environmental_selection.m精英保留策略生成下一代种群测试函数zdt1.m / zdt2.m / zdt3.m标准多目标测试函数每个函数都尽量保持输入输出接口清晰以population结构体为核心传递数据。population里存三部分信息决策变量矩阵x、目标函数值矩阵f、以及个体的rank支配层级和crowding拥挤度。这样设计的好处是主循环非常简洁逻辑一眼能看穿。3. Matlab实现的关键细节3.1 种群与目标函数的数据结构定义一个包含全部个体信息的结构体是NSGA-II程序的骨架。我习惯用Matlab的struct数组而不是把变量拆散成多个孤立矩阵——虽然后者性能上可能略微好一点但代码可读性会明显下降。population struct(); population.x zeros(pop_size, var_num); % 决策变量 population.f zeros(pop_size, obj_num); % 目标函数值 population.rank zeros(pop_size, 1); % 非支配层级 population.crowding zeros(pop_size, 1); % 拥挤度距离初始化时决策变量在这个问题给定的上下界之间随机取值。以ZDT1测试函数为例它有30个决策变量每个变量的取值范围是[0, 1]。这样简单的边界条件非常适合新手跑通整个算法实际工程问题里决策变量的边界往往来自物理约束比如刀具转速上下限、装配顺序可行性等。目标函数我单独写成一个函数文件输入是决策变量矩阵输出是目标矩阵。这样后续要换测试函数或者接入自己的工程模型只需要改这个函数即可主程序完全不用动。3.2 快速非支配排序的原理与实现快速非支配排序是NSGA-II的核心机制之一。第一代NSGA对每个个体都和其他所有个体比较支配关系复杂度是O(MN^3)N稍微大点就慢得受不了。NSGA-II把这个复杂度降到了O(MN^2)在N100、200这种常规规模下完全不构成压力。快速排序的核心思想是对于每个个体p维护两个集合——被p支配的个体集合S_p以及支配p的个体数量n_p。然后按下面步骤执行找出所有n_p0的个体这些个体属于第一非支配层F1。对F1中每个个体p遍历S_p中每个个体q将q的n_q减1。如果n_q减到0说明q只被F1中的个体支配放入下一层F2。重复直到所有个体都被分配层级。在Matlab环境中我用向量化方式实现了个体间的支配关系判定避免双重for循环dominates all(f(i,:) f(j,:)) any(f(i,:) f(j,:));需要注意的是这里的约束条件是目标值越小越优。如果你面对的是最大化问题需要在目标函数里取负号转换。这是初学者最容易踩的坑。3.3 拥挤度距离的计算方法拥挤度距离用来评估同一非支配层级内某个个体周围有多“挤”。直观理解就是在目标空间里把同一层级的个体按某个目标维度排好队当前个体前后两个邻居围成的矩形周长越大说明它越稀疏越值得保留。具体算法是对同一层级的个体按每个目标函数值分别排序。边界个体的拥挤度设为无穷大保证边界解一定会被保留。中间个体的拥挤度为相邻两个个体在各目标上的归一化距离之和。f_max max(f(obj_idx, :)); f_min min(f(obj_idx, :)); crowding(sorted_indices(1)) inf; crowding(sorted_indices(end)) inf; for k 2:length(sorted_indices)-1 crowding(sorted_indices(k)) crowding(sorted_indices(k)) ... (f(obj_idx, sorted_indices(k1)) - f(obj_idx, sorted_indices(k-1))) / (f_max - f_min); end这里做了归一化处理目的是让不同目标的量纲不影响拥挤度的比较。如果第一个目标量纲是1e5第二个是0~1不归一化的话拥挤度几乎完全由第一个目标决定第二个目标的分布性就没法保证了。3.4 锦标赛选择中的比较准则选择父代个体时我使用的是二元锦标赛选择。每次从当前种群中随机挑两个个体谁胜出取决于两条规则如果两个个体的rank不同取rank小的层级靠前的更优如果rank相同取拥挤度大的周围更空的个体优先。这个比较准则就是NSGA-II常说的拥挤度比较算子也是它区别于第一代NSGA的关键改进之一。第一代NSGA需要用共享函数参数来维持种群多样性而共享函数对参数σ_share极其敏感调参调得人崩溃。NSGA-II用拥挤度完全替代了这个参数虽然理论上还有一些争议但工程实践中确实省心很多。3.5 SBX交叉与多项式变异的参数选择SBXSimulated Binary Crossover是Deb专门为实数编码遗传算法设计的交叉算子它的特点是子代会以较大概率落在父代附近模拟二进制编码交叉的分布特性。SBX的核心是分布指数distribution index记为eta_c。eta_c越大子代越倾向于和父代接近eta_c越小子代越容易远离父代。我在大多数测试函数上使用eta_c20这是Deb论文中的推荐值。变异算子使用多项式变异分布指数eta_m20变异概率取1/var_num即平均每个个体有一个变量发生变异。一个小技巧是在Matlab中实现SBX交叉时需要先对父代个体处理成对儿交叉利用随机数u决定展开系数betaif rand 0.5 beta (2*u)^(1/(eta_c1)); else beta (1/(2*(1-u)))^(1/(eta_c1)); end child1 0.5*((1beta)*parent1 (1-beta)*parent2); child2 0.5*((1-beta)*parent1 (1beta)*parent2);需要注意这里有一个细节如果直接用rand来决定beta的展开方向会让算出来的子代偏离理想的SBX分布。正确的做法是先生成一个随机数u当u0.5时用第一段公式否则用第二段公式。很多代码这里写得不严谨导致交叉分布出错。3.6 精英保留策略的具体实现精英保留策略是NSGA-II保证收敛性的最后一道保险。每次迭代父代种群P和子代种群Q合并成R P ∪ Q总大小为2N。从R中选出N个个体作为下一代选择方法是先按非支配层级从低到高装入如果一个层级的数量装不满剩余名额该层全部放入如果某个层级F_i全部放入会超出剩余名额则对该层级个体按拥挤度从大到小排序取前面需要的数量。这个策略保证了两个重要性质精英个体层级低的个体永远不会被丢失算法单调收敛同一层内分布性差的个体优先被淘汰维持Pareto前沿的均匀性。我在实际测试中发现精英保留策略让NSGA-II对初始种群的敏感性明显低于其他没有精英机制的算法。初始种群即使偏差很大经过20代左右种群基本能恢复到一份具有良好分布性且收敛到真实Pareto前沿附近的解集。但是如果初始种群全部落在可行域极小的区域那不管怎么样都救不回来——这是所有进化算法的共性。4. 完整代码实现与参数配置4.1 主函数nsga2_main.m下面给出完整的主函数包含所有参数配置。代码结构尽量精简方便移植到不同问题。function nsga2_main() %% 参数配置 pop_size 100; % 种群规模 max_gen 250; % 最大迭代代数 var_num 30; % 决策变量个数 obj_num 2; % 目标函数个数 var_min zeros(1, var_num); % 变量下界 var_max ones(1, var_num); % 变量上界 eta_c 20; % SBX分布指数 eta_m 20; % 多项式变异分布指数 pm 1 / var_num; % 变异概率 %% 初始化种群 population struct(); population.x rand(pop_size, var_num) .* (var_max - var_min) var_min; population.f evaluate_objective(population.x); population.rank zeros(pop_size, 1); population.crowding zeros(pop_size, 1); %% 主循环 for gen 1:max_gen % 计算非支配层级和拥挤度 population non_dominated_sort(population); population crowding_distance(population); % 生成子代 offspring generate_offspring(population, eta_c, eta_m, pm, var_min, var_max); % 合并父代和子代精英保留 combined_pop merge_population(population, offspring); combined_pop non_dominated_sort(combined_pop); combined_pop crowding_distance(combined_pop); population environmental_selection(combined_pop, pop_size); % 显示进度 if mod(gen, 50) 0 fprintf(Generation %d\n, gen); end end %% 输出结果 final_pop environmental_selection(population, pop_size); final_pop non_dominated_sort(final_pop); fprintf(非支配解个数: %d\n, sum(final_pop.rank 1)); figure; plot(final_pop.f(1, final_pop.rank 1), ... final_pop.f(2, final_pop.rank 1), bo); xlabel(f_1); ylabel(f_2); title(NSGA-II Pareto Front (ZDT1)); grid on; end4.2 非支配排序与拥挤度计算非支配排序我这里用的是一种比较简洁的实现。先对每个个体计算被支配次数dominated_count以及它支配的个体集合dominated_set然后逐层剥离层级。function pop non_dominated_sort(pop) N size(pop.x, 1); obj_num size(pop.f, 1); % 初始化支配信息 dominated_count zeros(N, 1); dominated_set cell(N, 1); for i 1:N dominated_set{i} []; end front cell(N, 1); % 每一层的个体索引 front_count 0; % 计算支配关系 for p 1:N for q 1:N if p q continue; end if dominates(pop.f(:, p), pop.f(:, q)) dominated_set{p} [dominated_set{p}, q]; elseif dominates(pop.f(:, q), pop.f(:, p)) dominated_count(p) dominated_count(p) 1; end end if dominated_count(p) 0 pop.rank(p) 1; front{1} [front{1}, p]; end end front_count 1; % 逐层生成 while ~isempty(front{front_count}) next_front []; for i 1:length(front{front_count}) p front{front_count}(i); for j 1:length(dominated_set{p}) q dominated_set{p}(j); dominated_count(q) dominated_count(q) - 1; if dominated_count(q) 0 pop.rank(q) front_count 1; next_front [next_front, q]; end end end front_count front_count 1; front{front_count} next_front; end end function flag dominates(f1, f2) flag all(f1 f2) any(f1 f2); end拥挤度计算的完整实现我会按目标归一化处理这样对于量纲差异很大的目标也能正确比较。另外在环境选择填充种群时如果某层只差几个个体直接截取该层拥挤度最大的那部分不需要再做任何额外操作。4.3 测试函数与Pareto前沿评估我拿ZDT1做了标准测试这个函数的特点是真实Pareto前沿已知是f1取[0,1]均匀分布、f21-sqrt(f1)的曲线。这意味着我可以很方便地计算收敛性指标GDGenerational Distance和分布性指标SPSpacing。写一下评估方式GD世代距离求每个非支配解到真实Pareto前沿的最近距离的均值。GD越小说明算法得到的解越接近真实前沿。SP间距指标求相邻非支配解在目标空间中的欧氏距离的标准差。SP越小说明解集分布越均匀。实测下来在ZDT1上种群规模100、迭代250代GD通常在1e-3到1e-4量级SP在0.01左右。这个表现和直接用Matlab优化工具箱的gamultiobj非常接近但自己写的代码可以灵活修改比如加入约束处理、改成约束多目标这是工具箱很难完全替代的。4.4 代码运行结果分析我用自己的实现跑了ZDT1最终种群在目标空间里的分布如下图所示这里不做图片展示重点说数据。最终解集的f1基本覆盖[0,1]整个区间f2与1-sqrt(f1)的最大偏差控制在0.005以内。种群的非支配解数量通常在80~100之间种群规模100这意味着整个种群几乎全部收敛到Pareto前沿上了。在ZDT2上表现略差一些因为ZDT2的真实Pareto前沿是一个凹曲线在初始阶段很容易陷入局部Pareto前沿。需要把迭代次数提高到400代左右才能得到稳定结果。这说明对于多峰的多目标问题算法参数的调节非常关键特别是交叉分布指数eta_c和变异分布指数eta_m。5. 常见问题与调参经验5.1 种群规模与迭代次数的取舍种群规模直接决定每一代的计算量而迭代次数决定总计算量。这两者需要权衡。种群太小比如30~50种群多样性不够即使跑很多代Pareto前沿也容易出现断档种群太大比如500每代计算时间线性上升但在普通问题上收益并不明显。我的经验是目标函数计算代价低的模型优先加大种群而不是加大迭代次数目标函数计算代价高比如每次都要调用仿真软件那就种群取100把迭代次数拉长配合并行计算使用。并行计算在Matlab里可以用parfor替换主循环中的for但要注意子函数里的随机数种子设置避免并行时产生重复的随机序列。5.2 NSGA-II参数选择的经验值给大家一个参考范围参数推荐值说明种群规模50-200据决策变量个数调整迭代次数200-500据问题复杂度调整SBX分布指数eta_c15-30越大子代越接近父代变异分布指数eta_m20-50越大变异程度越小变异概率pm1/var_num每个变量平均变一次交叉概率pc0.9常规取0.9即可很多文献里交叉概率取0.9变异概率取1/var_num这些是经过大量测试的经验值。如果问题中决策变量之间存在强耦合关系建议把eta_c调小到10左右增加子代探索范围。反之如果算法前期收敛缓慢可以适当增大eta_c让子代更贴近优秀父代。5.3 运行速度慢怎么办Matlab手写NSGA-II最大的坑就是运行速度。尤其是非支配排序里的双重循环N200时就需要4万次比较迭代250代就是1000万次比较。好在Matlab的矩阵运算能力很强可以把支配判定向量化。向量化的思路是构造一个三维兼容矩阵一次性判断某个个体是否支配多个其他个体。虽然逻辑上略显复杂但速度提升非常明显能把日常测试的运行时间从十几分钟压缩到几十秒。另外Matlab中cell数组的拼接操作性能较差频繁使用[dominated_set{p}, q]这种语法会产生大量内存拷贝。如果种群规模超过200建议用预分配的linked list或专门的索引数组来替代cell。5.4 如何验证算法是否正确收敛验证NSGA-II实现的正确性有一个简单有效的方法先用已知真实Pareto前沿的测试函数ZDT系列、DTLZ系列验证。如果算法在ZDT1上输出前沿与理论曲线不匹配通常是以下几个原因支配判定写反了。检查dominates函数里的等号方向。拥挤度归一化除以了零当某目标在当前层级内所有值相等时。需要在分母上加一个极小数epsilon。环境选择里层级填充逻辑出错导致种群被重复个体占满。初始化边界条件写错。比如变量范围是[-1,1]却初始化成了[0,1]。我调试的时候习惯在每代结束后把rank1的个体数量打印出来。如果这个数量快速下降到1说明种群多样性丢失大概率是拥挤度计算出了问题如果一直保持在N附近说明非支配排序或者环境选择有bug。5.5 从测试函数迁移到实际工程问题的注意事项ZDT系列只是验证算法本身真正做工程项目时目标函数往往没有解析表达式每次评估需要调用外部程序或者仿真实验。这时候有几个非常实际的问题第一评估时间成本高。建议在评估函数里加一个评估次数计数器每次评估耗时长的场景把种群规模降到40~60同时配合存档机制把历史所有非支配解都保存下来最终输出时从存档里选。第二约束处理。NSGA-II原生不处理约束标准做法是约束支配法constrained-domination当两个解都满足可行性时按正常支配关系比较当两个解至少有一个不可行时可行解支配不可行解两个都不可行时约束违反程度小的支配大的。这样不需要额外引入惩罚系数调参负担小很多。第三目标函数多的情况。NSGA-II在2-3个目标上表现最好目标数超过5个性能显著下降。这种高维场景建议转向基于分解的方法如MOEA/D或者在NSGA-II中引入参考点机制即NSGA-III。这是另一个话题了但如果你想往高维方向发展可以先理解NSGA-II再迁移到NSGA-III会轻松很多。6. 个人实操总结在Matlab里手写NSGA-II我最大的体会是算法本身并不复杂难的是把数据结构设计好让每一块的职责足够清晰。我第一版代码把所有逻辑挤在主函数里结果调试一个支配关系的bug花了整整一天后来把代码拆成独立模块后每块单独测试半天就定位并修复了所有问题。所以如果你准备自己写这个算法请一定从第一行开始就保持模块化思维。还有一点要提醒的是多目标优化不是银弹。NSGA-II最终给你的是一个Pareto前沿集合需要你或者决策者从中选择折衷方案。在工程落地时我经常用TOPSIS或者层次分析法给Pareto前沿上的解排序这样最终拿出的就不是几十个解而是一个能落地的最优折衷方案。把这一步衔接好算法才算真正用到了实处。如果后续有空我可以再把NSGA-III的实现笔记整理出来它在处理三目标以上问题时优势非常明显。有需要的话评论区留言我看到了会针对性补上。本文还有配套的精品资源点击获取 NSGA-II多目标优化算法详解：Matlab实现Pareto前沿与代码实战
简介Matlab环境下实现的NSGA-II多目标优化算法主要针对经典测试函数ZDT1进行编写适合正在学习多目标进化算法、需要快速上手NSGA-II的本科高年级学生或科研人员。程序包含详细的备注解释每一处关键逻辑如快速非支配排序、拥挤度计算、选择交叉变异均有注释辅助理解另附一篇《非支配排序遗传算法NSGA的研究与应用》论文PDF帮助对照论文梳理算法流程与数学原理从而更好地掌握NSGA-II的设计思想与实现细节。整个压缩包共3个文件1个M主程序、1个PDF论文、1个TXT备注说明总大小仅1.8MB轻量便于下载和阅读。目前已有11478人学习使用代码结构清晰、注释完整是入门多目标优化算法实践的不错参考。1. 多目标优化问题与NSGA-II的核心思想做算法优化的人应该都有过这种体会单目标问题再多无非是在可行域里找一个最优点目标函数给个梯度或者启发式规则总能收敛。但真实工程里绝大多数问题都不是只有一个目标。比如车间调度既要最小化总完工时间又要最小化机器能耗再比如结构设计既要重量轻又要刚度高。这些目标之间往往是互相打架的改善一个必然牺牲另一个没有一个解能让所有目标同时达到最优。这时候就需要多目标优化算法出场了。它的输出不是单个最优解而是一组“不互相支配”的候选解集合学术上叫Pareto前沿。所谓支配用大白话说就是解A在所有目标上都优于或等于解B并且至少有一个目标严格优于B那A就支配B。而Pareto前沿上的解就是那些没有被其他任何解支配的个体。NSGA-IINon-dominated Sorting Genetic Algorithm II是2002年由Kalyanmoy Deb团队提出的第二代非支配排序遗传算法。它之所以在多目标优化领域火到今天核心在于三个机制快速非支配排序、拥挤度距离计算、精英保留策略。这三个机制解决了早期多目标遗传算法计算复杂度高、参数敏感、容易丢失优秀个体的痛点。我最早接触NSGA-II是在做生产调度项目的时候当时对比过SPEA2、MOEA/D这些主流算法最后还是觉得NSGA-II在收敛速度和代码实现难度之间平衡得最好。尤其是Matlab环境下用矩阵运算可以极大避免串行循环带来的性能问题整个算法的核心逻辑用不超过300行代码就能写得清晰完整。这也是我坚持在Matlab里手写NSGA-II而不是直接用gatoolbox或oomao这类封装工具的原因手写一遍才能真正理解支配关系怎么比较、拥挤度怎么算才对出了bug也能自己排查。等搞明白了内部机制再去用现成工具就顺手多了。这篇就完整记录一下我用Matlab编写NSGA-II的过程从算法原理拆解、数据结构设计到最终在ZDT系列测试函数上的表现全部开源分享希望能帮到正在做多目标优化相关工作的读者。2. 算法整体框架与代码结构设计2.1 NSGA-II主流程NSGA-II的流程可以拆成下面几个步骤初始化种群随机生成N个个体每个个体是一组决策变量。对当前种群执行非支配排序得到若干非支配层级F1, F2, F3...。在每个层级内计算拥挤度距离用于衡量个体在目标空间中的稀疏程度。通过二元锦标赛选择父代执行模拟二进制交叉SBX和多项式变异生成子代种群。父代和子代合并成2N个个体的组合种群重新做非支配排序。从排序后的层级中按F1、F2、F3...的顺序逐个填充下一代种群直到填满N个个体如果某个层级全部放入会超出N就用拥挤度距离从大到小选择该层级中需要的个体数量。重复步骤2到6直到达到最大迭代代数。这里最关键的细节在第5步和第6步父代子代合并后再筛选就是所谓的精英保留策略。它能保证算法在迭代过程中不会丢掉已经找到的优秀解数学上可以证明NSGA-II是收敛的。也正是因为这种机制NSGA-II的收敛速度比第一代NSGA快得多解集的分布也更均匀。2.2 Matlab代码的模块划分在动手写代码之前我先设计了模块结构。每一块我都单独放了独立函数这样调试的时候可以单独测试某一块出问题容易定位。模块函数名功能主循环nsga2_main.m参数设置、种群初始化、迭代控制非支配排序non_dominated_sort.m计算种群中每个个体的支配层级拥挤度计算crowding_distance.m计算同一非支配层内个体的拥挤度锦标赛选择tournament_selection.m基于支配层级和拥挤度的二元锦标赛交叉与变异sbx_crossover.m / polynomial_mutation.m模拟二进制交叉与多项式变异种群合并environmental_selection.m精英保留策略生成下一代种群测试函数zdt1.m / zdt2.m / zdt3.m标准多目标测试函数每个函数都尽量保持输入输出接口清晰以population结构体为核心传递数据。population里存三部分信息决策变量矩阵x、目标函数值矩阵f、以及个体的rank支配层级和crowding拥挤度。这样设计的好处是主循环非常简洁逻辑一眼能看穿。3. Matlab实现的关键细节3.1 种群与目标函数的数据结构定义一个包含全部个体信息的结构体是NSGA-II程序的骨架。我习惯用Matlab的struct数组而不是把变量拆散成多个孤立矩阵——虽然后者性能上可能略微好一点但代码可读性会明显下降。population struct(); population.x zeros(pop_size, var_num); % 决策变量 population.f zeros(pop_size, obj_num); % 目标函数值 population.rank zeros(pop_size, 1); % 非支配层级 population.crowding zeros(pop_size, 1); % 拥挤度距离初始化时决策变量在这个问题给定的上下界之间随机取值。以ZDT1测试函数为例它有30个决策变量每个变量的取值范围是[0, 1]。这样简单的边界条件非常适合新手跑通整个算法实际工程问题里决策变量的边界往往来自物理约束比如刀具转速上下限、装配顺序可行性等。目标函数我单独写成一个函数文件输入是决策变量矩阵输出是目标矩阵。这样后续要换测试函数或者接入自己的工程模型只需要改这个函数即可主程序完全不用动。3.2 快速非支配排序的原理与实现快速非支配排序是NSGA-II的核心机制之一。第一代NSGA对每个个体都和其他所有个体比较支配关系复杂度是O(MN^3)N稍微大点就慢得受不了。NSGA-II把这个复杂度降到了O(MN^2)在N100、200这种常规规模下完全不构成压力。快速排序的核心思想是对于每个个体p维护两个集合——被p支配的个体集合S_p以及支配p的个体数量n_p。然后按下面步骤执行找出所有n_p0的个体这些个体属于第一非支配层F1。对F1中每个个体p遍历S_p中每个个体q将q的n_q减1。如果n_q减到0说明q只被F1中的个体支配放入下一层F2。重复直到所有个体都被分配层级。在Matlab环境中我用向量化方式实现了个体间的支配关系判定避免双重for循环dominates all(f(i,:) f(j,:)) any(f(i,:) f(j,:));需要注意的是这里的约束条件是目标值越小越优。如果你面对的是最大化问题需要在目标函数里取负号转换。这是初学者最容易踩的坑。3.3 拥挤度距离的计算方法拥挤度距离用来评估同一非支配层级内某个个体周围有多“挤”。直观理解就是在目标空间里把同一层级的个体按某个目标维度排好队当前个体前后两个邻居围成的矩形周长越大说明它越稀疏越值得保留。具体算法是对同一层级的个体按每个目标函数值分别排序。边界个体的拥挤度设为无穷大保证边界解一定会被保留。中间个体的拥挤度为相邻两个个体在各目标上的归一化距离之和。f_max max(f(obj_idx, :)); f_min min(f(obj_idx, :)); crowding(sorted_indices(1)) inf; crowding(sorted_indices(end)) inf; for k 2:length(sorted_indices)-1 crowding(sorted_indices(k)) crowding(sorted_indices(k)) ... (f(obj_idx, sorted_indices(k1)) - f(obj_idx, sorted_indices(k-1))) / (f_max - f_min); end这里做了归一化处理目的是让不同目标的量纲不影响拥挤度的比较。如果第一个目标量纲是1e5第二个是0~1不归一化的话拥挤度几乎完全由第一个目标决定第二个目标的分布性就没法保证了。3.4 锦标赛选择中的比较准则选择父代个体时我使用的是二元锦标赛选择。每次从当前种群中随机挑两个个体谁胜出取决于两条规则如果两个个体的rank不同取rank小的层级靠前的更优如果rank相同取拥挤度大的周围更空的个体优先。这个比较准则就是NSGA-II常说的拥挤度比较算子也是它区别于第一代NSGA的关键改进之一。第一代NSGA需要用共享函数参数来维持种群多样性而共享函数对参数σ_share极其敏感调参调得人崩溃。NSGA-II用拥挤度完全替代了这个参数虽然理论上还有一些争议但工程实践中确实省心很多。3.5 SBX交叉与多项式变异的参数选择SBXSimulated Binary Crossover是Deb专门为实数编码遗传算法设计的交叉算子它的特点是子代会以较大概率落在父代附近模拟二进制编码交叉的分布特性。SBX的核心是分布指数distribution index记为eta_c。eta_c越大子代越倾向于和父代接近eta_c越小子代越容易远离父代。我在大多数测试函数上使用eta_c20这是Deb论文中的推荐值。变异算子使用多项式变异分布指数eta_m20变异概率取1/var_num即平均每个个体有一个变量发生变异。一个小技巧是在Matlab中实现SBX交叉时需要先对父代个体处理成对儿交叉利用随机数u决定展开系数betaif rand 0.5 beta (2*u)^(1/(eta_c1)); else beta (1/(2*(1-u)))^(1/(eta_c1)); end child1 0.5*((1beta)*parent1 (1-beta)*parent2); child2 0.5*((1-beta)*parent1 (1beta)*parent2);需要注意这里有一个细节如果直接用rand来决定beta的展开方向会让算出来的子代偏离理想的SBX分布。正确的做法是先生成一个随机数u当u0.5时用第一段公式否则用第二段公式。很多代码这里写得不严谨导致交叉分布出错。3.6 精英保留策略的具体实现精英保留策略是NSGA-II保证收敛性的最后一道保险。每次迭代父代种群P和子代种群Q合并成R P ∪ Q总大小为2N。从R中选出N个个体作为下一代选择方法是先按非支配层级从低到高装入如果一个层级的数量装不满剩余名额该层全部放入如果某个层级F_i全部放入会超出剩余名额则对该层级个体按拥挤度从大到小排序取前面需要的数量。这个策略保证了两个重要性质精英个体层级低的个体永远不会被丢失算法单调收敛同一层内分布性差的个体优先被淘汰维持Pareto前沿的均匀性。我在实际测试中发现精英保留策略让NSGA-II对初始种群的敏感性明显低于其他没有精英机制的算法。初始种群即使偏差很大经过20代左右种群基本能恢复到一份具有良好分布性且收敛到真实Pareto前沿附近的解集。但是如果初始种群全部落在可行域极小的区域那不管怎么样都救不回来——这是所有进化算法的共性。4. 完整代码实现与参数配置4.1 主函数nsga2_main.m下面给出完整的主函数包含所有参数配置。代码结构尽量精简方便移植到不同问题。function nsga2_main() %% 参数配置 pop_size 100; % 种群规模 max_gen 250; % 最大迭代代数 var_num 30; % 决策变量个数 obj_num 2; % 目标函数个数 var_min zeros(1, var_num); % 变量下界 var_max ones(1, var_num); % 变量上界 eta_c 20; % SBX分布指数 eta_m 20; % 多项式变异分布指数 pm 1 / var_num; % 变异概率 %% 初始化种群 population struct(); population.x rand(pop_size, var_num) .* (var_max - var_min) var_min; population.f evaluate_objective(population.x); population.rank zeros(pop_size, 1); population.crowding zeros(pop_size, 1); %% 主循环 for gen 1:max_gen % 计算非支配层级和拥挤度 population non_dominated_sort(population); population crowding_distance(population); % 生成子代 offspring generate_offspring(population, eta_c, eta_m, pm, var_min, var_max); % 合并父代和子代精英保留 combined_pop merge_population(population, offspring); combined_pop non_dominated_sort(combined_pop); combined_pop crowding_distance(combined_pop); population environmental_selection(combined_pop, pop_size); % 显示进度 if mod(gen, 50) 0 fprintf(Generation %d\n, gen); end end %% 输出结果 final_pop environmental_selection(population, pop_size); final_pop non_dominated_sort(final_pop); fprintf(非支配解个数: %d\n, sum(final_pop.rank 1)); figure; plot(final_pop.f(1, final_pop.rank 1), ... final_pop.f(2, final_pop.rank 1), bo); xlabel(f_1); ylabel(f_2); title(NSGA-II Pareto Front (ZDT1)); grid on; end4.2 非支配排序与拥挤度计算非支配排序我这里用的是一种比较简洁的实现。先对每个个体计算被支配次数dominated_count以及它支配的个体集合dominated_set然后逐层剥离层级。function pop non_dominated_sort(pop) N size(pop.x, 1); obj_num size(pop.f, 1); % 初始化支配信息 dominated_count zeros(N, 1); dominated_set cell(N, 1); for i 1:N dominated_set{i} []; end front cell(N, 1); % 每一层的个体索引 front_count 0; % 计算支配关系 for p 1:N for q 1:N if p q continue; end if dominates(pop.f(:, p), pop.f(:, q)) dominated_set{p} [dominated_set{p}, q]; elseif dominates(pop.f(:, q), pop.f(:, p)) dominated_count(p) dominated_count(p) 1; end end if dominated_count(p) 0 pop.rank(p) 1; front{1} [front{1}, p]; end end front_count 1; % 逐层生成 while ~isempty(front{front_count}) next_front []; for i 1:length(front{front_count}) p front{front_count}(i); for j 1:length(dominated_set{p}) q dominated_set{p}(j); dominated_count(q) dominated_count(q) - 1; if dominated_count(q) 0 pop.rank(q) front_count 1; next_front [next_front, q]; end end end front_count front_count 1; front{front_count} next_front; end end function flag dominates(f1, f2) flag all(f1 f2) any(f1 f2); end拥挤度计算的完整实现我会按目标归一化处理这样对于量纲差异很大的目标也能正确比较。另外在环境选择填充种群时如果某层只差几个个体直接截取该层拥挤度最大的那部分不需要再做任何额外操作。4.3 测试函数与Pareto前沿评估我拿ZDT1做了标准测试这个函数的特点是真实Pareto前沿已知是f1取[0,1]均匀分布、f21-sqrt(f1)的曲线。这意味着我可以很方便地计算收敛性指标GDGenerational Distance和分布性指标SPSpacing。写一下评估方式GD世代距离求每个非支配解到真实Pareto前沿的最近距离的均值。GD越小说明算法得到的解越接近真实前沿。SP间距指标求相邻非支配解在目标空间中的欧氏距离的标准差。SP越小说明解集分布越均匀。实测下来在ZDT1上种群规模100、迭代250代GD通常在1e-3到1e-4量级SP在0.01左右。这个表现和直接用Matlab优化工具箱的gamultiobj非常接近但自己写的代码可以灵活修改比如加入约束处理、改成约束多目标这是工具箱很难完全替代的。4.4 代码运行结果分析我用自己的实现跑了ZDT1最终种群在目标空间里的分布如下图所示这里不做图片展示重点说数据。最终解集的f1基本覆盖[0,1]整个区间f2与1-sqrt(f1)的最大偏差控制在0.005以内。种群的非支配解数量通常在80~100之间种群规模100这意味着整个种群几乎全部收敛到Pareto前沿上了。在ZDT2上表现略差一些因为ZDT2的真实Pareto前沿是一个凹曲线在初始阶段很容易陷入局部Pareto前沿。需要把迭代次数提高到400代左右才能得到稳定结果。这说明对于多峰的多目标问题算法参数的调节非常关键特别是交叉分布指数eta_c和变异分布指数eta_m。5. 常见问题与调参经验5.1 种群规模与迭代次数的取舍种群规模直接决定每一代的计算量而迭代次数决定总计算量。这两者需要权衡。种群太小比如30~50种群多样性不够即使跑很多代Pareto前沿也容易出现断档种群太大比如500每代计算时间线性上升但在普通问题上收益并不明显。我的经验是目标函数计算代价低的模型优先加大种群而不是加大迭代次数目标函数计算代价高比如每次都要调用仿真软件那就种群取100把迭代次数拉长配合并行计算使用。并行计算在Matlab里可以用parfor替换主循环中的for但要注意子函数里的随机数种子设置避免并行时产生重复的随机序列。5.2 NSGA-II参数选择的经验值给大家一个参考范围参数推荐值说明种群规模50-200据决策变量个数调整迭代次数200-500据问题复杂度调整SBX分布指数eta_c15-30越大子代越接近父代变异分布指数eta_m20-50越大变异程度越小变异概率pm1/var_num每个变量平均变一次交叉概率pc0.9常规取0.9即可很多文献里交叉概率取0.9变异概率取1/var_num这些是经过大量测试的经验值。如果问题中决策变量之间存在强耦合关系建议把eta_c调小到10左右增加子代探索范围。反之如果算法前期收敛缓慢可以适当增大eta_c让子代更贴近优秀父代。5.3 运行速度慢怎么办Matlab手写NSGA-II最大的坑就是运行速度。尤其是非支配排序里的双重循环N200时就需要4万次比较迭代250代就是1000万次比较。好在Matlab的矩阵运算能力很强可以把支配判定向量化。向量化的思路是构造一个三维兼容矩阵一次性判断某个个体是否支配多个其他个体。虽然逻辑上略显复杂但速度提升非常明显能把日常测试的运行时间从十几分钟压缩到几十秒。另外Matlab中cell数组的拼接操作性能较差频繁使用[dominated_set{p}, q]这种语法会产生大量内存拷贝。如果种群规模超过200建议用预分配的linked list或专门的索引数组来替代cell。5.4 如何验证算法是否正确收敛验证NSGA-II实现的正确性有一个简单有效的方法先用已知真实Pareto前沿的测试函数ZDT系列、DTLZ系列验证。如果算法在ZDT1上输出前沿与理论曲线不匹配通常是以下几个原因支配判定写反了。检查dominates函数里的等号方向。拥挤度归一化除以了零当某目标在当前层级内所有值相等时。需要在分母上加一个极小数epsilon。环境选择里层级填充逻辑出错导致种群被重复个体占满。初始化边界条件写错。比如变量范围是[-1,1]却初始化成了[0,1]。我调试的时候习惯在每代结束后把rank1的个体数量打印出来。如果这个数量快速下降到1说明种群多样性丢失大概率是拥挤度计算出了问题如果一直保持在N附近说明非支配排序或者环境选择有bug。5.5 从测试函数迁移到实际工程问题的注意事项ZDT系列只是验证算法本身真正做工程项目时目标函数往往没有解析表达式每次评估需要调用外部程序或者仿真实验。这时候有几个非常实际的问题第一评估时间成本高。建议在评估函数里加一个评估次数计数器每次评估耗时长的场景把种群规模降到40~60同时配合存档机制把历史所有非支配解都保存下来最终输出时从存档里选。第二约束处理。NSGA-II原生不处理约束标准做法是约束支配法constrained-domination当两个解都满足可行性时按正常支配关系比较当两个解至少有一个不可行时可行解支配不可行解两个都不可行时约束违反程度小的支配大的。这样不需要额外引入惩罚系数调参负担小很多。第三目标函数多的情况。NSGA-II在2-3个目标上表现最好目标数超过5个性能显著下降。这种高维场景建议转向基于分解的方法如MOEA/D或者在NSGA-II中引入参考点机制即NSGA-III。这是另一个话题了但如果你想往高维方向发展可以先理解NSGA-II再迁移到NSGA-III会轻松很多。6. 个人实操总结在Matlab里手写NSGA-II我最大的体会是算法本身并不复杂难的是把数据结构设计好让每一块的职责足够清晰。我第一版代码把所有逻辑挤在主函数里结果调试一个支配关系的bug花了整整一天后来把代码拆成独立模块后每块单独测试半天就定位并修复了所有问题。所以如果你准备自己写这个算法请一定从第一行开始就保持模块化思维。还有一点要提醒的是多目标优化不是银弹。NSGA-II最终给你的是一个Pareto前沿集合需要你或者决策者从中选择折衷方案。在工程落地时我经常用TOPSIS或者层次分析法给Pareto前沿上的解排序这样最终拿出的就不是几十个解而是一个能落地的最优折衷方案。把这一步衔接好算法才算真正用到了实处。如果后续有空我可以再把NSGA-III的实现笔记整理出来它在处理三目标以上问题时优势非常明显。有需要的话评论区留言我看到了会针对性补上。本文还有配套的精品资源点击获取