Выберите методы безусловной оптимизации ошибки обучения - Исправление ошибок и поиск оптимальных решений проблем

Обзор методов численной оптимизации. Безусловная оптимизация: метод линий

Время прочтения
24 мин

Просмотры 22K

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

Введение

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

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

Некоторые посты по оптимизации на Хабре от других авторов:

Обзор градиентных методов в задачах математической оптимизации
Обзор основных методов математической оптимизации для задач с ограничениями
Метод BFGS или один из самых эффективных методов оптимизации. Пример реализации на Python
Квазиньютоновские методы, или когда вторых производных для Атоса слишком много
Метод оптимизации Trust-Region DOGLEG. Пример реализации на Python

Раздел численных методов, изучающий оптимизацию общих нелинейных функций, называется нелинейным программированием, или NLP (конечно, сейчас на слуху у всех другое NLP — Natural Language Processing, но это не тот случай).

Кто хочет сразу перейти к сердцевине — можно перескочить к главе «Характеристики методов».

Типы оптимизаций

Безусловная vs условная

Почитать

Оптимизация может быть безусловной или условной. Последнее означает, что при нахождении минимума должны выполняться некоторые условия, например, равенства и/или неравенства. Безусловная оптимизация работает с целевой функцией как таковой, предполагая, что у нас есть полная свобода в поиске её минимума повсюду в её области определения.

Следующая оптимизация

$underset{x}{min} (x-1)^2$

является безусловной и имеет решение $x = 1$ , в то время как добавление неравенства как показано

$begin{cases}underset{x}{min} (x-1)^2 \ xleqslant 0end{cases}$

делает оптимизацию условной и порождает решение $x = 0$ .

В механике равенства соответствуют двусторонним ограничениям. Это может быть, например, условие, при котором материальная точка никогда не покидает определённую поверхность, пока её потенциальная энергия сводится к минимуму. Односторонние ограничения, с другой стороны, соответствуют неравенствам и распространены в негладкой механике, например: столкновения, непроникновение поверхностей сквозь друг друга, трение и т. д. Такие системы рассматриваются с помощью множителей Лагранжа, представляющих силы реакции со стороны поверхностей.

В данном посте мы будем рассматривать безусловную оптимизацию, а к условной, может быть, перейдём в другой раз.

Локальная vs глобальная

Гладкая vs негладкая

Почитать

Разделение «Гладкая vs негладкая» может быть двояко. В одном случае оно может означать, что мы априори работаем только с гладкими функциями, но сами методы оптимизации для них делятся на те, которые используют напрямую или не используют сами производные функций в своей работе. Последние называют безградиентными (derivative-free). Первые используют как минимум информацию первого порядка, т.е., кроме функций, ещё её производную, а вторые — информацию только нулевого порядка, т.е. саму функцию. Например, метод наискорейшего спуска использует производные, в то время как метод доверительной области с интерполяцией не использует, но при этом он предполагает, что наша функция достаточно гладкая. Однако ничто, конечно же, не запрещает пробовать безградиентные методы на негладких функциях, авось будет удачно. На практике люди не сильно придирчивы и часто так и поступают без лишних разбирательств.

В другом случае смысл «Гладкая vs негладкая» состоит в том, что идёт разделение между самими функциями: какие-то гладкие, а какие-то только непрерывные или даже разрывные. Понятно, что считать производные для последних не всегда будет корректным, или её расчёт крайне неточен. Например, метод Нелдера-Мида не требует какой-то информации о гладкости функции выше непрерывности для своей работы и вполне себе подходит для работы с функциями класса $C^0$ (само по себе это, конечно, не значит, что на более гладких функциях он не будет работать надёжнее).

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

Ниже мы рассмотрим методы, которые нацелены как минимум на непрерывно-дифференцируемые функции, т.е. из класса $C^1$ . Последнее условие требуется, чтобы иметь возможность положиться на обширную доступную теорию гладких функций, а также использовать как минимум градиент функции для ускорения процесса поиска. Как известно, если $mathbf{x}_*$ является локальным минимизатором, то $nabla f(mathbf{x}_*) = mathbf{0}$ . Обратное, к сожалению, неверно.

Трудоёмкость целевой функции

Почитать

Есть задачи оптимизации, в фокусе которых «тяжёлые» целевые функции: вычисление значения такой функции в одной точке может потребовать часы, если не дни! Например, оптимизация формы крыла самолёта при каждом вычислении целевой функции шарашит газодинамический решатель на сетке! Понятно, что в таких условиях на ускорение линейной алгебры в самом алгоритме оптимизации становится по боку — лишь бы поменьше вызывать целевую функцию. Это в противовес задачам, вычисление целевых функций в которых относительно дёшево, но количество переменных зашкаливает, — в таких задачах «линейка» должна быть прилизана по максимуму.

Дополнительные трудности могут возникать в задачах, в которых целевая функции вычислется неточно, а с некоторым шумом, вносящим неопределённость.

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

Разделений ещё больше

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

Например, есть огромная подобласть Mixed-Integer Programming, в которой рассматриваются дискретные сценарии. Есть недетерминированная (хаотическая) оптимизация, использующая вероятностные подходы в своей работе. Есть робастная оптимизация, в которую заложены параметры, не являющиеся переменными. Есть оптимизация динамических систем, в которых присутствует эволюция во времени. Есть всевозможные мета-эвристические методы: метод имитируемого отжига (simulated annealing), генетические алгоритмы, методы эволюции роя. Есть оптимизация неопределённой логики (fuzzy logic optimization).

Есть бог знает что! Я знаю лишь какую-то часть из представленного, и обо всём точно не расскажу никогда.

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

Наши реалии

Почитать

Все методы, представленные в этом посте, итеративны по своей природе. Это значит, что, начиная с предоставленного пользователем начального положения $mathbf{x}_0$ в области определения , любой метод последовательно генерирует специальное направление на основе предыдущих данных и самой функции. Метод полагает, что направление указывает на точку возможного локального минимума , и поэтому смещает положение в этом направлении.

Численная минимизация непрерывной невыпуклой функции не дает гарантии, что полученное решение действительно является точкой локального минимума. Критерий останова метода проверяет магнитуду градиента функции и прерывает алгоритм в случае, если эта магнитуда ниже некоторого малого значения . Теоретически доказанные результаты сходимости с таким критерием существуют для многих методов, но его можно обмануть, заменив минимизатор произвольной стационарной точкой, что приведёт к невозможности предоставить истинное решение.
Например, для функции начальное приближение уже будет ответом, ввиду того, что , но ответ этот фейковый, т.к. представляет собой локальный максимум, в то время как реальная точка минимума — это любая из $x = pmfrac{1}{sqrt{2}}$ .

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

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

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

Что в этом случае делать? Ничего. Скрестить пальцы и надеяться на лучшее.

Характеристики методов

Рассматриваемые ниже методы минимизации можно разделить по множеству актуальных характеристик. В этом посте нам могут быть интересны следующие:

Глобальная сходимость.
Скорость локальной сходимости.
Размеры задачи, на которую они нацелены.
Требуется ли хранить в памяти матрицы или нет.
Требуется ли матрица Гессе или нет.
Требуется ли масштабирование или нет.

Первые два пункта являются теоретическими, другие же более приземлены и непосредственно относятся к практической реализации. Пройдёмся по каждому из пунктов.

Говорят, что метод имеет свойство глобальной сходимости, если его итерации сходятся к локальному минимизатору независимо от начального положения (не путать с поиском глобального минимизатора!). Скорость локальной сходимости показывает, насколько быстро метод «настигает» минимизатор, когда он уже достаточно близко подошёл к нему. Удивительно, но эти две характеристики метода, будучи обе очень желательными в реальных вычислениях, могут конфликтовать между собой в многих методах оптимизации. Так называемые гибридные методы умело используют оба этих преимущества.

Размер задачи сужает выбор до тех методов, которые мы можем себе позволить. Для крупных задач интенсивное использование памяти может быть недопустимым, и для таких задач используют методы, возникшие из таких ограничений и основанные на экономном потреблении памяти. Параллельно рассматривается вопрос, нужны ли нам плотные матрицы в итерационных вычислениях или нет, поскольку объём памяти для матриц растёт квадратично с увеличением размера задачи.

Матрица Гессе $nabla^2f(mathbf{x})$ это такая матрица, у которой $ij$ -компонента равна $frac{partial^2 f}{partial x_i partial x_j}$ . Эта матрица использует информацию о второй производной функции (кривизне) и, при правильном использовании, обеспечивает более быструю сходимость за счёт лучшего определения направления шага. Получить такую матрицу для решателя может быть проблематично, в зависимости от типа и условий задачи. Это особенно верно для крупномасштабных задач или сценариев, использующих опосредованный расчёт функции (когда функция не предоставляется нам удобной аналитической формулой). Из матана известно, что если в стационарной точке матрица Гессе положительно определена (SPD), то эта точка является строгим локальным минимизатором. На практике это условие не проверяется даже при наличии матрицы Гессе, ибо слишком затратно проверять знаковую определённость матрицы.

Наконец, для некоторых методов может потребоваться масштабирование (скейлинг). Эффективность метода может зависеть от того, как именно поставлена задача. Задача плохо масштабируется, если изменения для в определённом направлении приводят к гораздо большим изменениям в значении , чем в другом направлении. Простым примером является , которая чувствительна к изменениям $x$ , но не чувствительна к изменениям $y$ . Некоторые методы, применяемые к таким задачам, могут страдать от медленной сходимости. Чтобы ускорить сходимость в этом примере, рекомендуется выполнить масштабирование, изменив переменные на новую пару $(z, y): x to 10^{-4}z, ~y to y$ . Некоторые методы нечувствительны к масштабированию, поэтому они предпочтительнее в практических ситуациях. Увы, на практике о целевой функции мало что бывает известно, и никто нам не предоставляет масштабирующих коэффициентов.

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

Методы линий

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

k := 0

while

рассчитать направление $mathbf{p}_k$
найти такое , что функция $g(alpha_k) := f(mathbf{x}_k + alpha_k mathbf{p}_k)$ минимизируется
$mathbf{x}_{k+1} := mathbf{x}_k + alpha_k mathbf{p}_k$
k := k + 1

Методы этого класса между собой различаются, как правило, тем, как они генерируют направление $mathbf{p}_k$ на шаге 1, и — реже — тем, как они реализуют шаг 2. Тем не менее, что объединяет эти методы при выполнении шага 1 — это то, что все они пытаются рассчитать так называемое направление спуска $mathbf{p}_k$ . Это такое направление, которое обеспечивает мгновенное уменьшение целевой функции, если движется вдоль него. Можно провести аналогию из реальной жизни: человек, идущий по холму, всегда знает, в каком направлении в данном месте идет подъём, а в каком спуск, и выбирает то, которое ему нужно.

Формально, направление спуска $mathbf{p}_k$ определяется так, что выполняется неравенство $mathbf{p}_kcdot nabla f(mathbf{x}_k)<0$ , означающее, что угол между градиентом функции, всегда указывающим в сторону роста, и направлением $mathbf{p}_k$ больше 90 градусов. Указание направления спуска гарантирует корректность одномерной подзадачи минимизации на шаге 2, если предполагать, конечно, что целевая функция начинает расти, когда становится достаточно большим. А это, в свою очередь, гарантируется ограниченностью функции снизу. Не будем же мы искать минимум у неограниченной снизу функции.

Шаг 2 почти никогда не вычисляется точно, несмотря на одномерность задачи. Дело в том, что для этого может потребоваться дорогостоящая оценка матрицы Гессе целевой функции. По этой причине эта подзадача заменяется неточным аналогом, запрашивающим такую позицию $mathbf{x}_k + alpha_k mathbf{p}_k$ , которая обеспечивает в каком-то смысле «хорошее» уменьшение функции. Существует несколько типов условий, которые используются как критерии решения этой неточной задачи, и все они формализованы как те или иные неравенства. Наиболее известны следующие:

Условие Армиджо (или Армихо, фиг знает), или условие достаточного убывания
Условия Вульфа двух типов: стандартное и сильное
Условия Голдштайна

Первое условие заботится о том, чтобы в направлении поиска достигалось достаточно существенное уменьшение целевой функции. Несоблюдение этого условия может привести к стагнации процесса редукции (например, уменьшение может происходить на положительное значение, уменьшающееся от итерации к итерации) во всём алгоритме минимизации. Условия Вульфа в оригинале состоят из двух, первое из которых — само условие Армиджо, а второе называется условием кривизны. Последнее в свою очередь делится на стандартное и сильное условия кривизны и означает, что число достаточно велико, чтобы обеспечивать разумный прогресс в пространстве переменных (простым языком: мы хотим, чтобы сдвигался достаточно далеко, а не топтался на месте).

Сильное условие Вульфа, как очевидно из названия, «усиливает» условие кривизны таким образом, что помещается в окрестность точного минимизатора. Наконец, условия Голдштайна имеют примерно ту же цель, что и условия Вульфа, хотя они могут исключить фактический одномерный минимизатор из области поиска по чистой случайности. Условия представлены в таблице ниже, где мы используем обозначения $phi(alpha) = f(mathbf{x}_k + alpha mathbf{p}_k)$ , $phi'(alpha) = mathbf{p}_kcdot nabla f(mathbf{x}_k + alphamathbf{p}_k)$ , $phi'(0) = mathbf{p}_kcdot nabla f(mathbf{x}_k ) < 0$ .

Ниже изображены четыре условия из таблицы, в указанном порядке: Армиджо, стандартный Вульф, сильный Вульф, Голдштайн. График относится к во всех случаях. Его красные части соответствуют тем значениям $alpha$ , которые удовлетворяют соответствующему условию. обозначает точный одномерный минимизатор: $alpha_* = underset{alpha>0}{text{argmin}}f(mathbf{x}_k + alpha mathbf{p}_k)$ .

Кстати, в машинном обучении коэффициент называют скоростью обучения — learning rate (в самой же оптимизации названия не помню, может что-то вроде «шаг поиска»?).

Условия Вульфа обладают некоторыми доказанными свойствами, связанными с глобальной сходимостью. Приведём ниже очень важную теорему о глобальной сходимости методов линий.

Теорема. Пусть целевая функция ограничена снизу и везде дифференцируема, причём её градиент непрерывен по Липшицу: существует константа такая, что для любой пары . Пусть $mathbf{p}_k$ есть направление спуска, а удовлетворяет стандратному условию Вульфа для всех . Обозначим положительный косинус угла между $mathbf{p}_k$ и как $cos theta_k = -frac{mathbf{p}_kcdotnabla f(mathbf{x}_k)}{|mathbf{p}_k||nabla f(mathbf{x}_k)|} > 0$ . Тогда ряд $sumlimits_{k = 0}^{infty}cos^2 theta_k |nabla f(mathbf{x}_k)|^2$ сходится.
Непрерывность по Липшицу не является слишком ограничивающим условием. Напомним из математического анализа, что это условие гладкости находится между равномерной непрерывностью и дифференцируемостью. Другими словами, если дважды дифференцируема, то липшицева непрерывность градиента гарантируется.

Хорошая новость в том, что сходимость ряда означает сходимость его последовательности к нулю, что также известно из математического анализа. Иными словами, заключаем, что $cos^2 theta_k |nabla f(mathbf{x}_k)|^2to 0$ . Плохая же новость в том, что это само по себе не гарантирует сходимости к стационарной точке $f(mathbf{x}_k)to 0$ из-за присутствия косинуса рядом. Если можно показать, что последовательность косинусов ограничена снизу положительной константой, то сходимость к стационарной точке подразумевается. Но если выяснится, что , т.е. что направление $mathbf{p}_k$ становится всё меньше и меньше направлением спуска из-за стремления к ортогональности градиенту, то факта стремления к стационарной точке может и не быть. Этот фундаментальный недостаток алгоритма может быть одной из основных причин, по которым методы доверительной области имеют такое же право на существование, как и методы линий.

Одномерная минимизация

Второй шаг в алгоритме методов линий выглядит так:

$alpha_k := underset{alpha>0}{text{argmin}} f(mathbf{x}_k + alphamathbf{p}_k)$

Как можно добиться приближённого решения для этой задачи одномерной минимизации?

Более конкретно, нас интересует удовлетворение того или иного условия Вульфа.

Есть относительно простой алогоритм для этих целей. Допустим, нам требуется сильное условие Вульфа. Начиная с некоторого начального приближения $alpha^{(0)}$ , мы каким-нибудь образом шагаем либо вперёд (т.е. последовательность $alpha^{(i)}$ возрастает), либо назад (т.е. убывает, по-английски это называется backtracking). А по ходу шагания, мы проверяем следующие условия в представленном порядке, при учёте невыполнения предыдущих:

Если выполнено второе условие, то цель достигнута. Если выполнено первое или третье, то мы проваливаемся в следующий алгоритм по названием «Зум», позволяющий локализовать нужную нам точку для выполнения условия Вульфа. Это напоминает работу под увеличительным стеклом, отсюда и название.
Если выполнено первое условие, то мы на вход «Зуму» подаём два параметра $alpha^{(i-1)}, alpha^{(i)}$ . Если же третье, то подаём их же, но в другом порядке.

Вот так графически могут выглядеть итерации алгоритма «Зум», где красный участок соответствует выполнению сильного условия Вульфа:

Давайте теперь посмотрим на самые известные алгоритмы безусловной оптимизации метода линий.

Метод наискорейшего спуска

Он же метод градиентного спуска. В алгоритме 1, шаг 1 цикла просто нам выдаёт $mathbf{p}_k = -nabla f(mathbf{x}_k)$ . Да, вот так просто. Для этого случая просто всегда имеем , и глобальная сходимость по теореме выше нам гарантирована. Скорость сходимости линейна, и метод требует масштабирования. Например, для целевой функции с начальным положением $left(sqrt{2}, frac{sqrt{2}}{10}right)$ наискорейший спуск может потребовать десятки итераций, чтобы сойтись к минимизатору , в то время как при масштабировании с коэффициентами $left(1, 10^{-2}right)$ , которые применяются к градиенту, уже одной итерации может быть достаточно, поскольку такое направление указывает непосредственно на начало координат. Метод не требует никакого матричного хранения в памяти.

На рисунке выше показаны первые три итерации, выполненные методом наискорейшего спуска в сторону минимизатора на возможной выпуклой функции, например, выпуклой квадратичной функции. Здесь продемонстированы точные одномерные минимизаторы вдоль направлений $mathbf{p}_k$ , что эквивалентно тому, что вектор антиградиента $-nabla f(mathbf{x}_k)$ касателен линии уровня функции, на которой лежит следующее приближение $mathbf{x}_{k+1}$ . Это также эквивалентно тому, что два последовательных градиента ортогональны друг другу. На практике, т.е. при неточной одномерной минимизации вдоль направления, такое случается крайне редко.

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

Антиградиент $-nabla f(mathbf{x}_0)$ сразу указывает на точный минимизатор! Обратите внимание, что в случае неточного поиска вдоль линий, например, с использованием условий Вульфа, даже этот сценарий может потребовать немного больше одной итерации, чтобы сойтись. Все направления в этом случае будут параллельны друг другу.

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

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

Метод Ньютона

Классический метод Ньютона полагается на решение на каждой итерации системы линейных уравнений $nabla^2 f(mathbf{x}_k)mathbf{p}_k = -nabla f(mathbf{x}_k)$ с матрицей Гессе в левой части. При возможности используется коэффициент в одномерном поиске. Любая матрица Гессе симметрична. Если дополнительно предположить, что она положительно определена (SPD), то созданное направление $mathbf{p}_k$ является наравлением спуска, что несложно показать напрямую:

$cos theta_k = -frac{mathbf{p}_kcdotnabla f(mathbf{x}_k)}{|mathbf{p}_k||nabla f(mathbf{x}_k)|}=frac{nabla f(mathbf{x}_k)cdotleft[left(nabla^2 f(mathbf{x}_k)right)^{-1}nabla f(mathbf{x}_k)right]}{|mathbf{p}_k||nabla f(mathbf{x}_k)|}>0.$

Если известно, что число обусловленности матрицы Гессе равноверно ограничено, т.е. существует константа $M$ такая, что $text{cond}(nabla^2 f(mathbf{x}_k))leqslant M$ для всех $k$ , то $cos theta_k geqslant frac{1}{M} > 0$ , и, cледовательно, глобальная сходимость гарантируется в случае матрицы Гессе из класса SPD. Однако матрица Гессе не является SPD для невыпуклой функции во всей области её определения. В таких случаях используют разнообразные методы для корретикровки собственных значений матрицы Гессе, до тех пор, пока участки положительной определённости в пространстве переменных не будут достигнуты, так чтобы потом итерации могли продолжать свою работу с настоящей матрице Гессе. Такие методы относятся к классу модифицированного метода Ньютона, и затрагивать здесь мы их не будем.

Несмотря на отсутствие глобальной сходимости стандартного метода Ньютона для общего класса функций, преимущество его использования заключается в квадратичной скорости сходимости вблизи минимизатора. Вблизи минимизатора выбор сразу же удовлетворяет условиям Вульфа. Метод также инвариантен по отношению к масштабированию задачи, что является безусловным плюсом. Естественно, необходимость вычислять и сохранять в памяти матрицу Гессе, а также решать линейные системы на каждом шаге, во многих случаях не делает метод Ньютона напрямую доступным на практике.

Откуда шаг в методе Ньютона берётся? Из полиномиальной аппроксимации второго порядка для функции в окрестности точки $mathbf{x}_k$ . Мы просто заменяем функцию её квадратичным многочленом Тейлора и ищем для него минимум, предполагая, конечно, что этот многочлен выпуклый:

$f(mathbf{x}_k + mathbf{p}) approx f(mathbf{x}_k) + nabla f(mathbf{x}_k)cdotmathbf{p} + frac{1}{2}mathbf{p}cdotleft(nabla^2 f(mathbf{x}_k)mathbf{p}right)implies nabla f(mathbf{x}_k) + nabla^2 f(mathbf{x}_k)mathbf{p} = 0.$

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

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

Методы квази-Ньютона

Этот класс методов пытается воспользоваться преимуществами высокой скорости сходимости исходного метода Ньютона, в то же время избегая прямого вычисления матрицы Гессе, а иногда даже решения линейной системы, показанной выше. Вычисление направления в этом случае может быть условно описано как $mathbf{p}_k = -mathbf{H}_knabla f(mathbf{x}_k)$ с некоторой приближённой обратной матрицей Гессе.

В методах семейства Бройдена используются некоторые эвристические условия (такие как, например, уравнение секущей) для аппроксимации фактической матрицы Гессе или её обратной таким способом, который желателен для сходимости. Наиболее известные методы из этого семейства — DFP, BFGS, SR1, при этом BFGS является среди них признанным лидером. Построение обратной матрицы Гессе в соответствии с алгоритмом BFGS сохраняет симметрию и положительную определённость матрицы при условии выполнения условий Вульфа на каждой итерации. Последнее требование важно, так как эти условия сильно привязаны к построению матрицы. Следовательно, BFGS гарантирует направление спуска. Однако мало что известно об ограниченности числа обусловленности построенной матрицы для общего класса целевых функций, и поэтому известные теоретические результаты о глобальной сходимости довольно узкие (в частности, они требуют выпуклости функции). Тем не менее, BFGS оказался очень надёжным методом на практике. Кроме того, он обладает сверхлинейной сходимостью, хотя квази-ньютоновские методы в этом отношении все ещё уступают чистому методу Ньютона.

Общее преимущество методов семейства Бройдена заключается в отсутствии необходимости вычислять фактическую матрицу Гессе. Матрица строится исключительно на основе информации о градиенте на предыдущем и текущем шагах. Метод не требует масштабирования, как и метод Ньютона.

SR1 не гарантирует, что его приближённые матрицы Гессе являются положительно определёнными, но, похоже, он хорошо работает в рамках методов доверительной области, где приближённые неопределённые матрицы Гессе могут быть даже предпочтительнее. Мы вернёмся к технике SR1 в случае, если я напишу пост о методах доверительной области.

Покажем наконец формулы для обновления приближённой обратной матрицы Гессе в соответствии с методом BFGS:

$mathbf{s}_k = mathbf{x}_{k+1} - mathbf{x}_k,quad mathbf{y}_k = nabla f(mathbf{x}_{k+1}) - nabla f(mathbf{x}_{k}),\ rho_k = frac{1}{mathbf{s}_kcdotmathbf{y}_k},\ mathbf{H}_{k+1} = (mathbf{I} - rho_kmathbf{s}_kmathbf{y}_k^T)mathbf{H}_k(mathbf{I} - rho_kmathbf{s}_kmathbf{y}_k^T)^T +rho_kmathbf{s}_kmathbf{s}_k^T.$

Обратите внимание на то, что член $mathbf{s}_kcdotmathbf{y}_k>0$ всегда положителен при выполнении условия кривизны в неравенствах Вульфа.

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

$f''(x_k)(x_{k+1} - x_{k}) = -f'(x_k).$

Обозначим . Мы, в результате, ищем корень функции $g(x)$ , т.е. решение уравнения , методом Ньютона

$g'(x_k)(x_{k+1} - x_{k}) = -g(x_k).$

Метод секущих заменяет производную тангенсом угла секущей, построенной на двух ближайших приближениях:

$g'(x_k)approxfrac{g(x_k) - g(x_{k-1})}{x_k - x_{k-1}}.$

См. рисунок с итерациями на графике $g(x)$ :

Это и порождает квази-ньютоновский алгоритм оптимизации (правда, в данном случае, с фиксированным шагом поиска ):

$frac{f'(x_k) - f'(x_{k-1})}{x_k - x_{k-1}}(x_{k+1}-x_k)=-f'(x_k).$

Приглушённый BFGS

В случае, если условия Вульфа не выполняются на каком-то шаге $k$ , мы не можем быть уверены, что наша следующая матрица Гессе, полученная с помощью формулы BFGS, будет положительно определённой. Это может нарушить глобальную сходимость метода линий. Чтобы избежать этой проблемы, можно расширить формулу BFGS на более общий случай, как показано ниже.

Обозначим $u_k = mathbf{s}_kcdotmathbf{y}_k,quad v_k = mathbf{s}_kcdot(mathbf{B}_kmathbf{s}_k)$ , где $mathbf{B}_k$ означает приближённую SPD-матрицу Гессе на итерации $k$ . Член $u_k = mathbf{s}_kcdotmathbf{y}_k$ теперь может иметь любой знак, в то время как член $v_k = mathbf{s}_kcdot(mathbf{B}_kmathbf{s}_k)$ строго положителен. Рассмотрим

$gamma_k = begin{cases}1, &text{если } u_kgeqslant frac{1}{5}v_k,\ frac{4}{5}frac{v_k}{v_k - u_k}, &text{в противном случае,}end{cases}$

и $mathbf{r}_k = gamma_kmathbf{y}_k + (1-gamma_k)mathbf{B}_kmathbf{s}_k$ .
Тогда формула обновления приближённой матрицы Гессе (не её обратной, как в случае с оригинальным BFGS) следующая:

$mathbf{B}_{k+1} = mathbf{B}_k - frac{1}{v_k}mathbf{B}_kmathbf{s}_kmathbf{s}_k^Tmathbf{B}_k+frac{1}{u_k}mathbf{r}_kmathbf{r}_k^T.$

Эта формула, называемая формулой приглушённого BFGS (damped BFGS), обеспечивает положительную определённость матрицы $mathbf{B}_{k+1}$ . Релевантность подобной аппроксимации — тема для дебатов, но, по крайней мере, направление спуска будет сохранено (Подобные трюкачества с изменением матрицы Гессе, искажающие аппроксимацию, происходят и в модифицированном методе Ньютона ради достижения направления спуска, см. выше ). Эту формулу можно пытаться использовать для поиска минимизатора самой общей невыпуклой целевой функции в рамках метода линий.

Метод нелинейных сопряжённых градиентов (NCG)

Методы класса NCG в некотором смысле уникальны. Они возникли из прямой адаптации традиционного метода сопряжённых градиентов (CG) для решения линейных SPD-систем к нелинейной оптимизации.

Чуть подробнее

Имеется в виду вот что. Допустим, требуется решить линейную систему с положительно определённой матрицей $A$ . Метод CG основан на минимизации квадратичной выпуклой функции

$frac{1}{2}mathbf{x}cdot(Amathbf{x}) - mathbf{b}cdotmathbf{x},$

единственная точка минимума которой совпадает с решением оригинальной линейной системы. В детали итераций CG вдаваться не будем, кроме одной: каждая итерация полностью определяется вектором невязки линейной системы $Amathbf{x}_k - mathbf{b}$ , которая есть ничто иное как градиент квадратичной функции в текущем приближении $mathbf{x}_k$ ! Так вот, в далёком 1964 году чуваки по имени Fletcher и Reeves взяли и поставили в итерации на место градиента квадратичной функции градиент самой общей нелинейной функции. Это и дало начало семейству методов NCG.

Методы NCG являются единственными в этом списке, которые строят текущее направление на основе предыдущего: $mathbf{p}_k = -nabla f(mathbf{x}_k) + beta_k mathbf{p}_{k-1}$ . Все методы в пределах этого класса в основном различаются тем, как задаётся скалярный параметр . Они так или иначе используют информацию о текущем и предыдущем градиентах при вычислении . Идея NCG принадлежит авторам Fletcher-Reeves (FR). Последующие модификации предложены следующими авторами: Polac-Rabier (PR), Hestenes-Stiefel (HS), Dai-Yuan (DY), Hager-Zhang (HZ). Если также модификации с обозначениями PR+, FR-PR.

Результаты по глобальной сходимости весьма неопределённы. Для FR наилучший результат по сходимости заключается в том, что последовательность магнитуд $|f(mathbf{x}_k)|$ не отделена от нуля (т.е. нижний предел последовательности $|f(mathbf{x}_k)|$ равен нулю), но ничего не говорится о том, существует ли сам предел, равный нулю. Что касается PR, даже если на практике он работает лучше, чем FR, теоретические результаты по нему ещё хуже, поскольку для него был построен расходящийся пример с специальной, но достаточно гладкой целевой функцией. FR гарантирует направление спуска до тех пор, пока выполняются строгие условия Вульфа с параметрами $0 < c_1 < c_2 < frac{1}{2}$ . Но с практической точки зрения FR является наихудшим методом, так как он может стагнировать из-за почти прямого угла между градиентом и направлением спуска и никогда не выйти из этой стагнации. Его можно улучшить, выполняя перезапуски (restarts), т.е. полагая всякий раз, когда обнаруживается близкий к прямому угол. Это заменяет направление FR на направление по антиградиенту.

Методы PR+ и HS гарантируют направление спуска, но требуют выполнения так называемого условия достаточного спуска во время одномерного поиска. Для DY и HZ сильных условий Вульфа достаточно для генерации направления спуска на следующей итерации, и эти два метода по сей день остаются одними из самых надежных.

Методы NCG часто рекомендуют использовать только для крупномасштабных задач ввиду их дешевизны, поскольку они менее устойчивы по сравнению с другими методами линий. Они не требуют хранения матрицы в памяти, и считается, что они в целом сходятся быстрее, чем метод наискорейшего спуска. Скорость сходимости NCG обычно линейна, но в некоторых работах была показана квадратичная или даже так называемая суперквадратичная $n$ -сходимость. Некоторое масштабирование для NCG может быть предпочтительным.

Ниже приведём формулы для для упомянутых методов из класса NCG:

В последней формуле $mathbf{y}_k = nabla f(mathbf{x}_{k+1}) - nabla f(mathbf{x}_{k})$ , как и для формулы BFGS.

Усечённый метод Ньютона (truncated Newton method)

Для крупных задач, даже если матрица Гессе у нас в наличии и используются методы для решения разреженных систем линейных уравнений, часто бывает слишком сложно вычислить шаг Ньютона $nabla^2 f(mathbf{x}_k)mathbf{p}_k = -nabla f(mathbf{x}_k)$ «в лоб». По этой причине на сцену выходит итерационный линейный решатель, такой как метод сопряженных градиентов (CG, см. упоминание выше). Это даёт начало семейству так называемых неточных, или усечённых, методов Ньютона (не путать с модифицированным). Критерием останова CG в этом случае будет

$|nabla^2 f(mathbf{x}_k)mathbf{p}_k + nabla f(mathbf{x}_k)| < epsilon$

с некоторой точностью .

Этот метод имеет сверхлинейную скорость сходимости при соблюдении определённых условий. Поскольку решатель CG сходится только для матриц SPD, неточный метод Ньютона адаптирует CG так, чтобы тот обрывал итерации сразу после обнаружения отрицательной кривизны (т.е. попадает в зону невыпуклости функции), и затем использует полученное, «обрезанное» направление спуска.

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

$nabla^2f(mathbf{x})mathbf{p}approx frac{1}{h}(nabla f(mathbf{x} + hmathbf{p}) - nabla f(mathbf{x})).$

Недостатком усечённого метода Ньютона является его плохая сходимость для почти сингулярных или плохо обусловленных матриц Гессе (что верно и для обычного метода Ньютона). Жетально ли для усеченного метода Ньютона масштабирование, я до конца не уверен.

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

BFGS c ограниченной памятью

BFGS c ограниченной памятью, также известный как L-BFGS, является адаптацией стандартного метода BFGS для задач с большим количеством переменных. Проблема стандартного BFGS, как и остальных квази-ньютоновских методов, в том, что обновление матрицы Гессе или её обратной на каждом шаге использует внешнее произведение плотных векторов. Имеется в виду произведение вертикального вектора на горизонтальный, например, вот так: $mathbf{s}_kmathbf{s}_k^T$ (см. формулу BFGS выше). Такое произведение генерирует плотную матрицу. По этой причине, из-за быстрого расхода памяти с ростом размерности задачи, была разработана экономичная версия BFGS.

В основе этой версии лежит запись последовательности разностей градиентов $mathbf{y}_k = nabla f(mathbf{x}_{k+1}) - nabla f(mathbf{x}_{k})$ и координат $mathbf{s}_k = mathbf{x}_{k+1} - mathbf{x}_{k}$ из последних $m$ итераций, где $m$ обычно принимает значение от 3 до 20. В процессе расчёта квази-ньютоновского направления $mathbf{p}_k = -mathbf{H}_knabla f(mathbf{x}_k)$ , обновление матрицы заменяется конкретной рекурсивной формулой произведения матрицы на вектор градиента с использованием заданной последовательности на каждом шаге. После этого последовательность обновляется путём удаления самой старой пары ${mathbf{s}_{k-m}, mathbf{y}_{k-m}}$ и добавления в очередь новой пары. По сути мы сохраняем память, немного жертвуя производительностью за счёт прикручивания двух циклов на каждой итерации метода.
L-BFGS точно воспроизводит процедуру BFGS, как если бы она применялась для последних $m$ шагов, начиная с заданного начального матричного приближения. Это приближение предоставляется специальным образом.

Вот скетч алгоритма L-BFGS для вычисления $H_knabla f(mathbf{x}_k)$ .

В алгоритме выше член $rho_i$ означает $frac{1}{mathbf{s}_icdotmathbf{y}_i}>0$ .

Несмотря на то, что L-BFGS плохо себя ведёт для плохо обусловленных задач (с широким распределением собственных значений в истинной матрице Гессе), он часто является удовлетворительным методом для задач безусловной минимизации с большим количеством переменных. Существуют и другие квази-ньютоновские методы с ограниченной памятью, но они мне кажутся менее популярными, чем L-BFGS.

Итог

Подытожим. В таблице ниже представлены описанные выше методы с точки зрения 6 характеристик, которые для нас представляли интерес.

Оговорюсь, не в 100% пунктов выше я уверен, например, что касается масштабирования, но в целом, думаю, таблица достаточно информативная.

Примеры

Рассмотрим два кейса для минимизации в $mathbb{R}^2$ :

В втором кейсе есть участки невыпуклости функции, что должно создать дополнительные трудности для оптимизатора. Точка в этом примере является стационарной (более конкретно — седловой), поэтому наше начальное приближение немного сдвинуто от неё, однако оно расположено в зоне невыпуклости.

Критерием остановки мы в обоих случаях полагаем $|nabla f(mathbf{x})|<10^{-4}$ .
Мы посмотрим, как сработают три алгоритма метода линий: наискорейший спуск, BFGS, FR. Последний мы используем с перезапуском: при нарушении условия , где есть угол между $-nabla f(mathbf{x}_k)$ и $-nabla f(mathbf{x}_k) + beta_kmathbf{p}_{k-1}$ , мы полагаем $mathbf{p}_k := -nabla f(mathbf{x}_k)$ .
Вот такие параметры мы использовали для одномерной минимизации:

Ниже покажем результаты. На картинках изображены линии уровня целевой функции, а также — самое интересное — итерации каждого из трёх методов.

Первый кейс

Наискорейший спуск, 9 итераций:

BFGS, 7 итераций (ожидаемо чуть увереннее, чем наискорейший спуск, на выпуклой функции):

FR, 14 итераций:

Второй кейс

Наискорейший спуск, 14 итераций:

BFGS, 15 итераций:

FR, 25 итераций:

В обоих случаях, FR потребовал больше итераций, чем даже наискорейший спуск, хотя необходимость перезапуска возникла только раз — на третьей итерации в втором кейсе. Также ни в одном из двух кейсов BFGS не потребовал приглушения, т.к. нам удалось на каждой итерации удовлетворить условию Вульфа.

В будущем, если ничто не помешает, хотел бы продолжить эту тему. Следующий кандидат на очереди — метод доверительной области (trust-region method).

А на горизонте маячит условная оптимизация.

Всем добра!

Список литературы

(Разумеется, список далеко не полный)

[1] Аттетков, А. В., Галкин, С. В., & Зарубин, В. С. (2001). Методы оптимизации. Изд. МГТУ им. Баумана.

[2] Nocedal, J., Wright, S. J., Numerical Optimization. Second Edition. Springer 2006.

[3] More, J. J., Newton’s method. Studies in Numerical Analysis, MAA Studies in Mathematics, 24 (1984): 29–82.

[4] Gill P.E., Murray W., Wright M. H., Practical Optimization, Academic Press, 1981.

[5] Powell M. J. D., Some global convergence properties of a variable metric algorithm for minimization without exact line searches, Nonlinear Programming, SIAM-AMS Proceedings, R. W. Cottle and C. E. Lemke, eds., SIAM Publications, IX(1976): 53–72.

[6] Fletcher R., Reeves C. M., Function minimization by conjugate gradients, Computer Journal, 7 (1964): 149–154.

[7] Dai Y., Yuan Y., A nonlinear conjugate gradient method with a strong global convergence property, SIAM Journal on Optimization, 10 (1999): 177–182.

[8] Hager W. W., Zhang H., A new conjugate gradient method with guaranteed descent and an efficient line search, SIAM Journal on Optimization, 16 (2005): 170–192.

[9] Al-Baali M., Descent property and global convergence of the Fletcher-Reeves method with inexact
line search, I.M.A. Journal on Numerical Analysis, 5 (1985): 121–124.

[10] Powell M. J. D., Nonconvex minimization calculations and the conjugate gradient method, Lecture Notes in Mathematics, 1066 (1984): 122–141.

[11] Gilbert J., Nocedal J., Global convergence properties of conjugate gradient methods for optimization, SIAM Journal on Optimization, 2 (1992): 21–42.

Источник

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

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

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

В настоящее время разработано множество численных методов для задач как безусловной, так и условной оптимизации. Качество численного метода характеризуется многими факторами: скоростью сходимости, временем выполнения одной итерации, объемом памяти ЭВМ, необходимым для реализации метода, классом решаемых задач и т. д. Решаемые задачи также весьма разнообразны: они могут иметь высокую и малую размерность, быть унимодальными и многоэкстремальными и т. д. Один и тот же метод, эффективный для решения задач одного типа, может оказаться совершенно неприемлемым для задач другого типа.

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

пример, в [12,14].

Начнем с рассмотрения прямых методов безусловной оптимизации, когда ограничения отсутствуют.

Смысл прямых методов безусловной оптимизации состоит в построении последовательности точек X[0], X[1], …, X[N], таких,

что f(X[0])>f(X[1])>… …>f(X[n]). В качестве начальной точки X[0] может быть выбрана произвольная точка, однако стремятся ее выбрать как можно ближе к точке минимума. Переход (итерация) от точки Х[k] к точке Х[k+1], k=0,1,2,… состоит из двух этапов:

–выбор направления движения из точки Х[k];

–определение шага вдоль этого направления.

Методы построения таких последовательностей часто называют методами спуска, так как осуществляется переход от бо́льших значений функции к меньшим.

Математически методы спуска описываются соотношением

X[k+1]=X[k]+ak p[k], k=0,1,2,…,

где p[k] – единичный вектор, определяющий направление спуска;

ak – длина шага.

Различные методы спуска отличаются друг от друга способами выбора p[k] и ak . На практике применяются только методы, обладающие сходимостью. Они позволяют за конечное число шагов получить точку минимума или подойти к ней достаточно близко. Качество сходящихся итерационных методов оценивают по скорости сходимости.

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

щения	аргумента		X[k]− X[k −1]		< ε	или	функции

f (X[k])− f (X[k −1]) < γ . Здесь k – номер итерации; ε, γ – задан-

ные величины точности решения задачи.

Методы поиска точки минимума называются детерминированными, если оба параметра перехода от X[k] к X[k+1] (направление движения и величина шага) выбираются однозначно по доступной в точке X[k] информации. Если же при переходе используется какой-либо случайный механизм, то алгоритм поиска называется случайным поиском минимума.

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

при необходимости дополнительного вычисления вторых производных – методы второго порядка.

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

Методы одномерного поиска. Перечень методов одномерного поиска – численного поиска экстремума функции одного аргумента f(x) – достаточно широк и хорошо освещен в литературе [12]. Поэтому здесь ограничимся рассмотрением только одного метода, который, по опыту авторов, является одним из наиболее эффективных, – метода «золотого сечения».

Идея метода состоит в последовательном сокращении интервала неопределенности – интервала значений аргумента x, содержащего искомую точку минимума, – до длины, не превышающей

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

Внутри любого интервала [a;b] содержатся две точки x=y0 и x=z0 , выполняющие его «золотое сечение» – разбиение на две неравные части такие, что отношение большей части к длине всего интервала совпадает с отношением меньшей части к большей. Очевидно, эти точки расположены симметрично относительно центра интервала (рис. 26). Координаты точек «золотого сечения» могут быть найдены из соответствующих пропорций:

b − y0	=	y0 − a	= δ ,	z0 − a	=	b − z0	= δ,

b − a		b − y		b − a		z	0	− a
		0

откуда нетрудно получить δ=(1–δ)/δ и прийти к уравнению: δ2+δ–1=0. В результате получим относительные доли, определяющие «золотое сечение» интервала: δ=0,618, 1–δ=0,382. «Золотое сечение» обладает важным свойством: точка y0 является одной из точек «золотого сечения» интервала [a;z0 ], точка z0 – одной из точек «золотого сечения» интервала [y0 ;b]. В этом убе-

ждает простой расчет: 0,382/0,618 = 0,618 и (0,618–0,382)/0,618 = = 0,382.

Алгоритм поиска минимума, построенный на основе метода «золотого сечения», предусматривает на каждой итерации выбор в качестве одной из границ сокращенного интервала левой или правой точки «золотого сечения» таким образом, чтобы искомый минимум сохранялся внутри него:

1.Задают k=0, исходный интервал неопределенности [ak;bk], допустимую погрешность результата ε.

2.Вычисляют координаты точек «золотого сечения»:

yk=ak+0,382(bk–ak), zk=ak+0,618(bk–ak).

3. Вычисляют значения целевой функции в найденных точках

f(yk) и f(zk).

4. Если f(yk)≤f(zk) (рис. 26, а), присваивают ak + 1 =ak , bk + 1 =zk , zk + 1 =yk , yk + 1 =ak +zk –yk , k=k+1. В противном случае (рис. 26, б) ak + 1 =yk , bk + 1 =bk , yk + 1 =zk , zk + 1 =yk +bk –zk , k=k+1.

5. Проверяют выполнение условия завершения поиска

bk+1 −ak+1 ≤ ε . В случае его выполнения в качестве решения выбирают точку x = (yk+1 + zk+1)2 . В противном случае переходят к шагу 2.

Рис. 26

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

Далее обратимся к методам многомерного поиска.

Метод прямого поиска (метод Хука-Дживса). Задаются неко-

торой начальной точкой Х[0]. Поочередно изменяя компоненты вектора Х, обследуют окрестность данной точки, в результате чего находят точку (новый базис), определяющую направление, в котором происходит уменьшение минимизируемой функции f(Х). В выбранном направлении осуществляют спуск, убеждаясь, что значение функции уменьшается. Процедура циклически повторяется, пока удается находить направление спуска с учетом принятого условия останова.

Алгоритм метода прямого поиска в самом общем виде можно сформулировать следующим образом:

1. Задаются значениями координат хi [0], i=1,2,…n, начальной точки (k=0), вектором начальных приращений координат

∆X= (∆х1 , ∆ х2 ,…, ∆ хn ) в процессе обследования окрестности, наименьшим допустимым значением ε компонент ∆X, ускоряющим множителем λ ≥1, определяющим скорость спуска, масштабным коэффициентом d>1.

2. Принимают Х[k] за «старый базис» [12]: Xб=Х[k]. Вычисляют

значение f(Xб).

3. Поочередно изменяют каждую координату хбi, i=1,2,…n,

точки Xб на величину ∆хi , то есть принимают хi [k]=хбi+∆хi , затем

хi [k]=хбi–∆хi . Вычисляют значения f(X[k]) в получаемых пробных точках и сравнивают их со значением f(Xб). Если f(X[k])< < f(Xб), то соответствующая координата хi [k] приобретает новое значение, вычисленное по одному из приведенных выражений. В противном случае значение этой координаты остается неизменным. Если после изменения последней n-й координаты f(X[k])<f(Xб), принимают X[k] за «новый базис» и переходят к п. 4. В противном случае – к п. 6.

4. Осуществляют спуск в направлении от «старого» к «новому» базису через последний, т. е. вычисляют координаты новой точки

X[k+1]: хi [k+1]=хi [k]+λ(хi [k]–хбi), i=1,2,…n. Вычисляют значение f(X[k+1]). Если выполняется условие f(X[k+1])<f(X[k]),

«новый» базис принимают «старым» (Xб=Х[k], f(Xб)=f(X[k])) и переходят к п. 5. В противном случае принимают хi [k+1]=хi [k], i=1,2,…n.

5.Как и в п. 3, поочередно изменяют каждую координату точки X[k+1], сравнивая соответствующие значения функции f(Х) со значением f(X[k+1]), полученным в п. 4. После изменения последней координаты сравнивают соответствующее значение

функции f(X[k+1]) со значением f(Xб), полученным в п. 4. Если f(X[k+1])<f(Xб), полученную точку X[k+1] принимают за «новый базис», увеличивают значение k на единицу и переходят к п. 4, в противном случае – к п. 3.

6.Если для всех i ∆хi <ε, вычисления прекращаются. В противном случае уменьшают значения ∆хi в d раз и переходят к п. 3.

Работа алгоритма проиллюстрирована рис. 27. Показаны линии

уровня минимизируемой функции f(x1 ,x2 ), т. е. линии, получаемые из условий f(x1 ,x2 )=f1 =const, f(x1 ,x2 )=f2 =const и так далее. Здесь f1 >f2 >f3 . Сплошные линии – результаты однократного выполнения пп. 3…5 (поиск направления уменьшения функции и спуск), пунктирная линия – следующий спуск.

Рис. 27

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

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

Метод деформируемого многогранника (метод Нелдера— Мида) состоит в том, что для минимизации функции n переменных f(X) в n-мерном пространстве строится многогранник, содержащий n+1 вершину. Очевидно, что каждая вершина соответствует некоторому вектору Xi. Вычисляют значения целевой функции f(Xi), i=1,2,…, n+1, в каждой из вершин многогранника, определяют максимальное из этих значений и соответствующую ему вершину Xh. Через эту вершину и центр тяжести остальных вершин проводят проецирующую прямую, на которой находится точка Xq с меньшим значением целевой функции, чем в вершине Xh (рис. 28, а). Затем исключают вершину Xh. Из оставшихся вершин и точки Xq строят новый многогранник, с которым повторяют описанную процедуру. В процессе выполнения таких операций многогранник изменяет свои размеры, что и обусловило название метода.

Рис. 28

Введем следующие обозначения: X[i,k] – вектор координат i-й вершины многогранника на k-м шаге поиска, i=1,2,…n+1, k=1,2,…; h – номер вершины, в которой значение целевой

функции максимально: f (X[h,k])= max f (X[i,k]); l	– номер вер-
i
шины, в которой значение целевой	функции	минимально:
f (X[h,k])= min f (X[i,k]); X[n+2,k] –	центр тяжести всех вер-
i

шин, за исключением X[h,k]. Координаты центра тяжести вычис-

	xj [n + 2,k]=	1	n+1
ляют по формуле	∑xj [i,k]− xj [h,k]	, j=1,2,…n.
		n
			j=1

Примерный алгоритм метода деформируемого многогранника состоит в следующем:

1.Задаются коэффициентами отражения α, растяжения γ>1, сжатия β<1, допустимой погрешностью определения координат

точки минимума ε. Выбирают координаты вершин исходного многогранника X[i,k], i=1,2,…n+1, k=1.

2.Вычисляют значения целевой функции во всех вершинах f(X[i,k]), i=1,2,…n+1, и находят точки X[h,k], X[l,k] (на рис. 28, б точки соответственно X2 и X1), а также X[n+2,k].

3.Осуществляют проецирование точки X[h,k] через центр тя-

жести: X[n+3,k]=X[n+2,k]+α(X[n+2,k]–X[h,k]).

4.Если f(X[n+3,k])≤ X[l,k], выполняют операцию растяже-

ния: X[n+4,k]=X[n+2,k]+γ(X[n+3,k]–X[n+2,k]). В противном случае переходят к п. 6.

5.Строят новый многогранник: если f(X[n+4,k])<X[l,k], то

заменой X[h,k] на X[n+4,k], в противном случае – заменой X[h,k] на X[n+3,k]. Продолжают вычисления с п. 2 при k=k+1.

6. Если X[h,k]>f(X[n+3,k])>X[i,k] для всех i, не равных h,

выполняют операцию сжатия: X[n+5,k]=X[n+2,k]+β(X[h,k]– X[n+2,k]). Строят новый многогранник заменой X[h,k] на X[n+5,k] и продолжают вычисления с п. 2 при k=k+1.

7. Если f(X[n+3,k])>X[h,k], то, сохраняя вершину X[l,k], строят новый многогранник, подобный текущему, уменьшением длин всех ребер в два раза: X[i,k]=X[l,k]+0,5(X[i,k]–X[l,k]) и продолжают вычисления с п. 2 при k=k+1.

В пп. 6, 7 перед переходом к п. 2 необходима проверка выполнения условия завершения поиска минимума, например, по усло-

вию max n∑+1(xj [i,k]− xj [n + 2,k])2 < ε2 .

ij=1

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

Большинство известных рекомендаций сводятся к выбору

α=1, 2≤ γ≤3, 0,4≤β≤0,6.

Метод вращающихся координат (метод Розенброка). Его суть состоит в последовательных поворотах системы координат в соответствии с изменением направления наиболее быстрого убывания целевой функции (рис. 29). Из начальной точки X[0] осуществляют спуск в точку X[1] по направлениям, параллельным координатным осям. На следующей итерации одна из осей должна проходить в направлении x’1=X[1]–X[0], остальные – в направлениях, перпендикулярных к x’1. Спуск вдоль этих осей пр и- водит в точку X[2], что дает возможность построить новый вектор x’’1=X[2]–X[1] и на его базе новую систему направлений поиска

точки минимума X .

Рис. 29

100

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

Вчастности, данный метод в отличие от многих других эффективен при минимизации так называемых «овражных» функций (с сильно вытянутыми поверхностями уровня), так как результирующее направление поиска стремится расположиться вдоль оси «оврага».

Метод параллельных касательных (метод Пауэлла). Его суть состоит в последовательном проведении одномерного поиска минимума целевой функции по n+1 направлению каким-либо из известных одномерных методов. На первой итерации в качестве первых n направлений выбираются координатные, в качестве (n+1)-го направления используется первое из них (рис. 30). На каждой последующей итерации поиск начинается со второго направления предшествующей итерации, соответственно номера направлений уменьшаются на единицу; (n+1)-е направление последующей итерации задается вектором X[1]–X[n+1] – из точки минимума, найденной на первом шаге предшествующей итерации, через точку минимума, найденную на последнем ее шаге.

Рис. 30

101

Источник

Безусловная оптимизация

Cтраница 2

Все рассмотренные ранее методы являются методами безусловной оптимизации.
[16]

С помощью естественного обобщения получим задачу безусловной оптимизации, охватывающую весь блок ВПБ.
[17]

Среди методов поиска локального экстремума методы безусловной оптимизации составляют наиболее многочисленную группу.
[18]

Как правило, при решении задач безусловной оптимизации для достаточно гладких выпуклых или вогнутых функций F () методы, использующие первые и вторые производные, сходятся быстрее, чем методы нулевого порядка. Однако при оптимизации схем в условиях сложного непредсказуемого рельефа целевых функций F ( X), их алгоритмического, неявного способа задания, например посредством решения системы дифференциальных уравнений, а также при сложной форме ограничений использование методов нулевого порядка часто предпочтительнее. Кроме того, при неявном задании / 7 ( Х) ее производные приходится определять численно, а возникающие при этом ошибки, особенно в окрестности экстремума, создают значительные трудности для точного определения точки оптимума.
[19]

Для нейронных сетей хорошо работают многие методы безусловной оптимизации, часто лучше, чем узкоспециальные, придуманные для обучения нейронных сетей.
[20]

В этом случае мы говорим о задаче безусловной оптимизации.
[21]

Квазиньютоновские методы являются эффективным средством решения задач безусловной оптимизации. Их отличает высокая скорость сходимости, в то же время при реализации квазиньютоновских алгоритмов не приходится выполнять такие трудоемкие операции, как вычисление матрицы вторых производных или обращение матрицы. Однако при большой размерности пространства необходимость хранения и пересчета на каждом шаге матриц Hk обусловливает высокие требования к объему занимаемой памяти ЭВМ. Этот недостаток не присущ изучаемому в следующем параграфе методу сопряженных градиентов.
[22]

Задача поиска Е ( W) является задачей безусловной оптимизации.
[23]

На рис. 12.12 приведены наиболее часто используемые методы безусловной оптимизации.
[24]

В настоящем разделе рассматриваются аспекты практического применения методов безусловной оптимизации к обучению нейросетевых моделей.
[26]

При решении задач условной оптимизации целесообразно использовать методы безусловной оптимизации, учитывая большое количество разработанных по этим методам программ. С этой целью задача условной оптимизации сводится к задаче безусловной оптимизации устранением ограничений путем преобразования параметра xt, на значения которого наложены ограничения, в неограничиваемый.
[27]

Теория необходимых и достаточных условий оптимальности в задачах безусловной оптимизации излагается в любом курсе математического анализа.
[28]

Там же были приведены соответствующие результаты для задачи безусловной оптимизации и классической задачи на условный экстремум. В данной главе излагается общая теория необходимых и достаточных условий оптимальности в задачах математического программирования, включающая указанные результаты как частные случаи.
[29]

Рассмотрим один из вариантов симплексного метода решения задачи безусловной оптимизации.
[30]

Страницы:

Источник