ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

多目标优化实战:用Python从零实现NSGA-II算法

多目标优化实战:用Python从零实现NSGA-II算法 简介面向数据分析、机器学习与工程设计人员这套Python多目标优化学习包提供了一条从概念到实战的清晰路径。压缩包共5个文件包含2个可直接运行的Python示例脚本、2份PDF理论文档和1段讲解视频整体约293.7MB结构紧凑、便于按需学习。讲解视频系统介绍了多目标问题的基本思想与适用场景重点演示如何基于pymoo、scipy.optimize等库定义目标函数、设置约束条件并运行NSGA-II、MOEA/D等主流算法两个Python脚本则分别以具体算例展示完整的优化流程包括初始化种群、执行迭代优化以及绘制帕累托前沿可直接修改参数用于自己的问题。PDF材料中一份是多目标优化方法综述一份是中文讲解讲义能帮助读者理解各算法的数学原理与适用边界。目前已有201人学习下载适合具备一定Python基础、希望快速上手多目标优化的开发者。通过代码、视频、文献三者结合读者既能掌握理论要点也能获得可复用的脚本模板未来在机器学习超参数调优、工程方案权衡等任务中可以有效实践。1. 多目标优化在现实工程里的位置不是找最优而是找一堆选项很多人刚接触多目标优化时第一反应是“多个目标是不是要加个权重变成单目标”但权重怎么定本身就是个决策问题。成本低和性能高往往冲突你没法提前告诉算法哪边重要。多目标优化换了个思路省掉权重直接通过帕累托支配关系找出一组同样优秀的解再让你在结果里按业务偏好做决定。这就像买房子面积大、总价低、通勤近三个目标很难同时满足多目标优化给你的是所有值得考虑的房源不是唯一答案。这套方法在工程配置、调度、超参数搜索里非常常见用Python写起来也比想象中容易。这篇文章会把多目标优化的核心算法NSGA-II拆开用Python代码从零实现一遍过程中解释支配、非支配排序、拥挤距离这些术语再讲评估、可视化和调参。配套的讲解视频演示的运行节奏和这里的代码顺序一致先定义问题再跑算法最后看Pareto前沿。适合想真正理解多目标优化与决策逻辑、而不是只会调包的开发者。2. 多目标优化的核心概念与NSGA-II的进化逻辑2.1 支配关系、非支配排序、拥挤距离这些术语到底在说什么多目标优化与决策里目标默认是求最小值。解A支配解B当且仅当A在所有目标上都不比B差并且至少在一个目标上严格优于B。假如A在第1个目标等于B在第2个目标比B小那A就支配B。一个解如果不被其他任何解支配它就是一个非支配解全部非支配解组成第一前沿也就是帕累托前沿。去掉第一前沿后剩余解中的非支配解组成第二前沿以此类推。这个过程叫非支配排序它给种群里的个体划分了优劣等级前沿编号越小的解越值得保留。但仅仅按等级排同一前沿内的解没有区别。进化中我们希望这些解尽可能铺满整个前沿而不是挤成一团。拥挤距离就是为此定义的对某个前沿内的解在每个目标方向上把相邻两个解的距离除以该前沿在这个目标上的跨度然后求和。边界解的拥挤距离设为无穷大保证它们一定被保留。在后文代码里你会看到排序和距离计算用的是同一套数据rank和crowding_distance决定了一个个体在精英保留时的优先级rank小的优先rank相同时crowding_distance大的优先。这里有个经常搞错的地方支配关系里说的“不差于”是指每个目标都满足、至少一个目标严格满足。如果你的目标是求最大值那就需要把评估函数取反或者把比较符号反过来否则代码会跑出一堆完全不对的结果。我自己第一次写的时候就因为忘了处理最大化问题Pareto前沿看起来像一个倒三角浪费了半天时间。2.2 从NSGA到NSGA-II为什么快速非支配排序有效NSGA-II全称是非支配排序遗传算法第二代它的核心贡献比第一代多两点快速非支配排序和拥挤距离比较算子。第一代NSGA的时间复杂度是O(MN^3)N是种群规模M是目标数种群一大就跑不动。NSGA-II把复杂度降到了O(MN^2)做法是不再每层都扫描全局比较所有解而是先对每个个体统计被谁支配n值和支配谁s集合然后从n0的个体开始把它们的支配者n值减一减到0的放进下一层就像剥洋葱一样逐层剥离。这个优化让种群在200、500、1000这种规模下都能跑得动。精英保留是另一个关键每一代把父代和子代合并成一个2倍规模的临时种群做完非支配排序后按“rank越小越好rank一样crowding distance越大越好”的规则选出N个进下一代。这样即使这一代进化的结果整体变差了父代里的好个体也不会丢避免了单目标遗传算法里那种“好解被破坏后找不回来”的问题。精英保留配合锦标赛选择让算法在收敛和多样性之间取得了很实用的平衡这也是它成为NSGA-II多目标优化经典入门算法的原因。2.3 其他常见算法与适用场景除了NSGA-II多目标优化领域还有几个主流算法做技术选型时值得了解。MOEA/D把多目标问题拆成多个单目标子问题通过邻域更新维持解的多样性在高维目标数时比NSGA-II表现好NSGA-III沿用非支配排序框架但用分布均匀的参考点替代拥挤距离专门解决3个以上目标时拥挤距离失效的问题SPEA2用强度值和最近邻距离做适应度早期也很流行。如果你做的是电力调度、路径规划这类两个目标的问题NSGA-II的代码和参数调节成本最低我一般会先用它跑通流程。目标数到4个或以上优先换MOEA/D或NSGA-III。可以用一张表对比一下算法核心机制目标数典型场景NSGA-II非支配排序 拥挤距离2-3工程优化、超参数搜索、调度NSGA-III非支配排序 参考点3-10高维多目标如组合设计MOEA/D分解为单目标子问题3-6高维连续优化SPEA2强度支配 最近邻密度2-3算法对比、教学实验注意表里的“目标数”是经验值不是硬限制。下面的实现仍然以NSGA-II为例子因为它的两个核心模块——非支配排序和拥挤距离——是很多多目标优化代码的基础写清楚比直接调包更有用。3. 用Python实现NSGA-II最小可运行代码的参数拆解3.1 定义多目标评估函数和问题参数运行代码之前先确认本地有Python 3.8以上环境不需要安装任何第三方库标准库的random和math就够用。我们以ZDT1测试函数为例它是多目标优化里的经典基准问题决策变量维度n_var30目标个数是2真实Pareto前沿是f2 1 - sqrt(f1)可以用来验证算法结果。import random import math n_var 30 def zdt1(x): f1 x[0] g 1 9 * sum(x[1:]) / (n_var - 1) f2 g * (1 - math.sqrt(f1 / g)) return f1, f2zdt1函数返回两个目标值默认都求最小。g的表达式里用了sum(x[1:])也就是从第2个决策变量到第30个作用是把决策变量对f2的贡献压到一个可控范围让Pareto前沿的形状固定不变。如果换成实际工程问题你只需要改这个函数比如第一目标返回成本第二目标返回耗时后面所有代码不需要动。3.2 种群初始化与编码选择多目标优化代码中种群是一个二维列表每个个体是长度为n_var的实数列表。初始化时每个决策变量在[0,1]区间内均匀随机取值。def init_pop(pop_size): return [[random.random() for _ in range(n_var)] for _ in range(pop_size)]这里的随机采用均匀分布没有对初始解做任何启发式约束。对于带约束的问题常见的做法是直接在初始化时只生成满足约束的个体或者在评估函数里对违反约束的个体返回一个很大的惩罚值。ZDT1无约束所以先不看这个。注意决策变量上下界在ZDT1里都是0和1如果你自己的问题变量范围不同初始化、交叉和变异里的上下界参数都要同步改。3.3 非支配排序与拥挤距离的Python实现这是多目标优化代码的核心。先实现支配判断再实现快速非支配排序。def dominates(x, y): 判断x是否支配y目标值都按最小化处理 one_strict False for a, b in zip(x, y): if a b: return False if a b: one_strict True return one_strict逻辑是只要x在任何一个目标上比y大x就不能支配y如果x在所有目标上都不大于y并且至少有一个目标严格小于y才算支配。注意如果x和y目标值完全相同one_strict保持False返回False这是合理的两个相同的解谁也不能支配谁。非支配排序的代码实现def fast_non_dominated_sort(values): values: 每个个体的目标值列表返回前沿索引列表 pop_size len(values) s [[] for _ in range(pop_size)] n [0] * pop_size fronts [[]] for p in range(pop_size): for q in range(pop_size): if dominates(values[p], values[q]): s[p].append(q) elif dominates(values[q], values[p]): n[p] 1 if n[p] 0: fronts[0].append(p) i 0 while fronts[i]: next_front [] for p in fronts[i]: for q in s[p]: n[q] - 1 if n[q] 0: next_front.append(q) i 1 fronts.append(next_front) return fronts[:-1]这里的s[p]存放被p支配的个体索引n[p]存放支配p的个体数量。先把n[p]0的个体放进第一前沿然后对第一前沿里的每个个体p把它支配的所有q的n值减1减到0的q就是下一前沿。类似拓扑排序全程每个个体只被处理一次所以复杂度是O(MN^2)。返回的fronts是一个列表例如[[3,7],[1,2]]表示第0层有两个个体索引3和7第1层有1和2。拥挤距离的实现如下def crowding_distance(front, values): front: 某一前沿的个体索引列表values: 所有个体的目标矩阵 dist {idx: 0.0 for idx in front} num_obj len(values[0]) for m in range(num_obj): sorted_front sorted(front, keylambda idx: values[idx][m]) dist[sorted_front[0]] float(inf) dist[sorted_front[-1]] float(inf) span values[sorted_front[-1]][m] - values[sorted_front[0]][m] if span 0: continue for i in range(1, len(sorted_front) - 1): dist[sorted_front[i]] ( values[sorted_front[i 1]][m] - values[sorted_front[i - 1]][m] ) / span return dist这段代码里边界个体的距离被设为无穷大确保前沿两端的解不会在精英保留中被淘汰。每个目标的距离贡献都除以该目标的跨度实现了归一化避免成本是10000、时间是0.1这种量纲差异导致时间维度被忽略。如果某个目标在所有前沿个体中的值完全相同span为0直接跳过不影响总距离。3.4 选择、交叉、变异与精英保留实数编码下交叉用模拟二进制交叉SBX变异用多项式变异这两个算子都是从遗传算法里沿用下来的对连续决策变量效果稳定。def sbx_crossover(p1, p2, eta_c20): SBX交叉返回两个子代 c1, c2 p1[:], p2[:] for i in range(n_var): if random.random() 0.5: continue if abs(p1[i] - p2[i]) 1e-14: c1[i] p1[i] c2[i] p2[i] continue u random.random() if u 0.5: beta (2 * u) ** (1 / (eta_c 1)) else: beta (1 / (2 * (1 - u))) ** (1 / (eta_c 1)) c1[i] 0.5 * ((1 beta) * p1[i] (1 - beta) * p2[i]) c2[i] 0.5 * ((1 - beta) * p1[i] (1 beta) * p2[i]) c1[i] max(0.0, min(1.0, c1[i])) c2[i] max(0.0, min(1.0, c2[i])) return c1, c2eta_c称为分布指数控制子代与父代相近的程度值越大子代越接近父代。一般取20搜索不充分时调小到10希望收敛更慢、搜索更广时可以调到30。这里对每个决策变量按0.5的概率决定是否交叉而不是整个向量一次交叉是为了不破坏变量之间的关联性实际使用中可以按你的问题特点调整。多项式变异def poly_mutation(ind, eta_m20, mutate_prob0.1): for i in range(n_var): if random.random() mutate_prob: continue u random.random() delta 0.0 if u 0.5: delta (2 * u) ** (1 / (eta_m 1)) - 1 else: delta 1 - (2 * (1 - u)) ** (1 / (eta_m 1)) ind[i] max(0.0, min(1.0, ind[i] delta)) return indmutate_prob可以理解为每个变量发生变异的概率ZDT1这种30维问题建议取1/n_var左右也就是约0.03太小容易早熟太大就成了随机搜索。锦标赛选择和精英保留def tournament_selection(fronts, dists, pop_size): selected [] layer {} for rank, front in enumerate(fronts): for idx in front: layer[idx] rank for _ in range(pop_size): a, b random.sample(range(pop_size), 2) if layer[a] layer[b]: winner a elif layer[a] layer[b]: winner b else: winner a if dists[a] dists[b] else b selected.append(winner) return selected注意这里的dists是某个前沿的拥挤距离但选择时需要跨前沿比较所以先构造了layer字典记录每个个体所在前沿编号。选择逻辑是先看rankrank小的胜出rank相同则比拥挤距离。这是NSGA-II里的二元锦标赛选择标准。3.5 运行主循环并输出Pareto前沿把上面的模块串起来就是一个最小可运行的NSGA-II多目标优化代码。def nsga2_main(pop_size100, max_gen200, pc0.9, pm0.03, eta_c20, eta_m20): pop init_pop(pop_size) values [zdt1(ind) for ind in pop] for gen in range(max_gen): fronts fast_non_dominated_sort(values) dists {} for front in fronts: dists.update(crowding_distance(front, values)) parents_idx tournament_selection(fronts, dists, pop_size) parents [pop[i] for i in parents_idx] offspring [] while len(offspring) pop_size: p1, p2 random.sample(parents, 2) if random.random() pc: c1, c2 sbx_crossover(p1, p2, eta_c) else: c1, c2 p1[:], p2[:] c1 poly_mutation(c1, eta_m, pm) c2 poly_mutation(c2, eta_m, pm) offspring.append(c1) if len(offspring) pop_size: offspring.append(c2) combined pop offspring combined_values [zdt1(ind) for ind in combined] fronts fast_non_dominated_sort(combined_values) new_pop [] new_values [] for front in fronts: if len(new_pop) len(front) pop_size: for idx in front: new_pop.append(combined[idx]) new_values.append(combined_values[idx]) else: dists crowding_distance(front, combined_values) front_sorted sorted(front, keylambda idx: dists[idx], reverseTrue) need pop_size - len(new_pop) for idx in front_sorted[:need]: new_pop.append(combined[idx]) new_values.append(combined_values[idx]) break pop, values new_pop, new_values fronts fast_non_dominated_sort(values) return [pop[i] for i in fronts[0]], [values[i] for i in fronts[0]]主循环里每代先做非支配排序和拥挤距离计算锦标赛选出父代经过交叉变异生成子代再把父代子代合并成2倍规模种群做精英保留。精英保留的关键在最后一段代码逐个前沿填充新种群如果当前前沿塞不下就按拥挤距离降序取前几个填满。这样保证下一代种群规模始终是pop_size。跑完后pareto_x, pareto_f nsga2_main(pop_size100, max_gen200) print(fPareto解数量: {len(pareto_x)}) for f in pareto_f[:5]: print(f)参数说明表方便对照参数含义常用范围对结果的影响pop_size种群规模50-200太小容易早熟太大计算慢max_gen迭代代数100-500越大越收敛但有上限pc交叉概率0.8-1.0控制探索能力pm变异概率1/n_var左右过大会变成随机搜索eta_cSBX分布指数10-30影响子代对父代的继承程度eta_m多项式变异分布指数10-30影响变异步长从20代开始你会看到Pareto前沿上点逐渐变多到100代后形状基本稳定。如果跑到200代仍然乱七八糟优先检查评估函数和支配关系是否写反。4. 多目标优化代码的调参、评估与可视化别只看调包4.1 种群规模和迭代次数的经验范围很多初次尝试NSGA-II的人会犯一个错误种群设的很大但迭代次数很少结果Pareto前沿稀疏或者种群很小迭代很多最后所有解都挤进了少数几个区域。这两个参数并不是越大越好我一般按问题维度来定维度低于10pop_size50就够维度在30左右像ZDT1这种pop_size100比较稳如果你的评估函数本身要跑好几分钟可以适当减少pop_size同时增加max_gen比如pop_size40, max_gen500以单次计算成本衡量是不是在傻等。迭代次数的判断没有一个固定值更好的办法是每隔若干代记录一次当前第一前沿的超体积指标当连续20代指标提升小于0.1%时就停止。下面会讲超体积怎么算它是比“看一眼点分布”更有量化意义的多目标优化评估方式。4.2 用超体积HV和IGD评估解的分布两个解集不是一个数组不能用均方误差去比。多目标优化里最常用的两个指标是HVHypervolume超体积和IGDInverted Generational Distance反向世代距离。HV测量Pareto前沿与某个参考点之间的体积数值越大说明解的收敛性和均匀性综合越好但需要确定参考点通常取每个目标方向上略差于最差解的点IGD测量真实Pareto前沿上每个点到当前解集的最近距离的平均值需要知道真实前沿的解析式或一组稠密采样点数值越小越好。假设当前前沿的目标值在f1方向从0到1f2方向从0到1参考点可以取(1.1, 1.1)def hypervolume(values, ref_point): 计算两个目标时的HV简单蒙特卡洛采样 count 0 samples 100000 for _ in range(samples): f1 random.uniform(0, ref_point[0]) f2 random.uniform(0, ref_point[1]) dominated False for v in values: if v[0] f1 and v[1] f2: dominated True break if dominated: count 1 return count / samples * ref_point[0] * ref_point[1]蒙特卡洛HV精度不高但代码简单用于判断结果是否改善足够了。如果你想精确计算可以找成熟的库比如pymoo的hv计算函数但自己实现一遍能更好地理解HV在衡量什么如果一个解离真实前沿很远它覆盖的右上角区域就小HV自然低。注意这行代码只适合两个目标三个目标时蒙特卡洛在高维空间里的采样效率会断崖式下降到时可以用WFG工具包计算。指标方向需要真实前沿适用场景HV越大越好否综合评估收敛与分布IGD越小越好是已知真实前沿时对比精度4.3 绘制Pareto前沿的常见问题用Matplotlib画多目标优化解集时最常见的坑是没做归一化就画图导致f1范围是0到1f2范围是0到100时图型被压扁看不出前沿形状。另外颜色映射不要把每个点都画成同一颜色可以用拥挤距离或迭代代数做颜色轴这样能看出哪些位置是算法后期才填充的import matplotlib.pyplot as plt def plot_pareto(fronts, colorsNone): fig, ax plt.subplots() for i, front in enumerate(fronts): xs [f[0] for f in front] ys [f[1] for f in front] if colors is None: ax.scatter(xs, ys, s30, alpha0.8) else: ax.scatter(xs, ys, s30, ccolors, alpha0.8) ax.set_xlabel(f1) ax.set_ylabel(f2) ax.set_title(Pareto Front) ax.grid(True) return ax这里fronts可以传不同代的Pareto前沿列表用来观察算法收敛过程。绘图时如果目标量纲差异大务必先把每个目标缩放到[0,1]再画不然坐标轴刻度和点的分布会误导你对解集质量的判断。4.4 常见坑重复解、过早收敛、目标冲突不明第一类是大量重复解挤在一起。出现原因多数是变异概率pm设得太小加上SBX交叉在多个变量上连续跳过子代几乎等于父代。解决办法是把pm从1/n_var调大一倍或者把eta_c从20调低到15增加搜索步长。第二类过早收敛。典型现象是不到50代Pareto前沿就不再变化但形状明显没达到真实前沿。这时优先检查是否某个目标的值域远大于另一个导致拥挤距离排序里被大值目标的跨度主导。可以考虑对目标值做归一化后再计算拥挤距离但要注意归一化系数恒定否则每次排序的尺度会漂移。第三类是目标之间实际不冲突。如果你的两个目标本质上线性相关比如f1和f2都是x的单调递增函数那Pareto前沿会退化成一条线甚至一个点这时多目标优化没有意义。做优化前先对采样数据算一下Spearman相关系数如果接近1或-1应该改用单目标或降维。5. 多目标优化在超参数搜索和工程配置中的落地技巧5.1 从两个目标到三个目标可视化与决策技巧两个目标可以把Pareto前沿画成散点图三个目标还能画3D图但四个以上目标就没办法直观看了。一个实用的做法是画平行坐标图每条折线代表一个解每个轴代表一个目标目标值都转为越小越好后你会在图中看到成簇的折线。配合颜色和高亮挑选那种在多个目标上都处于中低位、没有特别差值的解。如果实在需要硬选一个可以用TOPSIS或灰色关联分析对Pareto集做排序常见思路是先确定每个目标的正负理想点再算每个解到理想点的距离。5.2 用多目标优化替代网格搜索的代码骨架机器学习超参数搜索里网格搜索按固定步长枚举组合纬度一高就爆炸。多目标优化可以把多个指标作为目标比如把验证集准确率最大化、推理延迟最小化或参数量最小化直接搜Pareto前沿。你只需要把本章前面代码里的zdt1替换成训练评估函数注意评估函数需要返回1到2个目标值。为了省时间可以加上提前停止如果当前子代和父代的HV变化低于阈值就跳出一代。代码骨架大致是def evaluate_pipeline(params): model build_model(params) acc train_and_eval(model) return -acc, params[time_cost]这里返回的每个目标都按最小化处理所以准确率取负号。5.3 一个实用小技巧从Pareto集里自动挑“膝盖解”实际决策时用户通常希望在改善一个目标的同时尽量不损害另一个目标直觉上会去选Pareto前沿曲率最大的位置也就是工程里常说的“膝盖解”。可以用最大曲率或最小距离法来自动挑选先把前沿上的点按f1排序然后对每个点计算它跟前沿两端点连线之间的垂直距离距离最大的点就是膝盖解。def find_knee(pareto_values): # 假设f1升序排列 pts sorted(pareto_values, keylambda v: v[0]) if len(pts) 3: return pts[0] p0 pts[0] p1 pts[-1] dx p1[0] - p0[0] dy p1[1] - p0[1] denom math.hypot(dx, dy) max_dist, knee -1, pts[0] for p in pts[1:-1]: dist abs(dy * (p[0] - p0[0]) - dx * (p[1] - p0[1])) / denom if dist max_dist: max_dist dist knee p return knee这个函数假设两个目标都是求最小所以连线端点是前沿两端。返回的膝盖解是一个兼顾两头的折中方案可以在没有业务偏好时作为默认推荐。本文还有配套的精品资源点击获取
返回列表