Оптимизация
Contents
Оптимизация¶
Подмодуль scipy.optimize содержит в себе множество методов для решения задач оптимизации (минимизация, максимизация).
Минимизация скалярной функции одного аргумента¶
Функция scipy.optimize.minimize_scalar позволяет минимизировать функцию \(f: \mathbb{R} \to \mathbb{R}\):
import numpy as np
from scipy import optimize
from matplotlib import pyplot as plt
def f(x):
return (x - 2) * x * (x + 2)**2
sol = optimize.minimize_scalar(f)
print(sol)
x = np.linspace(-4, 3, 100)
plt.plot(x, f(x))
plt.scatter(sol.x, f(sol.x), color="r")
plt.xlabel("x")
plt.ylabel("y = f(x)")
plt.legend(["Минимизируемая функция", "Найденный минимум"])
fun: -9.914949590828147
message: '\nOptimization terminated successfully;\nThe returned value satisfies the termination criteria\n(using xtol = 1.48e-08 )'
nfev: 15
nit: 11
success: True
x: 1.2807764040333458
<matplotlib.legend.Legend at 0x225a09f8f10>

У функции \(f(x) = (x - 1)x(x + 2)^2\) два локальных минимума и функция minimize_scalar
вернула из них глобальный. Можно сузить поиск область поиска минимума используя метод bounded
.
a, b = -3, 1
sol = optimize.minimize_scalar(f, bounds=(a, b), method="bounded")
plt.plot(x, f(x))
plt.scatter(sol.x, f(sol.x), color="r")
plt.xticks([a, b])
upper = np.array([100, 100])
lower = np.array([-30, -30])
plt.fill_between([a, b], lower, upper, where=upper>0, facecolor='blue', alpha=0.10)
plt.xlabel("x")
plt.ylabel("y = f(x)")
plt.legend(["Минимизируемая функция", "Найденный минимум"])
<matplotlib.legend.Legend at 0x225a2b02280>

Минимизация функции многих переменных¶
Функция scipy.optimize.minimize предназначена для минимизации функции многих переменных. Рассмотрим её на примере функции Розенброка
Её минимум равен 0 при \(x_1 = 1\) и \(x_2 = 1\). Эта функция часто применяется для тестирования алгоритмов оптимизации и поэтому присутствует в SciPy
под именем scipy.optimize.rosen.
import numpy as np
from scipy import optimize
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
f = optimize.rosen
jac = optimize.rosen_der
hess = optimize.rosen_hess
x = np.linspace(0, 2, 200)
X, Y = np.meshgrid(x, x)
fig, ax = plt.subplots(subplot_kw={"projection": "3d"})
fig.set_size_inches((10, 9))
ax.plot_surface(X, Y, f([X, Y]), cmap="coolwarm")
plt.show()

res = optimize.minimize(f, [2, 2])
print(res)
fig, ax = plt.subplots()
fig.set_size_inches((10, 9))
ax.contourf(X, Y, f([X, Y]), cmap="coolwarm")
ax.scatter(res.x[0], res.x[1], color="r")
plt.show()
fun: 1.8932893809017893e-11
hess_inv: array([[0.51675994, 1.03186494],
[1.03186494, 2.0655726 ]])
jac: array([ 5.27380711e-06, -2.50575298e-06])
message: 'Optimization terminated successfully.'
nfev: 105
nit: 30
njev: 35
status: 0
success: True
x: array([0.99999565, 0.99999129])

Как видно из сообщения, метод успешно сошелся и алгоритму потребовалось 30 итераций и 105 вызовов функции f
, чтобы сойтись. Передадим этому методу функцию, вычисляющую матрицу градиент искомой функции.
res = optimize.minimize(f, [2, 2], jac=jac)
res
fun: 3.271435067604758e-14
hess_inv: array([[0.51671918, 1.03179112],
[1.03179112, 2.06544014]])
jac: array([ 5.26067293e-06, -2.49959031e-06])
message: 'Optimization terminated successfully.'
nfev: 35
nit: 30
njev: 35
status: 0
success: True
x: array([1.00000013, 1.00000025])
Количество вызовов функции снизилось до 35, но количество итераций осталось прежним. Это достигается за счет того, что метод минимизации вычисляет производную, вызывая переданную функцию jac
, а до этого он её оценивал численно, что требовало дополнительных вызовов f
.
Передадим ещё информацию о второй производной функции f
, т.е. матрицу Гессе. Метод минимизации по умолчанию (BFGS
) — метод первого порядка, чтобы задействовать матрицу Гессе, необходимо использовать метод второго порядка, например, dogleg
res = optimize.minimize(f, [2, 2], jac=optimize.rosen_der, hess=optimize.rosen_hess, method="dogleg")
print(res)
fun: 7.335116582097013e-13
hess: array([[ 802.00135535, -400.00033735],
[-400.00033735, 200. ]])
jac: array([ 7.65326055e-06, -2.98326426e-06])
message: 'Optimization terminated successfully.'
nfev: 16
nhev: 13
nit: 15
njev: 14
status: 0
success: True
x: array([1.00000084, 1.00000167])
Методу второго порядка потребовалось всего 15 итераций и 16 вызовов функции f
, чтобы сойтись.