Сципи: свести к минимуму количество нарушений

Question

Я пытаюсь использовать scipy.optimize.minimize с простыми a <= x <= b границами. Тем не менее, часто бывает, что моя целевая функция оценивается за пределами границ. Насколько я понимаю, это происходит, когда minimize пытается определить градиент целевой функции на границе.

---

Минимальный пример:

import math
import numpy as np
from scipy.optimize import Bounds, minimize

constraint = Bounds([-1, -1], [1, 1], True)
def fun(x):
    print(x)
    return -math.exp(-np.dot(x,x))
result = minimize(fun, [-1, -1], bounds=constraint)

Выходные данные показывают, что минимизатор переходит к точке [1,1] и затем пытается оценить на [1.00000001, 1]:

[-1. -1.]
[-0.99999999 -1.        ]
[-1.         -0.99999999]
[-0.72932943 -0.72932943]
[-0.72932942 -0.72932943]
[-0.72932943 -0.72932942]
[-0.22590689 -0.22590689]
[-0.22590688 -0.22590689]
[-0.22590689 -0.22590688]
[1. 1.]
[1.00000001 1.        ]
[1.         1.00000001]
[-0.03437328 -0.03437328]
...

Конечно, в этом примере нет проблем, так как fun можно оценить и там. Но это не всегда так ...

---

В моей актуальной задаче минимум не может быть на границе, и у меня есть простой способ добавить эпсилон к границам. Но можно ожидать, что должно быть простое решение этой проблемы, которое также работает, если минимум может быть на границе?

PS: Было бы странно, если бы я был первым, кто столкнулся с этой проблемой - извините, если этот вопрос где-то задавался ранее, но я нигде его не нашел.

Answer 1

Как обсуждалось здесь (спасибо @ "Welcome to Stack Overflow" за комментарий, направляющий меня туда), проблема в том, что процедура градиента не учитывает границы.Я написал новый, который выполняет эту работу:

import math
import numpy as np
from scipy.optimize import minimize

def gradient_respecting_bounds(bounds, fun, eps=1e-8):
    """bounds: list of tuples (lower, upper)"""
    def gradient(x):
        fx = fun(x)
        grad = np.zeros(len(x))
        for k in range(len(x)):
            d = np.zeros(len(x))
            d[k] = eps if x[k] + eps <= bounds[k][1] else -eps
            grad[k] = (fun(x + d) - fx) / d[k]
        return grad
    return gradient

bounds = ((-1, 1), (-1, 1))
def fun(x):
    print(x)
    return -math.exp(-np.dot(x,x))
result = minimize(fun, [-1, -1], bounds=bounds,
                  jac=gradient_respecting_bounds(bounds, fun))

---

Обратите внимание, что это может быть немного менее эффективно, потому что fun(x) теперь оценивается дважды в каждой точке.Это кажется неизбежным, соответствующий фрагмент из _minimize_lbfgsb в lbfgsb.py:

if jac is None:
    def func_and_grad(x):
        f = fun(x, *args)
        g = _approx_fprime_helper(x, fun, epsilon, args=args, f0=f)
        return f, g
else:
    def func_and_grad(x):
        f = fun(x, *args)
        g = jac(x, *args)
        return f, g

Как видите, значение f может быть повторно использовано только внутренней функцией _approx_fprime_helper.