Pymoo学习 (3):使用多目标优化找到最优解的集合

865 阅读4分钟

1 双目标优化

image.png

首先需要将其转换为优化问题的标准形式:
1)f2(x)-f_2(x)替换f2(x)f_2(x)
2)不等式约束均转换为0\leq0的形式,即g2(x)0-g_2(x)\leq0
3)标准化变量,使得其重要性等同,即g1(x)/(2(0.1)(0.9))g_1(x)/(2*(-0.1)*(-0.9))g2(x)/(20(0.4)(0.6))g_2(x)/(20*(-0.4)*(-0.6))。转换结果如下:

image.png

2 Python求解

  本小节使用element-wise求解,其是三种问题解决方案之一 :

import numpy as np
from pymoo.core.problem import ElementwiseProblem


class MyProblem(ElementwiseProblem):

    def __init__(self):
        super(MyProblem, self).__init__(
            n_var=2,  # 变量数量
            n_obj=2,  # 优化目标的数量
            n_constr=2,  # 约束的数量
            xl=np.array([-2, -2]),  # 变量的下限
            xu=np.array([2, 2])  # 变量的上限
        )

    def _evaluate(self, x, out, *args, **kwargs):
        """
        :param x:       1D向量,长度为n_var
        :param out:
        :param args:
        :param kwargs:
        :return:
        """
        # 目标函数
        f1 = 100 * (x[0]**2 + x[1]**2)
        f2 = (x[0] - 1)**2 + x[1]**2

        # 约束条件
        g1 = 2 * (x[0] - 0.1) * (x[0] - 0.9) / 0.18
        g2 = - 20 * (x[0] - 0.4) * (x[0] - 0.6) / 4.8

        # 存放目标值和约束值
        out["F"] = [f1, f2]
        out["G"] = [g1, g2]

  注:三种问题求解方式如下:
1)element-wise:每次为xx 调用_evaluate;
2)vectorizedxx 表示求解方案的集合
3)functional:将每个目标和约束作为函数提供。

2.1 初始化

  示意中的优化问题很简单,只有两个目标和两个约束,以下为使用著名多目标算法NSGA-II的初始化:

from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.factory import get_sampling, get_crossover, get_mutation


algorithm = NSGA2(
    pop_size=40,  # 种群数量
    n_offsprings=10,  # 后代数量
    sampling=get_sampling("real_random"),  # 采样的方式,如随机采样
    crossover=get_crossover("real_sbx", prob=0.9, eta=15),  # 交配方式,如模拟二进制
    mutation=get_mutation("real_pm", eta=20),  # 变异方式,多项式变异
    eliminate_duplicates=True   # 确保后代的目标值不同
)

123456789101112

2.2 终止条件

  在优化之前需要定义终止条件。常见方法包括:
1)限制函数评估的总数;
2)限制算法的迭代次数;
3)定制:一些算法实现了自己的终止条件,例如单纯形退化时的Nelder-Mead或使用供应商库的CMA-ES。
由于本例较简单,使用40次迭代作为终止条件:

from pymoo.factory import get_termination

termination = get_termination("n_gen", 40)

2.3 优化

  最小化方法如下:

from pymoo.optimize import minimize

res = minimize(
    problem=MyProblem(),
    algorithm=algorithm,
    termination=termination,
    seed=1,
    save_history=True,
    verbose=True    # 输出控制变量
)

# 存储可行解
X = res.X
# 存储目标值
F = res.F

  输出如下,每一列分别表示:
1)迭代次数;
2)评估次数;
3)总体最小约束违规;
4)总体平均约束违规;
5)非支配解决方案的数量;
6)和7)搜索空间中运动相关指标。

=====================================================================================
n_gen |  n_eval |   cv (min)   |   cv (avg)   |  n_nds  |     eps      |  indicator  
=====================================================================================
    1 |      40 |  0.00000E+00 |  2.36399E+01 |       1 |            - |            -
    2 |      50 |  0.00000E+00 |  1.15486E+01 |       1 |  0.00000E+00 |            f
    3 |      60 |  0.00000E+00 |  5.277918607 |       1 |  0.00000E+00 |            f
    4 |      70 |  0.00000E+00 |  2.406068542 |       2 |  1.000000000 |        ideal
    5 |      80 |  0.00000E+00 |  0.908316880 |       3 |  0.869706146 |        ideal
    6 |      90 |  0.00000E+00 |  0.264746300 |       3 |  0.00000E+00 |            f
    7 |     100 |  0.00000E+00 |  0.054063822 |       4 |  0.023775686 |        ideal
    8 |     110 |  0.00000E+00 |  0.003060876 |       5 |  0.127815454 |        ideal
    9 |     120 |  0.00000E+00 |  0.00000E+00 |       6 |  0.085921728 |        ideal
   10 |     130 |  0.00000E+00 |  0.00000E+00 |       7 |  0.015715204 |            f
   11 |     140 |  0.00000E+00 |  0.00000E+00 |       8 |  0.015076323 |            f
   12 |     150 |  0.00000E+00 |  0.00000E+00 |       7 |  0.026135665 |            f
   13 |     160 |  0.00000E+00 |  0.00000E+00 |      10 |  0.010026699 |            f
   14 |     170 |  0.00000E+00 |  0.00000E+00 |      11 |  0.011833783 |            f
   15 |     180 |  0.00000E+00 |  0.00000E+00 |      12 |  0.008294035 |            f
   16 |     190 |  0.00000E+00 |  0.00000E+00 |      14 |  0.006095993 |        ideal
   17 |     200 |  0.00000E+00 |  0.00000E+00 |      17 |  0.002510398 |        ideal
   18 |     210 |  0.00000E+00 |  0.00000E+00 |      20 |  0.003652660 |            f
   19 |     220 |  0.00000E+00 |  0.00000E+00 |      20 |  0.010131820 |        nadir
   20 |     230 |  0.00000E+00 |  0.00000E+00 |      21 |  0.005676014 |            f
   21 |     240 |  0.00000E+00 |  0.00000E+00 |      25 |  0.010464402 |            f
   22 |     250 |  0.00000E+00 |  0.00000E+00 |      25 |  0.000547515 |            f
   23 |     260 |  0.00000E+00 |  0.00000E+00 |      28 |  0.001050255 |            f
   24 |     270 |  0.00000E+00 |  0.00000E+00 |      33 |  0.003841298 |            f
   25 |     280 |  0.00000E+00 |  0.00000E+00 |      37 |  0.006664377 |        nadir
   26 |     290 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000963164 |            f
   27 |     300 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000678243 |            f
   28 |     310 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000815766 |            f
   29 |     320 |  0.00000E+00 |  0.00000E+00 |      40 |  0.001500814 |            f
   30 |     330 |  0.00000E+00 |  0.00000E+00 |      40 |  0.014706442 |        nadir
   31 |     340 |  0.00000E+00 |  0.00000E+00 |      40 |  0.003554320 |        ideal
   32 |     350 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000624123 |            f
   33 |     360 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000203925 |            f
   34 |     370 |  0.00000E+00 |  0.00000E+00 |      40 |  0.001048509 |            f
   35 |     380 |  0.00000E+00 |  0.00000E+00 |      40 |  0.001121103 |            f
   36 |     390 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000664461 |            f
   37 |     400 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000761066 |            f
   38 |     410 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000521906 |            f
   39 |     420 |  0.00000E+00 |  0.00000E+00 |      40 |  0.004652095 |        nadir
   40 |     430 |  0.00000E+00 |  0.00000E+00 |      40 |  0.000287847 |            f

2.4 可视化

  可视化以更好地展示求解结果:

import matplotlib.pyplot as plt
problem = MyProblem()
xl, xu = problem.bounds()
plt.figure(figsize=(7, 5))
plt.scatter(X[:, 0], X[:, 1], s=30, facecolors='none', edgecolors='r')
plt.xlim(xl[0], xu[0])
plt.ylim(xl[1], xu[1])
plt.title("Search space")
plt.show()

plt.figure(figsize=(7, 5))
plt.scatter(F[:, 0], F[:, 1], s=30, facecolors='none', edgecolors='blue')
plt.title("Objective Space")
plt.show()

1234567891011121314

  输出如下:

3 完整代码

import numpy as np
from pymoo.core.problem import ElementwiseProblem
from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.factory import get_sampling, get_crossover, get_mutation
from pymoo.factory import get_termination
from pymoo.optimize import minimize


class MyProblem(ElementwiseProblem):

    def __init__(self):
        super(MyProblem, self).__init__(
            n_var=2,  # 变量数量
            n_obj=2,  # 优化目标的数量
            n_constr=2,  # 约束的数量
            xl=np.array([-2, -2]),  # 变量的下限
            xu=np.array([2, 2])  # 变量的上限
        )

    def _evaluate(self, x, out, *args, **kwargs):
        """
        :param x:       1D向量,长度为n_var
        :param out:
        :param args:
        :param kwargs:
        :return:
        """
        # 目标函数
        f1 = 100 * (x[0]**2 + x[1]**2)
        f2 = (x[0] - 1)**2 + x[1]**2

        # 约束条件
        g1 = 2 * (x[0] - 0.1) * (x[0] - 0.9) / 0.18
        g2 = - 20 * (x[0] - 0.4) * (x[0] - 0.6) / 4.8

        # 存放目标值和约束值
        out["F"] = [f1, f2]
        out["G"] = [g1, g2]


algorithm = NSGA2(
    pop_size=40,  # 种群数量
    n_offsprings=10,  # 后代数量
    sampling=get_sampling("real_random"),  # 采样的方式,如随机采样
    crossover=get_crossover("real_sbx", prob=0.9, eta=15),  # 交配方式,如模拟二进制
    mutation=get_mutation("real_pm", eta=20),  # 变异方式,多项式变异
    eliminate_duplicates=True   # 确保后代的目标是不同
)

termination = get_termination("n_gen", 40)

res = minimize(
    problem=MyProblem(),
    algorithm=algorithm,
    termination=termination,
    seed=1,
    save_history=True,
    verbose=True    # 输出控制变量
)

# 存储可行解
X = res.X
# 存储目标值
F = res.F

import matplotlib.pyplot as plt
problem = MyProblem()
xl, xu = problem.bounds()
plt.figure(figsize=(7, 5))
plt.scatter(X[:, 0], X[:, 1], s=30, facecolors='none', edgecolors='r')
plt.xlim(xl[0], xu[0])
plt.ylim(xl[1], xu[1])
plt.title("Search space")
plt.show()

plt.figure(figsize=(7, 5))
plt.scatter(F[:, 0], F[:, 1], s=30, facecolors='none', edgecolors='blue')
plt.title("Objective Space")
plt.show()