MOO-02.多目标优化算法实现

上一篇 MOO-01.帕累托最优 建立了支配关系与帕累托前沿概念。本篇按**「怎么实现」展开:每类算法给出核心步骤、实现要点、pymoo 调用、偏好场景、优势与局限**,并附选型决策图。

段末注释MOEA(Multi-Objective Evolutionary Algorithm,多目标进化算法)通过种群迭代近似帕累托前沿;标量化(scalarization)把向量目标压成标量,复用单目标求解器。

前置MOO-01 帕累托最优 | Math-04 单目标优化 | 贝叶斯优化


一、算法族谱与选型

图 1 多目标优化算法选型

路线 代表算法 典型问题特征
数学规划 / 标量化 加权求和、$\varepsilon$-约束、Tchebycheff 目标可微或 black-box 但评估便宜;权重/约束已知
Pareto 进化算法 NSGA-II、NSGA-III、SPEA2 black-box、组合/混合变量、评估中等成本
分解进化 MOEA/D 目标数 $m \ge 3$,希望前沿分布均匀
代理模型 MOO ParEGO、EHVI、qEHVI 单次评估极贵(实验、大模型训练)

快速选型口诀

  • $m=2$,评估 $\lesssim 10^4$ 次 → NSGA-II
  • $m \ge 3$ → NSGA-IIIMOEA/D
  • 评估 $\lesssim 200$ 次、噪声小 → 多目标贝叶斯
  • 业务已定权重 → 加权求和(最快落地)
  • 需扫凹形/非凸前沿 → $\varepsilon$-约束Tchebycheff,勿单靠加权求和

二、标量化类:把 MOO 降维成单目标

2.1 加权求和(Weighted Sum)

形式(最小化):

$$
\min_{\mathbf{x} \in \mathcal{X}} \ S_w(\mathbf{x}) = \sum_{i=1}^{m} w_i f_i(\mathbf{x}), \quad w_i \ge 0,\ \sum_i w_i = 1
$$

实现步骤

  1. 各目标归一化到相近量级(否则大尺度目标主导):
    $\tilde{f}_i = (f_i - f_i^{\min}) / (f_i^{\max} - f_i^{\min})$
  2. 设定权重 $w$(业务打分卡或网格 ${0.1,0.2,\ldots,0.9}$)
  3. 调用任意单目标求解器(L-BFGS、CMA-ES、Optuna 等)
  4. 若要近似前沿:对多组 $w$ 重复运行,取所有解的非支配子集
1
2
3
4
5
6
7
8
9
10
11
12
13
import numpy as np
from scipy.optimize import differential_evolution

def weighted_sum(x, w, f_funcs):
"""f_funcs: list of callable, 均最小化。"""
return sum(wi * fi(x) for wi, fi in zip(w, f_funcs))

w = np.array([0.6, 0.4])
bounds = [(0, 1)] * 5
res = differential_evolution(
lambda x: weighted_sum(x, w, [f1, f2]),
bounds, seed=42, maxiter=200,
)
维度 说明
偏好场景 权重明确;凸帕累托前沿;快速 MVP(如 EVOLVEpro 线性加权)
优势 实现极简;可复用成熟单目标库;解释性强(「60% 看 A,40% 看 B」)
局限 无法得到凹前沿上的解;权重与解非线性对应,扫 $w$ 分布不均;多目标量纲未归一化时结果失真

段末注释:帕累托前沿时,加权求和可达边界上任意点;非凸时凹段上的折衷解永远不是任何 $w$ 的最优解。


2.2 $\varepsilon$-约束法(Epsilon-Constraint)

形式:选定主目标 $f_1$,其余转为约束:

$$
\min_{\mathbf{x} \in \mathcal{X}} f_1(\mathbf{x}) \quad \text{s.t.} \quad f_i(\mathbf{x}) \le \varepsilon_i,\ i = 2,\ldots,m
$$

实现步骤

  1. 估计各 $f_i$ 的可行范围(单目标极值或先验)
  2. 对 $f_2,\ldots,f_m$ 在范围内网格化 $\varepsilon_i$
  3. 对每个 $\boldsymbol{\varepsilon}$ 解一个约束单目标问题
  4. 合并所有可行解,做非支配筛选 → 近似前沿
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
from pymoo.core.problem import Problem
from pymoo.algorithms.soo.nonconvex.de import DE
from pymoo.optimize import minimize

class EpsilonConstraint(Problem):
def __init__(self, eps2):
super().__init__(n_var=2, n_obj=1, n_ieq_constr=1, xl=0, xu=1)
self.eps2 = eps2

def _evaluate(self, x, out, *args, **kwargs):
f1 = x[:, 0] ** 2
f2 = (x[:, 0] - 1) ** 2 + x[:, 1] ** 2
out["F"] = f1.reshape(-1, 1)
out["G"] = (f2 - self.eps2).reshape(-1, 1) # f2 <= eps2

front = []
for eps2 in np.linspace(0.05, 2.0, 30):
res = minimize(EpsilonConstraint(eps2), DE(), ("n_gen", 100), verbose=False)
front.append(res.F[0, 0])
维度 说明
偏好场景 需覆盖非凸前沿;有 SLA 型硬约束(如延迟 $\le 100$ ms 再 min loss)
优势 理论上可扫完整帕累托集(离散网格下);主目标语义清晰
局限 每点一次完整优化,$m$ 大时网格爆炸;$\varepsilon$ 选不好会得到不可行子问题

2.3 切比雪夫标量化(Weighted Tchebycheff)

形式($z_i^*$ 为理想点或各目标下界):

$$
\min_{\mathbf{x} \in \mathcal{X}} \ \max_{1 \le i \le m} \ w_i \left| f_i(\mathbf{x}) - z_i^* \right|
$$

实现要点

  • $w_i > 0$ 保证每个 Pareto 解可被某个 $w$ 命中(含非凸段)
  • 对 $\max$ 用 introduce auxiliary variable $t$:$\min t$ s.t. $w_i|f_i - z_i^*| \le t$
  • 常与归一化目标联用
维度 说明
偏好场景 非凸前沿;希望「最坏目标差距」最小化(均衡型决策)
优势 比加权求和更能覆盖非凸帕累托集
局限 对 $z^*$ 敏感;$\max$ 结构不可微,需进化或 smooth approximation

三、NSGA-II:二目标默认首选

前置:NSGA-II 是遗传算法(GA)的 Pareto 扩展。选择、交叉(SBX)、变异(PM)、种群迭代等基础见 MOO-03 遗传算法

图 2 NSGA-II 主循环

Non-dominated Sorting Genetic Algorithm II(Deb et al., 2002)——工程上 $m=2$ 时最常用。

3.1 核心机制

组件 作用
快速非支配排序 $O(MN^2)$ 将种群分为 Front 1(非支配)、Front 2…
拥挤距离 同 Front 内,边界解距离设 $\infty$;中间解为相邻目标值差之和
精英保留 合并 $P_t \cup Q_t$,按 Rank → Crowding 选前 $N$ 个
遗传算子 模拟二进制交叉 SBX + 多项式变异 PM

拥挤距离(第 $i$ 个目标上对解 $k$):

$$
d_k += \frac{f_i^{(k+1)} - f_i^{(k-1)}}{f_i^{\max} - f_i^{\min}}
$$

(边界解 $d = \infty$,保证极端折衷不被淘汰。)

3.2 伪代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
输入: 种群规模 N, 代数 G_max
P0 ← 随机初始化 N 个个体
for t = 0 .. G_max-1:
Qt ← TournamentSelect(Pt) → SBX + PM
Rt ← Pt ∪ Qt // 规模 2N
(F1, F2, ...) ← FastNonDominatedSort(Rt)
Pt+1 ← ∅
for each Front Fi:
if |Pt+1| + |Fi| <= N:
Pt+1 ← Pt+1 ∪ Fi
else:
CrowdingDistance(Fi)
Pt+1 ← Pt+1 ∪ TopByCrowding(Fi, N - |Pt+1|)
break
输出: Front 1 作为近似帕累托集

3.3 pymoo 实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.operators.crossover.sbx import SBX
from pymoo.operators.mutation.pm import PM
from pymoo.operators.sampling.rnd import FloatRandomSampling
from pymoo.optimize import minimize

algorithm = NSGA2(
pop_size=100,
sampling=FloatRandomSampling(),
crossover=SBX(prob=0.9, eta=15),
mutation=PM(eta=20),
eliminate_duplicates=True,
)
res = minimize(problem, algorithm, ("n_gen", 200), seed=1, verbose=False)
# res.X 决策变量, res.F 目标值, res.algorithm.opt 非支配 archive

离散/组合问题:将 FloatRandomSampling 换为 IntegerRandomSampling 或自定义 Sampling;交叉/变异按编码改(如排列用 OX)。

3.4 超参数建议

参数 典型值 说明
pop_size $50$–$200$ 前沿越长/pop 越大
n_gen $100$–$500$ 直到 HV 曲线平台
SBX eta $5$–$20$ 小 → 子代更接近父代
PM eta $10$–$30$ 变异步长控制
维度 说明
偏好场景 $m=2$ 或 $m=3$;black-box;混合/连续变量;评估 $10^3$–$10^5$ 次可接受
优势 实现成熟、文献/库丰富;拥挤距离在 2D/3D 前沿分布好;无需预设权重
局限 $m \ge 4$ 时拥挤距离失效(全解拥挤度趋同)→ 换 NSGA-III;高维决策空间收敛慢;不保证全局最优

段末注释SBX(Simulated Binary Crossover)模拟单点交叉的分布性质;PM(Polynomial Mutation)在变量边界内多项式扰动。


四、NSGA-III:Many-Objective($m \ge 3$)

图 3 NSGA-III 参考点机制

当目标 $\ge 4$ 时,几乎所有解互不支配,NSGA-II 的 Rank 失去区分力;NSGA-III结构化参考点维持多样性。

4.1 核心机制

  1. 参考点生成:在归一化目标单纯形上用 Das-Dennis 法均匀布点(如 $m=3$, $p=12$ 分割 → 91 个参考点)
  2. 目标归一化:用理想点 + 极值点估计,将 $F$ 映到超平面
  3. 关联:每个解关联到最近参考点(垂直距离或角度)
  4. Niche 保留:选解时优先「参考点覆盖数少」的区域,而非拥挤距离

4.2 实现要点

1
2
3
4
5
6
from pymoo.algorithms.moo.nsga3 import NSGA3
from pymoo.util.ref_dirs import get_reference_directions

ref_dirs = get_reference_directions("das-dennis", 3, n_partitions=12) # 3 目标
algorithm = NSGA3(ref_dirs=ref_dirs, pop_size=len(ref_dirs) * 10)
res = minimize(problem, algorithm, ("n_gen", 300), seed=1)

pop_size 建议 $\ge 4 \times$ 参考点数量,否则部分参考点无关联解。

维度 说明
偏好场景 $m \ge 3$,尤其 $m \ge 4$(many-objective);需均匀覆盖高维前沿
优势 高维目标下多样性可控;参考点可视化辅助理解决分布
局限 参考点数量随 $m,p$ 指数增长;$m=2$ 时不如 NSGA-II 简洁;对归一化质量敏感

五、SPEA2:外部 Archive + 强度适应度

Strength Pareto Evolutionary Algorithm 2 维护固定规模外部档案(archive),用强度适应度(被支配解数量)+ 密度($k$-NN 距离)选下一代。

5.1 实现流程

  1. 合并种群 $P_t$ 与档案 $A_t$
  2. 计算每个个体 $i$ 的强度 $S(i)$ = 支配的个体数
  3. 原始适应度 $R(i) = \sum_{j \in \text{支配 } i} S(j)$(被强者支配的强度和)
  4. 密度 $D(i)$ = 第 $k$ 近邻距离($k \approx \sqrt{N+N’}$)
  5. 适应度 $F(i) = R(i) + D(i)$;选 $F$ 最小的 $N$ 个进档案
  6. 档案满时:截断——反复删除密度最高(最拥挤)个体
1
2
3
4
from pymoo.algorithms.moo.spea2 import SPEA2

algorithm = SPEA2(pop_size=100)
res = minimize(problem, algorithm, ("n_gen", 200), seed=1)
维度 说明
偏好场景 需要稳定、有限规模的非支配 archive;$m=2,3$ 均可
优势 Archive 规模可控;理论性质清晰;某些问题上收敛性优于 NSGA-II
局限 每代 $k$-NN 计算 $O(N^2)$;截断策略对参数敏感;工程生态略小于 NSGA-II

六、MOEA/D:分解 + 邻域协作

图 4 MOEA/D 分解思想

Multi-Objective Evolutionary Algorithm based on Decomposition 将 MOO 拆成 $N$ 个标量子问题,每个子问题绑定一个权重向量 $\boldsymbol{\lambda}^{(j)}$,只与邻域内子问题共享信息。

6.1 子问题形式

常用 Tchebycheff 分解

$$
g^{\mathrm{te}}(\mathbf{x} \mid \boldsymbol{\lambda}, z^) = \max_{1 \le i \le m} \ \lambda_i \left| f_i(\mathbf{x}) - z_i^ \right|
$$

加权切比雪夫 / PBI(Penalty Boundary Intersection)。

6.2 实现流程

1
2
3
4
5
6
7
8
9
10
11
生成 N 个均匀权重向量 λ(1)..λ(N)
对每个 j: 定义邻域 B(j) = {与 λ(j) 距离最近的 T 个权重}
初始化: 对每个 j 一个解 x(j),评估 F(x(j))
for gen = 1 .. G_max:
for j = 1 .. N:
从 B(j) 中选父代,交叉变异得 y
更新理想点 z* = min(z*, F(y))
for k in B(j):
if g_te(F(y)|λ(k), z*) < g_te(F(x(k))|λ(k), z*):
x(k) ← y
输出: {x(j)} 的非支配子集
1
2
3
4
5
6
7
8
9
10
from pymoo.algorithms.moo.moead import MOEAD
from pymoo.util.ref_dirs import get_reference_directions

ref_dirs = get_reference_directions("das-dennis", 3, n_partitions=12)
algorithm = MOEAD(
ref_dirs=ref_dirs,
n_neighbors=20,
prob_neighbor_mating=0.7,
)
res = minimize(problem, algorithm, ("n_gen", 200), seed=1)
维度 说明
偏好场景 $m \ge 3$;前沿较长且希望均匀;可并行化各子问题
优势 单次迭代 $O(N)$ 邻域更新,可扩展;分解直观;many-objective 表现常优于 NSGA-II
局限 非凸前沿需特殊分解(如 Tchebycheff);权重固定,难动态改偏好;离散问题需定制交叉

段末注释PBI 分解在权重方向与距离方向之间加惩罚项,改善 MOEA/D 对非凸前沿的覆盖。


七、多目标贝叶斯优化:昂贵 black-box

单次评估耗时分钟~天时(湿实验、训练一轮模型),进化算法样本效率不够。多目标贝叶斯优化(Multi-Objective Bayesian Optimization,MOBO)用 高斯过程(Gaussian Process,GP)建模每个目标,按采集函数选下一批实验点。

7.1 ParEGO:随机标量化 + GP

思路:每轮随机抽权重 $w$,标量化 $S_w(\mathbf{x})$,对 $S_w$ 建单目标 GP,用 EI 选下一点;多轮后取所有观测的非支配集。

维度 说明
偏好场景 $m=2,3$;评估 $\lesssim 100$–$300$ 次;实现要快
优势 简单;复用单目标 BO 代码
局限 随机 $w$ 覆盖前沿不稳定;不直接优化 HV
1
2
3
4
5
6
7
8
9
10
11
12
13
# 概念示例;生产可用 botorch MultiObjectiveBayesianOptimization
import numpy as np
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import Matern

def parego_step(X_obs, F_obs, bounds, rng):
w = rng.dirichlet(np.ones(F_obs.shape[1]))
F_norm = (F_obs - F_obs.min(0)) / (F_obs.ptp(0) + 1e-8)
y = F_norm @ w
gp = GaussianProcessRegressor(kernel=Matern(nu=2.5), normalize_y=True)
gp.fit(X_obs, y)
# 在网格/随机候选上最大化 EI → 返回 x_next
...

7.2 EHVI / qEHVI:直接最大化超体积增益

Expected Hypervolume Improvement(EHVI)衡量「加入候选点 $\mathbf{x}$ 后,超体积期望增量」:

$$
\alpha_{\mathrm{EHVI}}(\mathbf{x}) = \mathbb{E}\left[ \mathrm{HV}(P \cup {\mathbf{f}(\mathbf{x})}) - \mathrm{HV}(P) \right]
$$

qEHVI 一次选 $q$ 个点(批量实验)。Botorch 实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
# pip install botorch gpytorch
from botorch.models.gp_regression import SingleTaskGP
from botorch.models.model_list_gp_regression import ModelListGP
from botorch.utils.multi_objective.box_decompositions.non_dominated import (
NondomDominatedPartitioning,
)
from botorch.acquisition.multi_objective.monte_carlo import (
qExpectedHypervolumeImprovement,
)

# train_data: (n, d), train_obj: (n, m), 均最小化
models = ModelListGP(
*[SingleTaskGP(train_X, train_obj[:, i:i+1]) for i in range(m)]
)
ref_point = train_obj.max(0).values + 0.1 # 参考点,需劣于所有点
partitioning = NondomDominatedPartitioning(ref_point, train_obj)
acq = qExpectedHypervolumeImprovement(models, ref_point, partitioning)
# 优化 acq 得下一批 candidate
维度 说明
偏好场景 评估极少($<200$);目标平滑、噪声可控;需 principled 探索-利用
优势 样本效率高;直接优化与 MOO 一致的 HV 指标;支持批量 q
局限 $m \ge 4$ 时 GP 与 HV 计算昂贵;对噪声/失败实验需 robust 变体;实现复杂度高

7.3 与单目标 BO 对比

单目标 BO MOBO (EHVI)
代理模型 1 个 GP $m$ 个 GP 或多输出 GP
采集函数 EI, UCB EHVI, NEHVI, ParEGO
输出 1 个最优 $\mathbf{x}$ 非支配集 + 近似前沿
典型库 scikit-optimize, Optuna+GP BoTorch, Emukit

八、算法横向对比

算法 实现复杂度 样本/评估效率 目标数 $m$ 非凸前沿 并行 主要局限
加权求和 中(多次单目标) 任意 ✗ 凹段缺失 ✓ 各 $w$ 独立 权重难选、非凸盲区
$\varepsilon$-约束 ★★ 2–4 ✓ 各 $\varepsilon$ 网格爆炸
Tchebycheff ★★ 任意 理想点敏感
NSGA-II ★★ 低–中 2–3 ✓ 种群 $m\ge4$ 多样性差
NSGA-III ★★★ 低–中 3+ 参考点指数增
SPEA2 ★★★ 低–中 2–3 部分 $O(N^2)$ 密度
MOEA/D ★★★ 低–中 3+ △ 依赖分解 ✓✓ 子问题 非凸需 PBI
ParEGO ★★ 2–3 序贯 前沿覆盖随机
qEHVI ★★★★ 最高 2–4 批量 q 实现/算力

九、工程落地 Checklist

  1. 明确 $m$ 与评估预算 $N_{\mathrm{eval}}$:决定进化 vs BO
  2. 目标归一化:否则算法被大数量纲目标绑架
  3. 约束处理:pymoo 用 n_ieq_constr / n_eq_constr;不可行解赋惩罚或约束支配
  4. 重复运行:MOEA 随机性大,建议 $\ge 5$ 次独立 seed,合并非支配集报告 HV 均值±std
  5. 终止条件:固定 n_gen HV 相对变化 $< \delta$ 连续 $K$ 代
  6. 验证:在 ZDT/DTLZ 基准上先跑通,再换业务 Problem
1
2
3
4
5
from pymoo.indicators.hv import Hypervolume

ref = np.max(res.F, axis=0) * 1.1 # 参考点劣于所有解
hv = Hypervolume(ref_point=ref).do(res.F)
print("Hypervolume:", hv)

十、与本仓库场景的映射

业务 推荐算法 理由
酶/抗体多指标(实验贵) qEHVI / ParEGO 每轮湿实验成本高
Embedding Recall–延迟 NSGA-II 两目标、评估中等(跑 benchmark)
超参(lr, batch, depth) NSGA-II 或 Optuna+加权 评估次数中等;权重已知时可标量化
调度(等待/能耗/公平) NSGA-III 或 MOEA/D $m=3$
已有 EVOLVEpro 加权分 加权求和 + 事后非支配筛选 兼容现有 pipeline

十一、小结

算法 一句话
加权求和 最快 MVP,非凸前沿不完整
$\varepsilon$-约束 扫非凸,网格成本换完备性
NSGA-II $m=2$ 默认;排序 + 拥挤距离
NSGA-III $m\ge3$;参考点 niche
SPEA2 外部 archive + 强度适应度
MOEA/D 分解子问题 + 邻域;高维均匀前沿
ParEGO / EHVI 昂贵评估;GP + 采集函数

推荐阅读顺序MOO-01MOO-03 遗传算法 → 本文 → 业务案例(规划中)。

段末注释HV(Hypervolume)衡量前沿与参考点围成超体积,是 MOO 最常用的质量指标;qEHVI 为其期望改进的批量版本,是 BoTorch 多目标 BO 默认强基线。

-------------本文结束感谢您的阅读-------------