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, rng=None, callback=None, disp=False, polish=True, init='latinhypercube', atol=0, updating='immediate', workers=1, constraints=(), x0=None, *, integrality=None, vectorized=False)[源代码]#
查找多元函数的全局最小值。
微分进化方法 [1] 本质上是随机的。它不使用梯度方法来寻找最小值,并且可以搜索候选空间的大片区域,但通常比传统的基于梯度的方法需要更多的函数评估。
该算法归功于 Storn 和 Price [2]。
- 参数:
- func可调用对象
要最小化的目标函数。 必须采用
f(x, *args)
的形式,其中x
是 1-D 数组形式的参数,args
是任何额外的固定参数的元组,这些参数用于完全指定该函数。参数的数量 N 等于len(x)
。- bounds序列或
Bounds
变量的边界。有两种指定边界的方法
Bounds
类的实例。用于
x
中每个元素的(min, max)
对,定义了 func 的优化参数的有限下界和上界。
边界的总数用于确定参数的数量 N。如果存在边界相等的参数,则自由参数的总数为
N - N_equal
。- argstuple, 可选
完全指定目标函数所需的任何其他固定参数。
- strategy{str, callable}, 可选
要使用的微分进化策略。应为以下之一
‘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, 可选
用于设置总种群大小的乘数。 种群有
popsize * (N - N_equal)
个个体。 如果通过 init 关键字提供了初始种群,则此关键字将被覆盖。当使用init='sobol'
时,种群大小计算为popsize * (N - N_equal)
之后的下一个 2 的幂。- tolfloat, 可选
收敛的相对公差,当
np.std(population_energies) <= atol + tol * np.abs(np.mean(population_energies))
时,求解停止,其中 atol 和 tol 分别是绝对公差和相对公差。- mutationfloat 或 tuple(float, float), 可选
突变常数。 在文献中,这也称为差分权重,用 \(F\) 表示。 如果指定为浮点数,则应在 [0, 2) 范围内。 如果指定为元组
(min, max)
,则采用抖动。 抖动在每一代的基础上随机更改突变常数。 该代的突变常数取自U[min, max)
。 抖动可以帮助显着加快收敛速度。 增加突变常数会增加搜索半径,但会减慢收敛速度。- recombinationfloat, 可选
重组常数,应在 [0, 1] 范围内。 在文献中,这也称为交叉概率,用 CR 表示。 增加此值允许更多的突变体进入下一代,但存在种群稳定性的风险。
- rng{None, int,
numpy.random.Generator
}, 可选 如果通过关键字传递了 rng,则将
numpy.random.Generator
以外的类型传递给numpy.random.default_rng
以实例化Generator
。 如果 rng 已经是Generator
实例,则使用提供的实例。 指定 rng 以实现可重复的函数行为。如果通过位置传递此参数或通过关键字传递了 seed,则应用参数 seed 的旧式行为
如果 seed 为 None (或
numpy.random
),则使用numpy.random.RandomState
单例。如果 seed 是一个整数,则会使用一个新的
RandomState
实例,并使用 seed 作为种子。如果 seed 已经是一个
Generator
或RandomState
实例,则会直接使用该实例。
在 1.15.0 版本中更改: 作为从使用
numpy.random.RandomState
过渡到使用numpy.random.Generator
的 SPEC-007 规范的一部分,此关键字从 seed 更改为 rng。在过渡期间,两个关键字都将继续有效,但一次只能指定一个。过渡期结束后,使用 seed 关键字的函数调用将发出警告。上面概述了 seed 和 rng 的行为,但在新代码中应仅使用 rng 关键字。- dispbool,可选
在每次迭代时打印评估的 func。
- callback可调用对象,可选
在每次迭代后调用的可调用对象。具有以下签名:
callback(intermediate_result: OptimizeResult)
其中
intermediate_result
是一个包含OptimizeResult
实例的关键字参数,该实例具有属性x
和fun
,分别表示目前找到的最佳解和目标函数值。请注意,参数的名称必须是intermediate_result
,才能将OptimizeResult
传递给回调函数。回调函数也支持类似以下的签名:
callback(x, convergence: float=val)
val
表示种群收敛的分数值。当val
大于1.0
时,函数将停止。使用内省来确定调用哪个签名。
如果回调函数引发
StopIteration
或返回True
,则全局最小化将停止;但仍会执行任何润色操作。在 1.12.0 版本中更改: 回调函数接受
intermediate_result
关键字。- polishbool,可选
如果为 True(默认值),则使用具有 L-BFGS-B 方法的
scipy.optimize.minimize
在最后对最佳种群成员进行润色,这可以稍微改善最小化效果。如果要研究约束问题,则会改为使用 trust-constr 方法。对于具有许多约束的大型问题,由于雅可比矩阵的计算,润色可能需要很长时间。在 1.15.0 版本中更改: 如果指定了 workers,则会将包装 func 的类似 map 的可调用对象提供给
minimize
,而不是直接使用 func。这允许调用者控制实际运行的调用方式和位置。- 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))
时停止求解,其中 atol 和 tol 分别是绝对容差和相对容差。- updating{‘immediate’, ‘deferred’},可选
如果为
'immediate'
,则在单个世代内不断更新最佳解向量 [4]。这可以加快收敛速度,因为试验向量可以利用最佳解的持续改进。如果为'deferred'
,则每个世代更新一次最佳解向量。只有'deferred'
与并行化或矢量化兼容,并且 workers 和 vectorized 关键字可以覆盖此选项。在 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 版本中添加。
- integrality一维数组,可选
对于每个决策变量,一个布尔值指示该决策变量是否被约束为整数值。该数组广播为
(N,)
。如果任何决策变量被约束为整数,则在润色期间不会更改它们。仅使用位于下限和上限之间的整数值。如果下限和上限之间没有整数值,则会引发 ValueError。在 1.9.0 版本中添加。
- vectorizedbool,可选
如果
vectorized 是 True
,则会向 func 发送一个 x 数组,其形状为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,并且通过 polishing 获得了更小的最小值,则 OptimizeResult 还包含jac
属性。如果最终解不满足应用的约束,则success
将为 False。
注释
差分进化是一种基于随机种群的方法,适用于全局优化问题。在每次遍历种群时,该算法都会通过与其他候选解混合来变异每个候选解,从而创建一个试验候选解。有几种创建试验候选解的策略 [3],这些策略更适合某些问题。“best1bin” 策略是许多系统的良好起点。在此策略中,随机选择种群的两个成员。它们的差异用于变异到目前为止的最佳成员(“best1bin” 中的“best”),\(x_0\),因此
\[b' = x_0 + F \cdot (x_{r_0} - x_{r_1})\]其中 \(F\) 是 mutation 参数。然后构建一个试验向量。从随机选择的第 i 个参数开始,试验向量会按顺序(以模数方式)填充来自
b'
或原始候选参数的参数。是否使用b'
或原始候选参数的选择是使用二项式分布(“best1bin” 中的 “bin”)进行的 - 生成 [0, 1) 中的随机数。如果此数字小于 recombination 常数,则从b'
加载参数,否则从原始候选参数加载。最终参数始终从b'
加载。构建试验候选参数后,将评估其适应度。如果试验候选参数优于原始候选参数,则它将取代其位置。如果它也优于最佳总体候选参数,则它也会取代该参数。其他可用策略在 Qiang 和 Mitchell (2014) 中概述 [3]。
rand1
: \(b' = x_{r_0} + F \cdot (x_{r_1} - x_{r_2})\)rand2
: \(b' = x_{r_0} + F \cdot (x_{r_1} + x_{r_2} - x_{r_3} - x_{r_4})\)best1
: \(b' = x_0 + F \cdot (x_{r_0} - x_{r_1})\)best2
: \(b' = x_0 + F \cdot (x_{r_0} + x_{r_1} - x_{r_2} - x_{r_3})\)currenttobest1
: \(b' = x_i + F \cdot (x_0 - x_i + x_{r_0} - x_{r_1})\)randtobest1
: \(b' = x_{r_0} + F \cdot (x_0 - x_{r_0} + x_{r_1} - x_{r_2})\)
其中整数 \(r_0, r_1, r_2, r_3, r_4\) 是从区间 [0, NP) 中随机选择的,其中 NP 是总种群大小,并且原始候选参数的索引为 i。用户可以通过向
strategy
提供可调用对象来完全自定义试验候选参数的生成。为了提高找到全局最小值的机会,请使用更高的 popsize 值,更高的 mutation 值(和抖动),但更低的 recombination 值。这具有扩大搜索半径的效果,但会减慢收敛速度。
默认情况下,最佳解向量在单个迭代中持续更新(
updating='immediate'
)。这是对原始差分进化算法的修改 [4],它可以加快收敛速度,因为试验向量可以立即从改进的解中受益。要使用原始的 Storn 和 Price 行为,每次迭代更新一次最佳解,请设置updating='deferred'
。'deferred'
方法与并行化和向量化('workers'
和'vectorized'
关键字)兼容。这些可以通过更有效地利用计算机资源来提高最小化速度。'workers'
将计算分布在多个处理器上。默认情况下,使用 Pythonmultiprocessing
模块,但也可以使用其他方法,例如在集群上使用的消息传递接口 (MPI) [6] [7]。这些方法(创建新进程等)的开销可能很大,这意味着计算速度不一定会随着使用的处理器数量而扩展。并行化最适合于计算量大的目标函数。如果目标函数的计算量较小,则'vectorized'
可能会通过每次迭代仅调用一次目标函数,而不是为所有种群成员调用多次目标函数来提供帮助;从而减少解释器开销。在 0.15.0 版本中添加。
参考文献
[1]差分进化,维基百科,http://en.wikipedia.org/wiki/Differential_evolution
[2]Storn, R 和 Price, K, Differential Evolution - a Simple and Efficient Heuristic for Global Optimization over Continuous Spaces, Journal of Global Optimization, 1997, 11, 341 - 359.
[3] (1,2)Qiang, J., Mitchell, C., A Unified Differential Evolution Algorithm for Global Optimization, 2014, https://www.osti.gov/servlets/purl/1163659
[4] (1,2)Wormington, M., Panaccione, C., Matney, K. M., Bowen, D. K., - Characterization of structures from X-ray scattering data using genetic algorithms, Phil. Trans. R. Soc. Lond. A, 1999, 357, 2827-2848
[5]Lampinen, J., A constraint handling approach for the differential evolution algorithm. Proceedings of the 2002 Congress on Evolutionary Computation. CEC’02 (Cat. No. 02TH8600). Vol. 2. IEEE, 2002.
示例
让我们考虑最小化 Rosenbrock 函数的问题。此函数在
rosen
inscipy.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, ... rng=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, rng=1) >>> result.x, result.fun (array([0., 0.]), 4.440892098500626e-16)
Ackley 函数以向量化的方式编写,因此可以使用
'vectorized'
关键字。请注意,函数评估的次数减少了。>>> result = differential_evolution( ... ackley, bounds, vectorized=True, updating='deferred', rng=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