Я долго готовился и собирал материал, надеюсь в этот раз получилось лучше. Эту статью посвящаю основным методам решения задач математической оптимизации с ограничениями, так что если вы слышали, что симплекс-метод — это какой-то очень важный метод, но до сих пор не знаете, что он делает, то возможно эта статья вам поможет.

P. S. Статья содержит математические формулы, добавленные макросами хабраредактора. Говорят, что они иногда не отображаются. Также есть много анимаций в формате gif.

Преамбула


Задача математической оптимизации — это задача вида «Найти в множестве $\mathcal{K}$ элемент $x^*$ такой, что для всех $x$ из $\mathcal{K}$ выполняется $f(x^*)\leq f(x)$», что в научной литературе скорее будет записано как-то так

$\begin{array}{rl} \mbox{минимизировать } & f(x) \\ \mbox{при условии } & x\in \mathcal{K}. \end{array} $


Исторически так сложилось, что популярные методы такие как градиентный спуск или метод Ньютона работают только в линейных пространствах (причем желательно простых, например $\mathbb{R}^n$). На практике же часто встречаются задачи, где нужно найти минимум не в линейном пространстве. Например нужно найти минимум некоторой функции на таких векторах $x=(x_1, \ldots, x_n)$, для которых $x_i\geq 0$, это может быть обусловлено тем, что $x_i$ обозначают длины каких-либо объектов. Или же например, если $x$ представляют координаты точки, которая должна быть на расстоянии не больше $r$ от $y$, т. е. $\|x-y\|\leq r$. Для таких задач градиентный спуск или метод Ньютона уже напрямую не применить. Оказалось, что очень большой класс задач оптимизации удобно покрывается «ограничениям», подобными тем, что я описал выше. Иначе говоря, удобно представлять множество $\mathcal{K}$ в виде системы равенств и неравенств

$ \begin{array}{l} g_i(x)\leq 0,~1\leq i \leq m, \\ h_i(x)=0,~1\leq i\leq k. \end{array} $


Задачи минимизации над пространством вида $\mathbb{R}^n$ таким образом стали условно называть «задачами без ограничений» (unconstrained problem), а задачи над множествами, заданными наборами равенств и неравенств — «задачами с ограничениями» (constrained problem).

Технически, совершенно любое множество $\mathcal{K}$ можно представить в виде одного равенства или неравенство с помощью индикатор-функции, которая определяется как

$I_\mathcal{K}(x)= \begin{cases} 0, & x\notin \mathcal{K} \\ 1, & x\in \mathcal{K}, \end{cases} $


однако такая функция не обладает разными полезными свойствами (выпуклость, дифференцируемость и т. п.). Тем не менее, часто можно представить $\mathcal{K}$ в виде нескольких равенств и неравенств, каждое из которых такими свойствами обладает. Основная теория подведена под случай

$\begin{array}{rl} \mbox{минимизировать } & f(x) \\ \mbox{при условии } & g_i(x)\leq 0,~1\leq i\leq m \\ & Ax=b, \end{array} $


где $f, g_i$ — выпуклые (но не обязательно дифференцируемые) функции, $A$ — матрица. Для демонстрации работы методов я буду использовать два примера:

  1. Задача линейного программирования

    $\begin{array}{rl} \mbox{минимизировать } & -2&x~~~- &y \\ \mbox{при условии } &-1.0 & ~x -0.1 & ~y \leq -1.0 \\ & -1.0 & ~x + ~0.6 &~ y \leq -1.0 \\ & -0.2 & ~x + ~1.5 & ~y \leq -0.2\\ & ~0.7 &~x + ~0.7 & ~y \leq 0.7\\ &~2.0 & ~x -0.2 & ~y \leq 2.0\\ &~0.5 & ~x -1.0 & ~y \leq 0.5\\ &-1.0 & ~x -1.5& ~ y \leq -1.0\\ \end{array} $


    По сути эта задача состоит в том, чтобы найти самую дальнюю точку многоугольника в направлении (2, 1), решение задачи — точка (4.7, 3.5) — самая «правая» в многоугольнике). А вот собственно и сам многоугольник

  2. Минимизация квадратичной функции с одним квадратичным ограничением

    $\begin{array}{rl} \mbox{минимизировать } & 0.7 (x - y)^2 + 0.1 (x + y)^2\\ \mbox{при условии } & (x-4)^2+(y-6)^2\leq 9 \end{array} $




Симплекс-метод


Из всех методов, которые я покрываю этим обзором, симплекс-метод наверно является самым известным. Метод был разработан специально для линейного программирования и единственный из представленных достигает точного решения за конечное число шагов (при условии, что для вычислений используется точная арифметика, на практике это обычно не так, но в теории возможно). Идея симплекс-метода состоит из двух частей:

  1. Системы линейных неравенств и равенств задают многомерные выпуклые многогранники (политопы). В одномерном случае это точка, луч или отрезок, в двумерном — выпуклый многоугольник, в трехмерном — выпуклый многогранник. Минимизация линейной функции — это по сути нахождение самой «дальней» точки в определенном направлении. Думаю, интуиция должна подсказывать, что этой самой дальней точкой должна быть некая вершина, и это действительно так. В общем случае, для системы из $m$ неравенств в $n$-мерном пространстве вершина — это точка, удовлетворяющая системе, для которой ровно $n$ из этих неравенств обращаются в равенства (при условии, что среди неравенств нет эквивалентных). Таких точек всегда конечное число, хоть их и может быть очень много.
  2. Теперь у нас есть конечный набор точек, вообще говоря можно их просто взять и перебрать, то есть сделать что-то такое: для каждого подмножества из $n$ неравенств решить систему линейных уравнений, составленных на выбранных неравенствах, проверить, что решение подходит в исходную систему неравенств и сравнить с другими такими точками. Это довольно простой неэффективный, но рабочий метод. Симплекс-метод вместо перебора двигается от вершины к вершине по ребрам таким образом, чтобы значений целевой функции улучшалось. Оказывается, если у вершины нет «соседей», в которых значений функции лучше, то она оптимальна.

Симплекс-метод является итеративным, то есть он последовательно по чуть-чуть улучшает решение. Для таких методов нужно с чего-то начинать, в общем случае это делается с помощью решения вспомогательной задачи

$\begin{array}{rl} \mbox{минимизировать } & s \\ \mbox{при условии } & g_i(x)\leq s,~1\leq i\leq m \\ \end{array} $


Если для решения этой задачи $x^*, s^*$ такое, что $s^*\leq 0$, то выполняется $g_i(x^*)\leq s\leq 0$, иначе исходная задача вообще задана на пустом множестве. Чтобы решить вспомогательную задачу, можно также использовать симплекс-метод, начальной же точкой можно взять $s=\max_ig_i(x)$ с произвольным $x$. Нахождение начальной точки можно условно назвать первой фазой метода, нахождение решение исходной задачи можно условно назвать второй фазой метода.

Траектория двухфазового симплекс-метода
Траектория была сгенерирована с помощью scipy.optimize.linprog.


Проективный градиентный спуск


Про градиентный спуск я недавно писал отдельную статью, в которой в том числе кратко описал и этот метод. Сейчас этот метод вполне себе живой, но изучается как часть более общего проксимального градиентного спуска. Сама идея метода совсем банальна: если мы применяем градиентный спуск к выпуклой функции $f$, то при правильном выборе параметров получаем глобальный минимум $f$. Если же после каждого шага градиентного спуска корректировать полученную точку, взяв вместо нее её проекцию на замкнутое выпуклое множество $\mathcal{K}$, то в результате мы получим минимум функции $f$ на $\mathcal{K}$. Ну или более формально, проективный градиентный спуск — это алгоритм, который последовательно вычисляет

$ \begin{cases} x_{k+1} = y_k - \alpha_k\nabla f(y_k) \\ y_{k+1} = P_{\mathcal{K}}(x_{k+1}), \end{cases} $


где

$ P_{\mathcal{K}}(x)=\mbox{argmin}_{y\in\mathcal{K}}\|x-y\|. $


Последнее равенство определяет стандартный оператор проекции на множество, по сути это функция, которая по точке $x$ вычисляет ближайшую к ней точку множества $\mathcal{K}$. Роль расстояния здесь играет $\|\ldots\|$, стоит отметить, что здесь можно использовать любую норму, тем не менее проекции с разными нормами могут отличаться!

На практике проективный градиентный спуск используется только в особых случаях. Основная его проблема состоит в том, что вычисление проекции может быть еще более сложной задачей, чем исходная, а её нужно вычислять много раз. Самый распространенный случай, для которого удобно применять проективный градиентный спуск — это «коробочные ограничения», которые имеют вид

$ \ell_i\leq x_i\leq r_i,~~1\leq i\leq n. $


В этом случае проекция вычисляется очень просто, по каждой координате получается

$ [P_{\mathcal{K}}(x)]_i= \begin{cases} r_i, & x_i>r_i \\ \ell_i, & x_i < \ell_i \\ x_i, & x_i\in [\ell_i, r_i]. \end{cases} $


Применение проективного градиентного спуска для задач линейного программирования совершенно бессмысленно, тем не менее если это все-таки сделать, то выглядеть будет как-то так

Траектория проективного градиентного спуска для задачи линейного программирования


А вот как выглядит траектория проективного градиентного спуска для второй задачи, если

выбирать большой размер шага


и если

выбирать небольшой размер шага


Метод эллипсоидов


Этот метод примечателен тем, что является первым полиномиальным алгоритмом для задач линейного программирования, его можно считать многомерным обобщением метода бисекции. Я начну с более общего метода разделяющей гиперплоскости:

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

Для задач оптимизации построение «разделяющей гиперплоскости» основано на следующем неравенстве для выпуклых функций

$ f(y)\geq f(x)+\nabla f(x)^T(y-x). $


Если зафиксировать $x$, то для выпуклой функции $f$ полупространство $\nabla f(x)^T(y-x)\geq 0$ содержит только точки со значением не меньше, чем в точке $x$, а значит их можно отсечь, так как эти точки не лучше, чем та, что мы уже нашли. Для задач с ограничениями можно аналогичным образом избавиться от точек, которые гарантированно нарушают какое-то из ограничений.

Самый простой вариант метода разделяющей гиперплоскости — это просто отсекать полупространства без добавления каких-либо точек. В результате на каждом шаге у нас будет некий многогранник. Проблема этого метода в том, что количество граней многогранника скорее всего будет возрастать от шага к шагу. Более того, оно может расти экспоненциально.

Метод эллипсоидов собственно на каждом шаге хранит эллипсоид. Точнее, после проведения гиперплоскости строится эллипсоид минимального объема, который содержит одну из частей исходного. Этого получается добиться за счет добавления новых точек. Эллипсоид всегда можно задать положительно определенной матрицей и вектором (центром эллипсоида) следующим образом

$ \mathcal{E}(P, x)=\{z~|~(z-x)^TP(z-x)\leq 1\}. $


Построение минимального по объему эллипсоида, содержащего пересечение полупространства и другого эллипсоида, можно осуществить с помощью в меру громоздких формул. К сожалению на практике этот метод оказался все еще на так хорош, как симплекс-метод или метод внутренней точки.

А вот собственно как он работает для

линейного программирования


и для

квадратичного программирования


Метод внутренней точки


Этот метод имеет долгую историю развития, одни из первых предпосылок появились примерно в то же время, когда был разработан симплекс-метод. Но в то время он был еще недостаточно эффективен, чтобы использоваться на практике. Позднее в 1984 был разработан вариант метода специально для линейного программирования, который был хорош как в теории, так и на практике. Более того, метод внутренней точки не ограничен только линейным программированием в отличие от симплекс-метода, и сейчас он является основным алгоритмом для задач выпуклой оптимизации с ограничениями.

Базовая идея метода — замена ограничений на штраф в виде так называемой барьерной функции. Функция $F:Int\mathcal{K}\rightarrow \mathbb{R}$ называется барьерной функцией для множества $\mathcal{K}$, если

$ F(x)\rightarrow +\infty~\mbox{при } x\rightarrow\partial \mathcal{K}. $


Здесь $Int\mathcal{K}$ — внутренность $\mathcal{K}$, $\partial\mathcal{K}$ — граница $\mathcal{K}$. Вместо исходной задачи предлагается решать задачу

$ \mbox{минимизировать по }x~~\varphi(x, t)=tf(x)+F(x). $


$F$ и $\varphi$ заданы только на внутренности $\mathcal{K}$ (по сути отсюда и название), свойство барьера гарантирует, что у $\varphi$ минимум по $x$ существует. Более того, чем больше $t$, тем больше влияние $f$. При достаточно разумных условиях можно добиться того, что если устремить $t$ к бесконечности, то минимум $\varphi$ будет сходиться к решению исходной задачи.

Если множество $\mathcal{K}$ задано в виде набора неравенств $g_i(x)\leq 0,~1\leq i\leq m$, то стандартным выбором барьерной функции является логарифмический барьер

$ F(x)=-\sum_{i=1}^m\ln(-g_i(x)). $


Точки минимума $x^*(t)$ функции $\varphi(x, t)$ для разных $t$ образует кривую, которую обычно называют центральный путь, метод внутренний точки как бы пытается следовать этому пути. Вот так он выглядит для

Примера с линейным программированием

Аналитический центр — это просто $x^*(0)$

Наконец сам метод внутренней точки имеет следующий вид

  1. Выбрать начальное приближение $x_0$, $t_0>0$
  2. Выбрать новое приближение методом Ньютона

    $ x_{k+1}=x_k-[\nabla^2_x \varphi(x_k, t_k)]^{-1}\nabla_x \varphi(x_k, t_k) $


  3. Увеличить $t$

    $ t_{k+1}=\alpha t_k,~~\alpha>1 $



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

Задача линейного программирования

Прыгающая черная точка — это $x^*(t_k)$, т.е. точка, к которой мы пытаемся приблизиться шагом метода Ньютона на текущем шаге.

Задача квадратичного программирования