scipy.optimize.

differential_evolution#

scipy.optimize.differential_evolution(func, bounds, args=(), strategy='best1bin', maxiter=1000, popsize=15, tol=0.01, mutation=(0.5, 1), recombination=0.7, seed=None, callback=None, disp=False, polish=True, init='latinhypercube', atol=0, updating='immediate', workers=1, constraints=(), x0=None, *, integrality=None, vectorized=False)[source]#

查找多元函数的全局最小值。

微分进化法 [1]本质上是随机的。它不使用梯度法查找最小值,并且可以在候选空间的较大区域中进行搜索,但通常需要比传统基于梯度的技术更多的函数评估。

该算法归功于 Storn 和 Price [2]

参数:
func可调用

要最小化的目标函数。必须采用 f(x, *args) 形式,其中 x 是 1-D 数组形式的参数,args 是完全指定函数所需的所有附加固定参数的元组。形参的数量 N 等于 len(x)

bounds序列或 Bounds

变量的界限。有两种方法来指定界限

  1. Bounds 类的实例。

  2. 对于 x 中的每个元素,(min, max) 对,定义指定func 的优化参数的有限下界和上界。

边界总数用于确定参数数量 N。如果存在边界相等的参数,则自由参数总数为 N - N_equal

args元组,可选

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

策略(str、可调用函数),可选

要使用的差分进化策略。应为以下之一

  • “best1bin”

  • “best1exp”

  • “rand1bin”

  • “rand1exp”

  • “rand2bin”

  • “rand2exp”

  • “randtobest1bin”

  • “randtobest1exp”

  • “currenttobest1bin”

  • “currenttobest1exp”

  • “best2exp”

  • “best2bin”

默认值为“best1bin”。“注释”中概述了可能实施的策略。或者,可以通过提供用于构造试验向量的可调用函数来自定义差分进化策略。可调用函数必须具有 strategy(candidate: int, population: np.ndarray, rng=None) 的格式,其中 candidate 是一个整数,指定正在演化的种群中哪一个条目,population 是形状为 (S, N) 的数组,其中包含所有种群成员(其中 S 是总种群大小),rng 是求解器中使用的随机数生成器。 candidate 将在范围 [0, S) 中。 strategy 必须返回形状为 (N,) 的试验向量。该试验向量的适应度将与 population[candidate] 的适应度进行比较。

1.12.0 版本中已更改: 通过可调用函数自定义进化策略。

maxiterint,可选

整个种群演化的最大世代数。没有优化项的最大函数评估次数为: (maxiter + 1) * popsize * (N - N_equal)

popsizeint,可选

用来设置总体数量的乘数。总体有 群落大小 * (N - N_对等) 个体。如果通过 init 关键字提供了初始种群,则会被覆盖此关键字。当使用 init='sobol' 时,群体规模被计算为 群落大小 * (N - N_对等) 之后的下一个 2 的幂。

tol浮点数,可选

收敛的相对容差,求解在以下情况停止 np.std(pop) <= atol + tol * np.abs(np.mean(pop_energy)),其中 atoltol 分别为绝对容差和相对容差。

mutation浮点数或元组(浮点数, 浮点数),可选

变异常数。文献中也称其为差分权重,记为 F。若指定为浮点数,则应在 [0, 2] 范围内。若指定为元组 (min, max),则采用抖动。抖动随机地逐代改变变异常数。该代的变异常数取自 U[min, max)。抖动可以显著加快收敛速度。增加变异常数会增加搜索半径,但会减慢收敛速度。

recombination浮点数,可选

重组常数,应在 [0, 1] 范围内。文献中也称其为交叉概率,记为 CR。增加此值可以让更多的突变体进入下一代,但有损于群体稳定性。

seed{None、int、numpy.random.Generatornumpy.random.RandomState},可选

如果 seed 为 None(或 np.random),将使用 numpy.random.RandomState 单例。如果 seed 是一个整数,将使用一个新的 RandomState 实例,并用 seed 为种子。如果 seed 已是一个 GeneratorRandomState 实例,则使用该实例。指定 seed 以进行可重复最小化。

disp布尔值,可选

打印每次迭代中经评估的 func

callback可调用对象,可选

每次迭代后调用的一个可调用对象。具有签名

callback(intermediate_result: OptimizeResult)

其中 intermediate_result 是一个关键字参数,包含一个 OptimizeResult,具有属性 xfun(迄今为止找到的最佳解和目标函数)。请注意,参数名称必须是 intermediate_result,以便将 OptimizeResult 传递给回调函数。

回调还支持如下类型的签名

callback(x, convergence: float=val)

val 表示人口收敛的分数。当 val 大于 1.0 时,该函数将终止。

使用内省来确定调用哪个签名。

如果回调引发 StopIteration 或返回 True,则全局最小化将停止;仍会执行任何抛光操作。

在版本 1.12.0 中更改:callback 接受 intermediate_result 关键字。

polish布尔值,可选

如果为 True(默认),则在最后使用带有 L-BFGS-B 方法的 scipy.optimize.minimize 来优化种群中最优的成员,这可以略微改善优化。如果正在研究一个约束问题,则相反,将使用 trust-constr 方法。对于具有许多约束的大型问题,由于雅各比计算,优化可能需要很长时间。

initstr 或类数组,可选

指定将执行哪种类型的种群初始化。应为下列选项之一

  • ‘latinhypercube’

  • ‘sobol’

  • ‘halton’

  • ‘random’

  • 指定初始种群的数组。该数组应具有形状 (S, N),其中 S 是种群总数,N 是参数数。在使用之前,init 会被剪辑至 bounds

默认值为 ‘latinhypercube’。拉丁超立方体采样尝试最大化可用参数空间的覆盖率。

‘sobol’ 和 ‘halton’ 属于更高级的替代,以更大的程度最大化参数空间。‘sobol’ 将强制执行初始种群大小,该大小被计算为 popsize * (N - N_equal) 之后的下一个 2 的次方数。‘halton’ 没有要求,但效率稍低。参见 scipy.stats.qmc 了解更多详细信息。

‘random’ 随机初始化种群,其缺点在于可能发生聚类,导致无法覆盖整个参数空间。使用数组来指定种群可用于实现某种目的,例如在已知存在解的位置创建一组紧密的初始猜测,从而减少收敛时间。

atolfloat,可选

收敛绝对公差,求解在 np.std(pop) <= atol + tol * np.abs(np.mean(population_energies)) 时停止,其中 atoltol 分别是绝对公差和相对公差。

updating{'immediate', 'deferred'},可选

如果为 'immediate',最佳解向量会在单个世代中持续更新 [4]。这可以带来更快的收敛,因为试验向量可以利用最佳解的连续改进。如果为 'deferred',最佳解向量每代更新一次。只有 'deferred' 与并行化或向量化兼容,而且 workersvectorized 关键字会覆盖此选项。

在版本 1.2.0 中添加。

workersint 或类似 map 的可调用对象,可选

如果 workers 是一个整数,群体会被细分为 workers 个部分,并行评估(使用 multiprocessing.Pool)。提供 -1 以使用所有可用的 CPU 核心。此外,提供类似 map 的可调用对象,如 multiprocessing.Pool.map,用于并行评估群体。此评估执行为 workers(func, iterable)。如果 workers != 1,此选项会将 updating 关键字覆盖为 updating='deferred'。如果 workers != 1,此选项会覆盖 vectorized 关键字。需要 func 可被序列化。

在版本 1.2.0 中添加。

constraints{NonLinearConstraint, LinearConstraint, Bounds}

除了 bounds kwd 应用的约束之外,对求解器的约束。使用 Lampinen 的方法 [5]

在版本 1.4.0 中添加。

x0None 或类似数组的对象,可选

提供对最小化的初始估计。一旦群体被初始化,此向量会替换第一个(最佳)成员。即使为 init 给出初始群体,也会完成此替换。 x0.shape == (N,)

在版本 1.7.0 中添加。

integrality1-D 数组,可选

对于每一个决策变量,一个布尔值用于指示决策变量是否为受限的整数。数组广播为 (N,)。如果任何决策变量受到约束,则在优化期间不会对这些变量进行更改。仅使用在上下限之间的整数。如果没有在上下限之间的整数,则会引发 ValueError

添加到版本 1.9.0 中。

vectorized布尔值,可选

如果 vectorized True,则会将 x 数组发送给 func,其 x.shape == (N, S),并预期返回形状为 (S,) 的数组,其中 S 是要计算的解决方案向量的个数。如果应用了约束,则用于构建 Constraint 对象的每一个函数都应接受形状为 x.shape == (N, S)x 数组,并返回形状为 (M, S) 的数组,其中 M 是约束分量的个数。此选项是 workers 提供并行化的替代方案,通过减少多次函数调用的解释器开销,在优化速度方面有所帮助。如果 workers != 1,则此关键字会被忽略。此选项将会把 updating 关键字覆盖为 updating='deferred'。请参见附录以进一步讨论何时使用 'vectorized' 和何时使用 'workers'

添加到版本 1.9.0 中。

返回:
resOptimizeResult

OptimizeResult对象表示的优化结果。重要的属性有:x 解决方案数组,success 布尔标志,指示优化器是否成功退出,message 描述终止原因,population 群体中的解决方案向量,population_energies population中每个条目的目标函数值。有关其他属性的说明,请参阅OptimizeResult。如果采用了polish,并且通过抛光获得了较低最小值,则优化结果还包含jac属性。如果最终解决方案不满足应用的约束,success 将是False

注释

差分进化是一种基于随机群体的方法,适用于全局优化问题。在每次通过群体时,算法都会通过与其他候选解决方案混合来改变每个候选解决方案,以创建一个试验候选解决方案。有几种创建试验候选者的策略[3],这些策略有的适用于一些问题,有的适用于其他问题。对于许多系统,“best1bin”策略是一个良好的起点。在此策略中,随机选择种群中的两个成员。利用它们之间的差异改变迄今为止的最佳成员(“best1bin”中的“best”)\(x_0\)

\[b' = x_0 + mutation * (x_{r_0} - x_{r_1})\]

然后构造一个试验矢量。从随机选择的首个参数开始,使用来自 b' 或原始候选者的参数按(取模)顺序填写试验。使用 b' 或原始候选者的选择由二项分布(“best1bin” 中的“bin” - 生成为 [0, 1) 中的随机数)确定。如果此数字小于重组常量,那么参数会从 b' 加载,否则会从原始候选者加载。最后一个参数始终从 b' 加载。创建试验候选者后,会对其适应度进行评估。如果试验比原始候选者更佳,则它会取代原始候选者。如果它也比最佳整体候选者更佳,则它还会取代该候选者。

其他可用策略在 Qiang and Mitchell (2014) [3] 中概述。

\[ \begin{align}\begin{aligned}rand1* : b' = x_{r_0} + mutation*(x_{r_1} - x_{r_2})\\rand2* : b' = x_{r_0} + mutation*(x_{r_1} + x_{r_2} - x_{r_3} - x_{r_4})\\best1* : b' = x_0 + mutation*(x_{r_0} - x_{r_1})\\best2* : b' = x_0 + mutation*(x_{r_0} + x_{r_1} - x_{r_2} - x_{r_3})\\currenttobest1* : b' = x_i + mutation*(x_0 - x_i + x_{r_0} - x_{r_1})\\randtobest1* : b' = x_{r_0} + mutation*(x_0 - x_{r_0} + x_{r_1} - x_{r_2})\end{aligned}\end{align} \]

其中整数 \(r_0, r_1, r_2, r_3, r_4\) 从区间 [0, NP) 中随机选择,其中 NP 是总种群大小,而原始候选者的索引为 i。用户可以通过向 strategy提供可调用项来完全自定义试验候选者的生成。

若要提高找到全局最小值的可能性,请使用较高的 popsize 值,较高的 mutation 和(抖动),但较低的 recombination 值。这会产生扩大搜索半径但减慢收敛速度的效果。

默认情况下,在单个迭代中,最佳解向量会持续更新 (更新='立即')。这是原始差分进化算法的修改 [4],可以通过使试验向量直接受益于改进的解,实现更快的收敛。要使用原始 Storn 和 Price 行为(每个迭代更新一次最佳解),请设置 更新='延迟''延迟' 方法与并行化和矢量化 ('工作人员''矢量化' 关键字) 都兼容。通过更有效率地使用计算机资源,它们可以提高最小化速度。'工作人员' 将计算任务分配给多个处理器。默认情况下,将使用 Python multiprocessing 模块,但也可以使用其他方法,例如群集上使用的消息传递接口 (MPI) [6] [7]。这些方法的开销(创建新进程等)可能很大,这意味着计算速度不一定随着所用处理器数量的增加而增加。并行化最适合计算开销大的目标函数。如果目标函数开销较小,则 '矢量化' 可以通过每个迭代仅调用一次目标函数(而不是对所有种群成员多次调用目标函数)来提供帮助;减少了解释器开销。

v0.15.0 中新增。

参考

[1]

差分进化,维基百科,http://en.wikipedia.org/wiki/Differential_evolution

[2]

Storn,R 和 Price,K,《差分进化 - 连续空间全局优化的一种简单而有效的启发式算法》,《全球优化杂志》,1997 年,11,341 - 359。

[3] (1,2)

祁ang,J.,Mitchell,C.,面向全局优化的统一差分进化算法,2014,https://www.osti.gov/servlets/purl/1163659

[4] (1,2)

Wormington,M.,Panaccione,C.,Matney,K. M.,Bowen,D. K. - 使用遗传算法对 X 射线散射数据中的结构进行表征,《英国皇家学会哲学评论 A》,1999,357,2827-2848

[5]

Lampinen,J.,差分进化算法的约束处理方法。《2002 年进化计算大会论文集。CEC'02 (目录号 02TH8600)》。卷 2。IEEE,2002。

示例

让我们考虑极小化 Rosenbrock 函数的问题。此函数在 rosen 中以 scipy.optimize 内实现。

>>> import numpy as np
>>> from scipy.optimize import rosen, differential_evolution
>>> bounds = [(0,2), (0, 2), (0, 2), (0, 2), (0, 2)]
>>> result = differential_evolution(rosen, bounds)
>>> result.x, result.fun
(array([1., 1., 1., 1., 1.]), 1.9216496320061384e-19)

现在重复,但要进行并行化。

>>> result = differential_evolution(rosen, bounds, updating='deferred',
...                                 workers=2)
>>> result.x, result.fun
(array([1., 1., 1., 1., 1.]), 1.9216496320061384e-19)

让我们进行受限极小化。

>>> from scipy.optimize import LinearConstraint, Bounds

我们添加 x[0]x[1] 之和必须小于等于 1.9 的限制。这是一个线性限制,可以写为 A @ x <= 1.9,其中 A = array([[1, 1]])。这可以编码为 LinearConstraint 实例

>>> lc = LinearConstraint([[1, 1]], -np.inf, 1.9)

使用 Bounds 对象指定限制。

>>> bounds = Bounds([0., 0.], [2., 2.])
>>> result = differential_evolution(rosen, bounds, constraints=lc,
...                                 seed=1)
>>> result.x, result.fun
(array([0.96632622, 0.93367155]), 0.0011352416852625719)

接下来找到 Ackley 函数的极小值 (https://en.wikipedia.org/wiki/Test_functions_for_optimization)。

>>> def ackley(x):
...     arg1 = -0.2 * np.sqrt(0.5 * (x[0] ** 2 + x[1] ** 2))
...     arg2 = 0.5 * (np.cos(2. * np.pi * x[0]) + np.cos(2. * np.pi * x[1]))
...     return -20. * np.exp(arg1) - np.exp(arg2) + 20. + np.e
>>> bounds = [(-5, 5), (-5, 5)]
>>> result = differential_evolution(ackley, bounds, seed=1)
>>> result.x, result.fun
(array([0., 0.]), 4.440892098500626e-16)

Ackley 函数以矢量化方式编写,因此可以使用 'vectorized' 关键字。注意,减少的函数评估次数。

>>> result = differential_evolution(
...     ackley, bounds, vectorized=True, updating='deferred', seed=1
... )
>>> result.x, result.fun
(array([0., 0.]), 4.440892098500626e-16)

下面的自定义策略函数模仿‘best1bin’

>>> def custom_strategy_fn(candidate, population, rng=None):
...     parameter_count = population.shape(-1)
...     mutation, recombination = 0.7, 0.9
...     trial = np.copy(population[candidate])
...     fill_point = rng.choice(parameter_count)
...
...     pool = np.arange(len(population))
...     rng.shuffle(pool)
...
...     # two unique random numbers that aren't the same, and
...     # aren't equal to candidate.
...     idxs = []
...     while len(idxs) < 2 and len(pool) > 0:
...         idx = pool[0]
...         pool = pool[1:]
...         if idx != candidate:
...             idxs.append(idx)
...
...     r0, r1 = idxs[:2]
...
...     bprime = (population[0] + mutation *
...               (population[r0] - population[r1]))
...
...     crossovers = rng.uniform(size=parameter_count)
...     crossovers = crossovers < recombination
...     crossovers[fill_point] = True
...     trial = np.where(crossovers, bprime, trial)
...     return trial