scipy.optimize.

shgo#

scipy.optimize.shgo(func, bounds, args=(), constraints=None, n=100, iters=1, callback=None, minimizer_kwargs=None, options=None, sampling_method='simplicial', *, workers=1)[source]#

使用 SHG 优化查找函数的全局最小值。

SHGO 代表“单纯形同源全局优化”。

参数::
func可调用对象

要最小化的目标函数。必须采用 f(x, *args) 的形式,其中 x 是参数,形式为一维数组,args 是元组,包含完全指定函数所需的其他固定参数。

bounds序列或 Bounds

变量的边界。有两种方法可以指定边界

  1. Bounds 类实例。

  2. (min, max) 对的序列,用于 x 中的每个元素。

args元组,可选

完全指定目标函数所需的任何其他固定参数。

constraints{Constraint, dict} 或 {Constraint, dict} 的列表,可选

约束定义。仅适用于 COBYLA、COBYQA、SLSQP 和 trust-constr。有关指定约束的更多详细信息,请参见教程 [5]

注意

目前,只有 COBYLA、COBYQA、SLSQP 和 trust-constr 局部最小化方法支持约束参数。如果局部优化问题中使用的 constraints 序列未在 minimizer_kwargs 中定义,并且使用了受限方法,那么将使用全局 constraints。(在 minimizer_kwargs 中定义 constraints 序列意味着不会添加 constraints,因此如果需要添加等式约束等,则需要将 constraints 中的不等式函数也添加到 minimizer_kwargs 中)。COBYLA 仅支持不等式约束。

版本 1.11.0 中变更: constraints 接受 NonlinearConstraintLinearConstraint

nint,可选

用于构建单纯形复形的采样点数。对于默认的 simplicial 采样方法,将生成 2**dim + 1 个采样点,而不是默认的 n=100。对于所有其他指定的 n 值,将生成 n 个采样点。对于 sobolhalton 和其他任意 sampling_methods,将生成 n=100 或其他指定的采样点数。

itersint,可选

用于构建单纯形复形的迭代次数。默认值为 1。

callback可调用对象,可选

在每次迭代后调用,如 callback(xk),其中 xk 是当前参数向量。

minimizer_kwargsdict,可选

要传递给最小化器 scipy.optimize.minimize 的额外关键字参数。一些重要的选项可能是

  • methodstr

    最小化方法。如果未给出,则选择为 BFGS、L-BFGS-B、SLSQP 之一,具体取决于问题是否有约束或边界。

  • args元组

    传递给目标函数 (func) 及其导数(雅可比矩阵、海森矩阵)的额外参数。

  • optionsdict,可选

    请注意,默认情况下,容差被指定为 {ftol: 1e-12}

optionsdict,可选

求解器选项的字典。全局例程中指定的大多数选项也会传递给 scipy.optimize.minimize 例程。也会传递给局部例程的选项标记为“(L)”。

停止标准,如果满足任何指定的标准,算法将终止。但是,默认算法不需要指定任何标准

  • maxfevint (L)

    可行域中函数评估的最大次数。(请注意,只有支持此选项的方法才能在精确指定的准确值处终止例程。否则,标准只会终止全局迭代期间)

  • f_min

    指定最小目标函数值,如果已知。

  • f_tolfloat

    停止标准中 f 值的精度目标。请注意,如果全局例程中的采样点在该容差范围内,全局例程也会终止。

  • maxiterint

    要执行的最大迭代次数。

  • maxevint

    要执行的最大采样评估次数(包括在不可行点进行搜索)。

  • maxtimefloat

    允许的最大处理运行时间

  • minhgrdint

    最小同源群秩微分。同源群的目标函数在每次迭代期间计算(近似)。该群的秩与目标函数中局部凸子域的数量有一一对应关系(在每个子域包含唯一的全局最小值后进行足够的采样点)。如果 maxhgrd 指定的迭代次数中,迭代之间的 hgr 差异为 0,则算法将终止。

目标函数知识

  • symmetrylist 或 bool

    指定目标函数是否包含对称变量。在完全对称的情况下,搜索空间(因此性能)最多会减少 O(n!) 倍。如果指定了 True,则所有变量都将设置为与第一个变量对称。默认设置为 False。

    例如 f(x) = (x_1 + x_2 + x_3) + (x_4)**2 + (x_5)**2 + (x_6)**2

    在这个等式中,x_2 和 x_3 与 x_1 对称,而 x_5 和 x_6 与 x_4 对称,这可以指定为求解器为

    symmetry = [0, # 变量 1

    0, # 与变量 1 对称 0, # 与变量 1 对称 3, # 变量 4 3, # 与变量 4 对称 3, # 与变量 4 对称 ]

  • jacbool 或可调用对象,可选

    目标函数的雅可比矩阵(梯度)。仅适用于 CG、BFGS、Newton-CG、L-BFGS-B、TNC、SLSQP、dogleg、trust-ncg。如果 jac 是布尔值并且为 True,则假设 fun 返回梯度以及目标函数。如果为 False,则将以数值方式估计梯度。 jac 也可以是可调用对象,返回目标函数的梯度。在这种情况下,它必须接受与 fun 相同的参数。(自动传递给 scipy.optimize.minimize

  • hess、hessp可调用对象,可选

    目标函数的海森矩阵(二阶导数矩阵)或目标函数乘以任意向量 p 的海森矩阵。仅适用于 Newton-CG、dogleg、trust-ncg。只需要给出 hessphess 之一。如果提供了 hess,则会忽略 hessp。如果既没有提供 hess 也没有提供 hessp,则将使用对 jac 的有限差分来近似海森矩阵乘积。 hessp 必须计算海森矩阵乘以任意向量。(自动传递给 scipy.optimize.minimize

算法设置

  • minimize_every_iterbool

    如果为 True,则每次迭代都会将有希望的全局采样点传递给局部最小化例程。如果为 True,则只会运行最终的最小化器池。默认为 True。

  • local_iterint

    每次迭代只评估几个最佳最小化器池候选。如果为 False,则所有潜在点都会传递给局部最小化例程。

  • infty_constraintsbool

    如果为 True,则会保存所有生成的超出可行域的采样点,并为其赋予 inf 的目标函数值。如果为 False,则会丢弃这些点。使用此功能可以提高在找到全局最小值之前函数评估的性能,指定 False 会使用更少的内存,但会略微降低性能。默认为 True。

反馈

  • dispbool (L)

    设置为 True 以打印收敛消息。

sampling_methodstr 或 function,可选

当前内置的采样方法选项包括 haltonsobolsimplicial。默认的 simplicial 提供了在有限时间内收敛到全局最小值的理论保证。 haltonsobol 方法在采样点生成速度方面更快,但代价是失去了收敛保证。它更适合大多数“更容易”的问题,这些问题的收敛速度相对较快。用户定义的采样函数必须接受两个参数,每个调用分别为维度为 dimn 个采样点,并输出一个形状为 n x dim 的采样点数组。

workersint 或 map-like 可调用对象,可选

并行采样并运行本地串行最小化。提供 -1 以使用所有可用的 CPU 内核,或提供一个整数以使用相应的进程数量(使用 multiprocessing.Pool)。

或者,提供一个 map-like 可调用对象,例如 multiprocessing.Pool.map 以进行并行评估。此评估将以 workers(func, iterable) 的形式进行。要求 func 可被 pickle 化。

在版本 1.11.0 中添加。

返回值:
resOptimizeResult

将优化结果表示为 OptimizeResult 对象。重要的属性包括:x 对应于全局最小值的解数组,fun 全局解处的函数输出,xl 本地最小值解的有序列表,funl 对应于本地解的函数输出,success 一个布尔标志,指示优化器是否成功退出,message 描述终止原因,nfev 目标函数评估的总次数(包括采样调用),nlfev 来自所有局部搜索优化的目标函数评估的总次数,nit 全局例程执行的迭代次数。

注释

使用单纯形同源全局优化 [1] 进行全局优化。适用于解决通用的 NLP 和黑盒优化问题以达到全局最优(低维问题)。

一般来说,优化问题具有以下形式

minimize f(x) subject to

g_i(x) >= 0,  i = 1,...,m
h_j(x)  = 0,  j = 1,...,p

其中 x 是一个或多个变量的向量。 f(x) 是目标函数 R^n -> Rg_i(x) 是不等式约束,h_j(x) 是等式约束。

可选地,可以使用 bounds 参数指定 x 中每个元素的下界和上界。

虽然 SHGO 的大多数理论优势仅在 f(x) 是 Lipschitz 光滑函数时得到证明,但该算法也被证明收敛到全局最优,即使 f(x) 是非连续、非凸和非光滑的,只要使用默认采样方法 [1]

可以使用 minimizer_kwargs 参数指定局部搜索方法,该参数会传递给 scipy.optimize.minimize。默认情况下,使用 SLSQP 方法。一般来说,如果为问题定义了不等式约束,建议使用 SLSQPCOBYLACOBYQA 局部最小化,因为其他方法不使用约束。

使用 scipy.stats.qmc 生成 haltonsobol 方法的点。可以使用任何其他 QMC 方法。

参考文献

[1] (1,2)

Endres, SC, Sandrock, C, Focke, WW (2018) “A simplicial homology algorithm for lipschitz optimisation”, Journal of Global Optimization.

[2]

Joe, SW and Kuo, FY (2008) “Constructing Sobol’ sequences with better two-dimensional projections”, SIAM J. Sci. Comput. 30, 2635-2654.

[3] (1,2)

Hock, W and Schittkowski, K (1981) “Test examples for nonlinear programming codes”, Lecture Notes in Economics and Mathematical Systems, 187. Springer-Verlag, New York. http://www.ai7.uni-bayreuth.de/test_problem_coll.pdf

[4]

Wales, DJ (2015) “Perspective: Insight into reaction coordinates and dynamics from the potential energy landscape”, Journal of Chemical Physics, 142(13), 2015.

示例

首先考虑最小化 Rosenbrock 函数 rosen 的问题。

>>> from scipy.optimize import rosen, shgo
>>> bounds = [(0,2), (0, 2), (0, 2), (0, 2), (0, 2)]
>>> result = shgo(rosen, bounds)
>>> result.x, result.fun
(array([1., 1., 1., 1., 1.]), 2.920392374190081e-18)

请注意,边界决定了目标函数的维数,因此是必需的输入,但是您可以使用 Nonenp.inf 之类的对象指定空边界,这些对象将被转换为较大的浮点数。

>>> bounds = [(None, None), ]*4
>>> result = shgo(rosen, bounds)
>>> result.x
array([0.99999851, 0.99999704, 0.99999411, 0.9999882 ])

接下来,我们考虑 Eggholder 函数,这是一个具有多个局部最小值和一个全局最小值的问题。我们将演示参数的使用以及 shgo 的功能。(https://en.wikipedia.org/wiki/Test_functions_for_optimization

>>> import numpy as np
>>> def eggholder(x):
...     return (-(x[1] + 47.0)
...             * np.sin(np.sqrt(abs(x[0]/2.0 + (x[1] + 47.0))))
...             - x[0] * np.sin(np.sqrt(abs(x[0] - (x[1] + 47.0))))
...             )
...
>>> bounds = [(-512, 512), (-512, 512)]

shgo 内置了低差异采样序列。首先,我们将输入 64 个 Sobol’ 序列的初始采样点

>>> result = shgo(eggholder, bounds, n=64, sampling_method='sobol')
>>> result.x, result.fun
(array([512.        , 404.23180824]), -959.6406627208397)

shgo 还返回了找到的任何其他局部最小值,可以使用以下方法调用这些值

>>> result.xl
array([[ 512.        ,  404.23180824],
       [ 283.0759062 , -487.12565635],
       [-294.66820039, -462.01964031],
       [-105.87688911,  423.15323845],
       [-242.97926   ,  274.38030925],
       [-506.25823477,    6.3131022 ],
       [-408.71980731, -156.10116949],
       [ 150.23207937,  301.31376595],
       [  91.00920901, -391.283763  ],
       [ 202.89662724, -269.38043241],
       [ 361.66623976, -106.96493868],
       [-219.40612786, -244.06020508]])
>>> result.funl
array([-959.64066272, -718.16745962, -704.80659592, -565.99778097,
       -559.78685655, -557.36868733, -507.87385942, -493.9605115 ,
       -426.48799655, -421.15571437, -419.31194957, -410.98477763])

这些结果在存在许多全局最小值且需要其他全局最小值的价值或局部最小值可以提供对系统洞察力的应用中很有用(例如物理化学中的形态 [4])。

如果我们想要找到更多数量的局部最小值,可以增加采样点的数量或迭代次数。我们将采样点的数量增加到 64,并将迭代次数从默认的 1 增加到 3。使用 simplicial,这将给我们 64 x 3 = 192 个初始采样点。

>>> result_2 = shgo(eggholder,
...                 bounds, n=64, iters=3, sampling_method='sobol')
>>> len(result.xl), len(result_2.xl)
(12, 23)

请注意,例如,n=192, iters=1n=64, iters=3 之间的区别。在第一种情况下,最小化池中包含的有希望的点仅处理一次。在后一种情况下,每 64 个采样点处理一次,总共处理 3 次。

为了演示如何解决具有非线性约束的问题,请考虑以下来自 Hock 和 Schittkowski 问题 73(牛饲料)的示例 [3]

minimize: f = 24.55 * x_1 + 26.75 * x_2 + 39 * x_3 + 40.50 * x_4

subject to: 2.3 * x_1 + 5.6 * x_2 + 11.1 * x_3 + 1.3 * x_4 - 5    >= 0,

            12 * x_1 + 11.9 * x_2 + 41.8 * x_3 + 52.1 * x_4 - 21
                -1.645 * sqrt(0.28 * x_1**2 + 0.19 * x_2**2 +
                              20.5 * x_3**2 + 0.62 * x_4**2)      >= 0,

            x_1 + x_2 + x_3 + x_4 - 1                             == 0,

            1 >= x_i >= 0 for all i

[3] 中给出的近似答案是

f([0.6355216, -0.12e-11, 0.3127019, 0.05177655]) = 29.894378
>>> def f(x):  # (cattle-feed)
...     return 24.55*x[0] + 26.75*x[1] + 39*x[2] + 40.50*x[3]
...
>>> def g1(x):
...     return 2.3*x[0] + 5.6*x[1] + 11.1*x[2] + 1.3*x[3] - 5  # >=0
...
>>> def g2(x):
...     return (12*x[0] + 11.9*x[1] +41.8*x[2] + 52.1*x[3] - 21
...             - 1.645 * np.sqrt(0.28*x[0]**2 + 0.19*x[1]**2
...                             + 20.5*x[2]**2 + 0.62*x[3]**2)
...             ) # >=0
...
>>> def h1(x):
...     return x[0] + x[1] + x[2] + x[3] - 1  # == 0
...
>>> cons = ({'type': 'ineq', 'fun': g1},
...         {'type': 'ineq', 'fun': g2},
...         {'type': 'eq', 'fun': h1})
>>> bounds = [(0, 1.0),]*4
>>> res = shgo(f, bounds, n=150, constraints=cons)
>>> res
 message: Optimization terminated successfully.
 success: True
     fun: 29.894378159142136
    funl: [ 2.989e+01]
       x: [ 6.355e-01  1.137e-13  3.127e-01  5.178e-02] # may vary
      xl: [[ 6.355e-01  1.137e-13  3.127e-01  5.178e-02]] # may vary
     nit: 1
    nfev: 142 # may vary
   nlfev: 35 # may vary
   nljev: 5
   nlhev: 0
>>> g1(res.x), g2(res.x), h1(res.x)
(-5.062616992290714e-14, -2.9594104944408173e-12, 0.0)