1 引入
收敛性分析可以帮助我们在再次解决同一个问题时，更好地决定终止标准。收敛性分析应当考虑两种情况：
1）Pareto前沿未知；
2）Pareto前沿通过解析推导而来，或者存在合理的近似值。
一个双目标优化示意如https://inkiyinji.blog.csdn.net//article/details/125914306，其Pareto前沿和解如下：
代码如下：
from pymoo.util.misc import stack
class MyTestProblem(MyProblem):
def _calc_pareto_front(self, flatten=True, *args, **kwargs):
f2 = lambda f1: ((f1/100) ** 0.5 - 1)**2
F1_a, F1_b = np.linspace(1, 16, 300), np.linspace(36, 81, 300)
F2_a, F2_b = f2(F1_a), f2(F1_b)
pf_a = np.column_stack([F1_a, F2_a])
pf_b = np.column_stack([F1_b, F2_b])
return stack(pf_a, pf_b, flatten=flatten)
def _calc_pareto_set(self, *args, **kwargs):
x1_a = np.linspace(0.1, 0.4, 50)
x1_b = np.linspace(0.6, 0.9, 50)
x2 = np.zeros(50)
a, b = np.column_stack([x1_a, x2]), np.column_stack([x1_b, x2])
return stack(a,b, flatten=flatten)
import matplotlib.pyplot as plt
problem = MyTestProblem()
pf_a, pf_b = problem.pareto_front(use_cache=False, flatten=False)
pf = problem.pareto_front(use_cache=False, flatten=True)
plt.figure(figsize=(7, 5))
plt.scatter(F[:, 0], F[:, 1], s=30, facecolors='none', edgecolors='b', label="Solutions")
plt.plot(pf_a[:, 0], pf_a[:, 1], alpha=0.5, linewidth=2.0, color="red", label="Pareto-front")
plt.plot(pf_b[:, 0], pf_b[:, 1], alpha=0.5, linewidth=2.0, color="red")
plt.title("Objective Space")
plt.legend()
plt.show()
示意中开启了一个历史记录选项save_history，这正是收敛性分析所需要的。当开启save_history后，其将记录每次迭代的算法对象并存储于Result对象中。该方法更占内存 (特别是对于多次迭代)，但具有可以事后分析任何算法相关变量的优点。
为了进行分析，算法的一些中间步骤需要被考虑。它们已经被以下实现：
1）Callback类存储了算法每次迭代的必要信息；
2）当使用最小化方法时开启save_history以记录。
此外使用Hypervolume和IGD来衡量算法的性能：
from pymoo.optimize import minimize
res = minimize(problem,
algorithm,
("n_gen", 40),
seed=1,
save_history=True,
verbose=False)
X, F = res.opt.get("X", "F")
hist = res.history
print(len(hist))
从history可以导出我们所需要的信息：
n_evals = [] # 函数评估次数
hist_F = [] # 每次迭代的目标值
hist_cv = [] # 每次迭代的约束违反
hist_cv_avg = [] # 所有种群的平均约束违反
for algo in hist:
# 存储评估数
n_evals.append(algo.evaluator.n_eval)
# 从算法中检索最优
opt = algo.opt
# 存储约束违反
hist_cv.append(opt.get("CV").min())
hist_cv_avg.append(algo.pop.get("CV").mean())
# 过滤出可行解，并添加至目标空间
feas = np.where(opt.get("feasible"))[0]
hist_F.append(opt.get("F")[feas])
2 约束满足
首先，可以按如下查询每一代可行解的数量：
k = np.where(np.array(hist_cv) <= 0.0)[0].min()
print(f"经过{n_evals[k]}次评估，第{k}代至少有1个可行解")
输出如下：
经过40次评估，第0代至少有1个可行解
因为这个问题的复杂性低，所以很快就找到了可行的解决方案。然而，对应于自己的优化问题可能完全不同，因此需要先行分析：
代码如下：
# 分析最不可行的最优解而非总体时，将’hist_cv_avg‘替换为’hist_cv‘
vals = hist_cv_avg
k = np.where(np.array(vals) <= 0.0)[0].min()
plt.figure(figsize=(7, 5))
plt.plot(n_evals, vals, color='black', lw=0.7, label="Avg. CV of Pop")
plt.scatter(n_evals, vals, facecolor="none", edgecolor='black', marker="p")
plt.axvline(n_evals[k], color="red", label="All Feasible", linestyle="--")
plt.title("Convergence")
plt.xlabel("Function Evaluations")
plt.ylabel("Constraint Violation")
plt.legend()
plt.show()
3 Pareto前沿未知
如果Pareto前沿未知，则无法知道算法是否已经收敛到真正的最优值。当然也不是没有任何进一步的信息。因为依然可以看到算法在优化过程中何时取得了大部分进展，并确定迭代次数是否增减。此外，以下指标用于将两种算法相互比较：1）Hypervolume和2）IGD。
3.1 Hypervolume (HV)
Hypervolume是多目标优化中的一个很重要的指标，是Pareto兼容的，需要基于预定义参考点和提供的解决方案之间的volume。因此，首先定义参考点ref_point，其应当大于Pareto前沿的最大值：
approx_ideal = F.min(axis=0)
approx_nadir = F.max(axis=0)
from pymoo.indicators.hv import Hypervolume
metric = Hypervolume(ref_point= np.array([1.1, 1.1]),
norm_ref_point=False,
zero_to_one=True,
ideal=approx_ideal,
nadir=approx_nadir)
hv = [metric.do(_F) for _F in hist_F]
plt.figure(figsize=(7, 5))
plt.plot(n_evals, hv, color='black', lw=0.7, label="Avg. CV of Pop")
plt.scatter(n_evals, hv, facecolor="none", edgecolor='black', marker="p")
plt.title("Convergence")
plt.xlabel("Function Evaluations")
plt.ylabel("Hypervolume")
plt.show()
输出如下：
Hypervolume可以针对2–3个目标进行精确有效的计算，然而Hypervolume的时间成本随着维数的增加而增高。因此对于更高的维度，一些使用近似方法，这在Pymoo中尚不可用。
3.2 运行指标
当真正的Pareto前沿未知时，另一种分析运行的方法是新近的运行指标。 运行指标显示了每一代目标空间的差异。如果没有定义默认的终止标准，这个指标也被用于Pymoo来确定终止条件。
例如，该分析表明算法从第4代到第5代存在显着改进：
代码如下：
from pymoo.util.running_metric import RunningMetric
running = RunningMetric(delta_gen=5,
n_plots=3,
only_if_n_plots=True,
key_press=False,
do_show=True)
for algorithm in res.history[:15]:
running.notify(algorithm)
设置更多的迭代次数时：
代码如下：
from pymoo.util.running_metric import RunningMetric
running = RunningMetric(delta_gen=10,
n_plots=4,
only_if_n_plots=True,
key_press=False,
do_show=True)
for algorithm in res.history:
running.notify(algorithm)
4 Pareto前沿已知或者近似
4.1 IGD/GD/IGD+/GD+
对于一个问题，其Pareto前沿可以手动提供，也可以直接在问题定义中实现，以动态分析运行。这里展示了一个使用算法历史作为附加后处理步骤的示例。
对于现实世界问题，必须使用近似值。可以通过运行算法几次并从所有解集中提取非支配解来获得近似值。如果只有一次运行，另一种方法是使用获得的非支配解集作为近似值。但是，结果仅表明算法在收敛到最终集合方面的进展程度。
代码如下：
from pymoo.indicators.igd import IGD
metric = IGD(pf, zero_to_one=True)
igd = [metric.do(_F) for _F in hist_F]
plt.plot(n_evals, igd, color='black', lw=0.7, label="Avg. CV of Pop")
plt.scatter(n_evals, igd, facecolor="none", edgecolor='black', marker="p")
plt.axhline(10**-2, color="red", label="10^-2", linestyle="--")
plt.title("Convergence")
plt.xlabel("Function Evaluations")
plt.ylabel("IGD")
plt.yscale("log")
plt.legend()
plt.show()
输出如下：
代码如下：
from pymoo.indicators.igd_plus import IGDPlus
metric = IGDPlus(pf, zero_to_one=True)
igd = [metric.do(_F) for _F in hist_F]
plt.plot(n_evals, igd, color='black', lw=0.7, label="Avg. CV of Pop")
plt.scatter(n_evals, igd, facecolor="none", edgecolor='black', marker="p")
plt.axhline(10**-2, color="red", label="10^-2", linestyle="--")
plt.title("Convergence")
plt.xlabel("Function Evaluations")
plt.ylabel("IGD+")
plt.yscale("log")
plt.legend()
plt.show()
输出如下：
扫一扫
1387
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
