scipy.optimize.minimize(method='trust-constr') 不会在 xtol 条件下终止

Posted 2023-03-12

【中文标题】scipy.optimize.minimize(method=\'trust-constr\') 不会在 xtol 条件下终止

【英文标题】：scipy.optimize.minimize(method=’trust-constr’) doesn't terminate on xtol conditionscipy.optimize.minimize(method='trust-constr') 不会在 xtol 条件下终止

【发布时间】：2019-11-28 03:27:59

【问题描述】：

我已经建立了一个线性等式约束的优化问题，如下所示

sol0 = minimize(objective, x0, args=mock_df, method='trust-constr',
                bounds=bnds, constraints=cons,
                options='maxiter': 250, 'verbose': 3)

objective 是一个加权求和函数，其系数/权重将被优化以使其最小化。因为我对系数和约束都有边界，所以我在scipy.optimize.minimize 中使用了trust-constr 方法。

最小化可行，但我不明白终止标准。根据trust-constr documentation，它应该终止于xtol

算法将在tr_radius < xtol 时终止，其中tr_radius 是算法中使用的信任区域的半径。默认为 1e-8。

但是，verbose 输出显示，终止确实是由 barrier_tol 参数触发的，如下面的清单所示

| niter |f evals|CG iter|  obj func   |tr radius |   opt    |  c viol  | penalty  |barrier param|CG stop|
|-------|-------|-------|-------------|----------|----------|----------|----------|-------------|-------|
C:\ProgramData\Anaconda3\lib\site-packages\scipy\optimize\_trustregion_constr\projections.py:182: UserWarning: Singular Jacobian matrix. Using SVD decomposition to perform the factorizations.
  warn('Singular Jacobian matrix. Using SVD decomposition to ' +
|   1   |  31   |   0   | -4.4450e+02 | 1.00e+00 | 7.61e+02 | 5.00e-01 | 1.00e+00 |  1.00e-01   |   0   |
C:\ProgramData\Anaconda3\lib\site-packages\scipy\optimize\_hessian_update_strategy.py:187: UserWarning: delta_grad == 0.0. Check if the approximated function is linear. If the function is linear better results can be obtained by defining the Hessian as zero instead of using quasi-Newton approximations.
  'approximations.', UserWarning)
|   2   |  62   |   1   | -2.2830e+03 | 6.99e+00 | 3.64e+02 | 7.28e-01 | 1.00e+00 |  1.00e-01   |   2   |
|   3   |  93   |   2   | -9.7651e+03 | 3.42e+01 | 5.52e+01 | 5.33e+00 | 1.00e+00 |  1.00e-01   |   2   |
|   4   |  124  |  26   | -4.9999e+03 | 3.42e+01 | 8.23e+01 | 9.29e-01 | 3.48e+16 |  1.00e-01   |   1   |
|   5   |  155  |  50   | -4.1486e+03 | 3.42e+01 | 5.04e+01 | 2.08e-01 | 3.48e+16 |  1.00e-01   |   1   |
...
|  56   | 1674  | 1127  | -1.6146e+03 | 1.77e-08 | 4.49e+00 | 3.55e-15 | 3.66e+33 |  1.00e-01   |   1   |
|  57   | 1705  | 1151  | -1.6146e+03 | 1.77e-09 | 4.49e+00 | 3.55e-15 | 3.66e+33 |  1.00e-01   |   1   |
|  58   | 1736  | 1151  | -1.6146e+03 | 1.00e+00 | 4.42e+00 | 3.55e-15 | 1.00e+00 |  2.00e-02   |   0   |
|  59   | 1767  | 1175  | -1.6146e+03 | 1.00e-01 | 4.42e+00 | 3.55e-15 | 1.00e+00 |  2.00e-02   |   1   |
|  60   | 1798  | 1199  | -1.6146e+03 | 1.00e-02 | 4.42e+00 | 3.55e-15 | 1.00e+00 |  2.00e-02   |   1   |
...
|  66   | 1984  | 1343  | -1.6146e+03 | 1.00e-08 | 4.42e+00 | 3.55e-15 | 1.00e+00 |  2.00e-02   |   1   |
|  67   | 2015  | 1367  | -1.6146e+03 | 1.00e-09 | 4.42e+00 | 3.55e-15 | 1.00e+00 |  2.00e-02   |   1   |
|  68   | 2046  | 1367  | -1.6146e+03 | 1.00e+00 | 4.36e+00 | 3.55e-15 | 1.00e+00 |  4.00e-03   |   0   |
|  69   | 2077  | 1391  | -1.6146e+03 | 1.00e-01 | 4.36e+00 | 3.55e-15 | 1.00e+00 |  4.00e-03   |   1   |
...
|  77   | 2325  | 1583  | -1.6146e+03 | 1.00e-09 | 4.36e+00 | 3.55e-15 | 1.00e+00 |  4.00e-03   |   1   |
|  78   | 2356  | 1583  | -1.6146e+03 | 1.00e+00 | 4.35e+00 | 3.55e-15 | 1.00e+00 |  8.00e-04   |   0   |
|  79   | 2387  | 1607  | -1.6146e+03 | 1.00e-01 | 4.35e+00 | 3.55e-15 | 1.00e+00 |  8.00e-04   |   1   |
...
|  87   | 2635  | 1799  | -1.6146e+03 | 1.00e-09 | 4.35e+00 | 3.55e-15 | 1.00e+00 |  8.00e-04   |   1   |
|  88   | 2666  | 1799  | -1.6146e+03 | 1.00e+00 | 4.34e+00 | 3.55e-15 | 1.00e+00 |  1.60e-04   |   0   |
|  89   | 2697  | 1823  | -1.6146e+03 | 1.00e-01 | 4.34e+00 | 3.55e-15 | 1.00e+00 |  1.60e-04   |   1   |
...
|  97   | 2945  | 2015  | -1.6146e+03 | 1.00e-09 | 4.34e+00 | 3.55e-15 | 1.00e+00 |  1.60e-04   |   1   |
|  98   | 2976  | 2015  | -1.6146e+03 | 1.00e+00 | 4.34e+00 | 3.55e-15 | 1.00e+00 |  3.20e-05   |   0   |
|  99   | 3007  | 2039  | -1.6146e+03 | 1.00e-01 | 4.34e+00 | 3.55e-15 | 1.00e+00 |  3.20e-05   |   1   |
...
|  167  | 5053  | 3527  | -1.6146e+03 | 1.00e-07 | 1.35e+01 | 2.12e-11 | 1.00e+00 |  2.05e-09   |   1   |
|  168  | 5084  | 3551  | -1.6146e+03 | 1.00e-08 | 1.35e+01 | 2.12e-11 | 1.00e+00 |  2.05e-09   |   1   |
|  169  | 5115  | 3575  | -1.6146e+03 | 1.00e-09 | 1.35e+01 | 2.12e-11 | 1.00e+00 |  2.05e-09   |   1   |
`xtol` termination condition is satisfied.
Number of iterations: 169, function evaluations: 5115, CG iterations: 3575, optimality: 1.35e+01, constraint violation: 2.12e-11, execution time: 3.8e+02 s.

很明显，一旦tr_radius < xtol，tr_radius 被重置为其默认值1，barrier param 被减少。一旦barrier param < barrier_tol（即1e-8）和tr_radius < xtol，优化成功终止。文档说关于barrier_tol

当存在不等式约束时，算法将仅在障碍参数小于barrier_tol 时终止。

这将解释不等式约束情况下的行为，但我所有的约束都是定义为字典的等式约束

con0 = 'type': 'eq', 'fun': constraint0

是否有人深入了解trust-constr 可以向我解释这一点？

【问题讨论】：

能否提供objective的实际实现。否则 SO 用户将能够重现您的问题

@tstanisl 感谢您的评论。 objective 是带有 66 个参数的加权和 objective 本身的定义，所有约束都是 80 行代码。我正在考虑在问题中发布目标，但是由于我使用无法共享的数据，因此无法重现该问题，因此我退出了它。但是，我的问题不是解决我的优化问题，而是了解minimize(method='trust-constr') 算法如何处理其中断条件，我认为这对社区来说也比我的具体问题更有趣

【参考方案1】：

您是否有具有上限的变量？也许求解器正在将这些作为约束来实现，例如 var < UPPER_BOUND。

（如果我有声誉分数，我会将此作为评论）

【讨论】：

好点。事实上，我对所有变量（下限和上限）都有界限，我通过bounds 键在问题中实现了这些界限。

如果它会在内部将所有边界转换为不等式约束，我可以想象最终会出现在我上面描述的场景中。我会检查source code

你是对的。源代码显示了它（请参阅我的答案）。感谢您的提示，您获得了赏金

【参考方案2】：

它通过PreparedConstraints 类和minimize(method='trust-constr') 内函数_minimize_trustregion_constr 中的initial_constraints_as_canonical 函数与变量边界到不等式约束的内部转换相关联。

定义的源代码可以在scipy/scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py中找到

负责的代码行是

if bounds is not None:
    if sparse_jacobian is None:
        sparse_jacobian = True
    prepared_constraints.append(PreparedConstraint(bounds, x0,
                                                   sparse_jacobian))

算法将定义的变量边界bounds 作为PreparedConstraint 附加到已经在prepared_constraints 中准备好的最初定义的约束列表中。后续线路

# Concatenate initial constraints to the canonical form.
c_eq0, c_ineq0, J_eq0, J_ineq0 = initial_constraints_as_canonical(
    n_vars, prepared_constraints, sparse_jacobian)

将每个边界转换为两个不等式约束（x > lb 和 x < ub）并因此返回两倍于边界数量的附加约束。

_minimize_trustregion_constr 然后检测那些不等式约束并正确选择算法tr_interior_point

# Choose appropriate method
if canonical.n_ineq == 0:
    method = 'equality_constrained_sqp'
else:
    method = 'tr_interior_point'

在下文中，该问题被视为最初包含不等式约束的问题，因此正确终止于问题中描述的xtol 条件和barrier_parameter 条件。

感谢@Dylan Black 的提示，他的回答获得了赏金。