Канал по Машинному обучению
@ai_machinelearning_big_data
Перейти
Введение в Оптимизацию с ограничениями на SciPy.
1. Введение
Оптимизация на SciPy
Мы выполняем множество задач по оптимизации в нашей повседневной жизни: находим кратчайший или быстрейший маршрут, чтобы добраться до пункта назначения, готовим список дел с ежедневными заданиями, упорядоченными по приоритету, пишем список покупок и многое другое.
Давайте выразим оптимизацию математически.
Давайте представим, что мы организуем поездку в другой город и пытаемся определить подходящее время отправления. В этом примере целевая функция f(x) представляет собой продолжительность поездки в зависимости от времени отправления x.
Мы можем сформулировать задачу оптимизации как определение минимального или максимального значения целевой функции. В нашем примере мы хотим определить время отправления, которое сведёт к минимуму продолжительность поездки:
В других сценариях мы можем захотеть максимизировать f(x). Например, когда цель представляет собой вероятность возврата инвестиций. Тем не менее, максимизация функции эквивалентна минимизации её отрицательного значения. Таким образом, можно сосредоточиться только на проблемах минимизации:
В реальных приложениях нам может потребоваться применить ограничения к нашей задаче оптимизации. Например, мы можем захотеть найти самый быстрый маршрут, но не хотим платить за проезд или путешествовать ночью. Мы определяем ограниченную оптимизацию как процесс минимизации целевой функции при некоторых логических условиях, которые могут отражать:
- ограничения реального мира;
- физический смысл входных переменных;
- контекстуальные обстоятельства.
В этом посте мы делимся примером оптимизации с использованием SciPy, популярной библиотеки Python для научных вычислений. В частности, мы исследуем наиболее распространённые типы ограничений: границы, линейные и нелинейные ограничения.
2. Реализация
2.1 Неограниченная оптимизация
Мы начинаем с простой задачи неограниченной оптимизации и позже добавляем ограничения к входным переменным.
Импортируйте необходимые библиотеки:
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import minimize, Bounds, LinearConstraint, NonlinearConstraintПредставьте себе следующую многомерную целевую функцию:
Её уклон относительно x₀ и x₁ равен…
def f(x):
'''Objective function'''
return 0.5*x[0]**2 + 2.5*x[1]**2 + 4 * x[0] * np.sin(np.pi * x[0]) + 5
def df(x):
'''Gradient of the objective function'''
return np.array([x[0] + 4 * np.sin(np.pi * x[0]) + 4 * np.pi * x[0] * np.cos(np.pi * x[0]), 5*x[1]])Давайте сгенерируем данные и понаблюдаем за значениями функции для x₀, x₁ ∈ [-1, 1]:
# Generate data
X0, X1 = np.meshgrid(np.linspace(-1, 1, 100), np.linspace(-1, 1, 100))
Z = f(np.stack([X0, X1]))
# Plot
fig = plt.figure(figsize=(30, 20))
# First subplot
ax = fig.add_subplot(1, 3, 1, projection='3d')
ax.contour3D(X0, X1, Z, 70, cmap='plasma')
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
# Second subplot
ax = fig.add_subplot(1, 3, 2, projection='3d')
ax.contour3D(X0, X1, Z, 70, cmap='plasma')
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.yaxis.set_ticklabels([])
ax.view_init(0, 80)
# Third subplot
ax = fig.add_subplot(1, 3, 3, projection='3d')
ax.contour3D(X0, X1, Z, 70, cmap='plasma')
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.zaxis.set_ticklabels([])
ax.view_init(89, -90);
В нашем примере целевая функция невыпуклая и обладает несколькими минимумами.
Это подразумевает, что, в зависимости от начальной точки, проблема может сходиться к другому минимизатору.
Мы можем решить проблему оптимизации, используя полезную функцию scipy.optimize.minimize следующим образом:
# Starting point
x_start = np.array([0.5, 0.5])
# Optimization
result = minimize(
f, x_start, method='trust-constr', jac=df)
result.x.round(3)Примечательно, что мы применяем метод trust-constr. Он позволяет оптимизировать функцию с учётом ограничений. Более подробная информация о методе доступна в документации к пакету и в “Trust-region methods” (Conn, Gould and Toint; 2000).
Приведённый выше фрагмент кода возвращает найденный минимизатор:
array([-0., 0.])Давайте построим его:
# Minimum from unconstrained optimization
min_x0, min_x1 = np.meshgrid(result.x[0], result.x[1])
min_z = f(np.stack([min_x0, min_x1]))
# Plot
fig = plt.figure(figsize=(30, 20))
# First subplot
ax = fig.add_subplot(1, 3, 1, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
# Second subplot
ax = fig.add_subplot(1, 3, 2, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.yaxis.set_ticklabels([])
ax.view_init(0, 80)
# Third subplot
ax = fig.add_subplot(1, 3, 3, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.zaxis.set_ticklabels([])
ax.view_init(89, -90);
Теперь мы можем поэкспериментировать с добавлением ограничений.
2.2 Границы
Давайте рассмотрим наш предыдущий пример о поиске самой быстрой поездки между двумя городами с использованием времени отправления в качестве входной переменной. Мы можем ожидать большего или меньшего трафика в зависимости от часа дня. Стремясь свести к минимуму продолжительность поездки, модель может также предложить, например, путешествовать ночью.
Хотя это может привести к самой короткой поездке, мы можем предпочесть путешествовать в дневное время. Таким образом, мы могли бы попросить модель найти кратчайшую поездку, учитывая время отправления только с 7.00 утра до 6.00 вечера.
Вот тут-то и появляются границы. Границы – это просто ограничения на равенство или неравенство входных переменных. Они позволяют оценивать целевую функцию исключительно в пределах указанных диапазонов.
В нашем случае предположим, что у нас есть следующие приемлемые значения для x₀ и x₁:
Мы можем легко передать эти значения объекту Bounds и выполнить новый эксперимент по оптимизации следующим образом:
lim = [0.25, 0.30, 0.75, 0.8]
bounds = Bounds([lim[0], lim[1]], # [min x0, min x1]
[lim[2], lim[3]]) # [max x0, max x1]
result = minimize(
f, x_start, method='trust-constr', jac=df, bounds=bounds)
result.x.round(3)Задача оптимизации теперь приводит к другому решению, поскольку предыдущий массив точек ([0, 0]) не попадает в допустимую область:
array([0.25, 0.301])Наконец, мы можем построить новый минимум и допустимую область и наблюдать область, в которой оценивалось f (x):
# Feasible region (bounds)
X0_bound, X1_bound = np.meshgrid(np.linspace(lim[0], lim[2], 20), np.linspace(lim[1], lim[3], 20))
Z_bound = f(np.stack([X0_bound, X1_bound]))
# New minimum within bounds
min_x0_bounds, min_x1_bounds = np.meshgrid(result.x[0], result.x[1])
min_z_bounds = f(np.stack([min_x0_bounds, min_x0_bounds]))
# Plot
fig = plt.figure(figsize=(30, 20))
# First subplot
ax = fig.add_subplot(1, 3, 1, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_bounds, min_x1_bounds, min_z_bounds, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_bound, X1_bound, Z_bound, color='black', alpha=0.6)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
# Second subplot
ax = fig.add_subplot(1, 3, 2, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_bounds, min_x1_bounds, min_z_bounds, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_bound, X1_bound, Z_bound, color='black', alpha=0.6)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.yaxis.set_ticklabels([])
ax.view_init(0, 80)
# Third subplot
ax = fig.add_subplot(1, 3, 3, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_bounds, min_x1_bounds, min_z_bounds, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_bound, X1_bound, Z_bound, color='black', alpha=0.6)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.zaxis.set_ticklabels([])
ax.view_init(89, -90);
2.3 Линейные ограничения
Линейные ограничения определяют линейные отношения между переменными оптимизации. Например, давайте представим, что это x₀ и x₁:
Мы можем легко переписать эти условия как линейную систему и передать их объекту LinearConstraint перед выполнением задачи оптимизации:
linear_constraint = LinearConstraint(
[[1, 0], [1, -1]], [0.25, -np.inf], [np.inf, 0.1])
result = minimize(
f, x_start, method='trust-constr', jac=df, constraints=linear_constraint)
result.x.round(3)Новое решение заключается в…
array([0.25, 0.15])Допустимая область f(x) соответствует части пространства, ограниченной пересечением гиперплоскостей. Давайте наметим эти границы:
# Linear constraints: first hyperplane
X0_lin_1 = np.repeat(0.25, 20)
X1_lin_1, Z_lin_1 = np.meshgrid(np.linspace(-1, 1, 20), np.linspace(4, 10, 20))
# Linear constraints: second hyperplane
X1_lin_2 = np.linspace(-1, 1, 20)
X0_lin_2 = X1_lin_2 + 0.1
# New minimum with linear constraints
min_x0_lin_constr, min_x1_lin_constr = np.meshgrid(result.x[0], result.x[1])
min_z_lin_constr = f(np.stack([min_x0_lin_constr, min_x0_lin_constr]))
# Plot
fig = plt.figure(figsize=(30, 20))
# First subplot
ax = fig.add_subplot(1, 3, 1, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_lin_constr, min_x1_lin_constr, min_z_lin_constr, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_lin_1, X1_lin_1, Z_lin_1, color='green', alpha=0.3)
ax.plot_surface(X0_lin_2, X1_lin_2, Z_lin_1, color='yellow', alpha=0.2)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
# Second subplot
ax = fig.add_subplot(1, 3, 2, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_lin_constr, min_x1_lin_constr, min_z_lin_constr, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_lin_1, X1_lin_1, Z_lin_1, color='green', alpha=0.2)
ax.plot_surface(X0_lin_2, X1_lin_2, Z_lin_1, color='yellow', alpha=0.2)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.yaxis.set_ticklabels([])
ax.view_init(0, 80)
# Third subplot
ax = fig.add_subplot(1, 3, 3, projection='3d')
ax.contour3D(X0, X1, Z, 40, cmap='plasma')
ax.scatter(min_x0, min_x1, min_z, marker='o', color='red', linewidth=10)
ax.scatter(min_x0_lin_constr, min_x1_lin_constr, min_z_lin_constr, marker='o', color='blue', linewidth=10)
ax.plot_surface(X0_lin_1, X1_lin_1, Z_lin_1, color='green', alpha=1)
ax.plot_surface(X0_lin_2, X1_lin_2, Z_lin_1, color='yellow', alpha=1)
ax.set_xlabel('$x_{0}$')
ax.set_ylabel('$x_{1}$')
ax.set_zlabel('$f(x)$')
ax.axes.zaxis.set_ticklabels([])
ax.view_init(89, -90);
2.4 Нелинейные ограничения
Мы также можем исследовать целевую функцию в пределах области, определяемой нелинейными ограничениями, используя объект NonlinearConstraint. Предположим, что это x₀ и x₁:
Мы оптимизируем f(x) следующим образом:
non_linear_eq= lambda x: x[0]**2 + x[1]**2
non_linear_constr = NonlinearConstraint(
non_linear_eq, 0.2, np.inf)
result = minimize(
f, np.array([0.5, 1]), method='trust-constr', jac=df, constraints=non_linear_constr)
result.x.round(3)array([-0., 0.447])Аналогично предыдущим примерам, мы могли бы наблюдать цель и найденный минимизатор, учитывая текущее ограничение:
2.5 Совместное применение различных типов ограничений
Мы можем комбинировать границы, а также линейные и нелинейные ограничения следующим образом:
result = minimize(
f,
x_start,
method='trust-constr',
jac=df,
bounds=bounds,
constraints=[linear_constraint, non_linear_constr]
)
result.x.round(3)array([0.25, 0.371])Отметим, что не все методы оптимизации поддерживают границы и/или ограничения. Дополнительную информацию можно найти в документации к пакету.
3. Заключение
В этом посте мы рассмотрели различные типы Оптимизации с ограничениями. В частности, мы поделились практическими примерами Python с использованием библиотеки SciPy. Примеры сопровождаются графиками, которые позволяют визуально проверить различные ограничения.
@data_analysis_ml – анализ данных
+1
3
+1
6
+1
2
+1
+1
