Методы оптимизации

Назад § 5. Обзор методов безусловной оптимизации Вперед

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

О. Генри. Джефф Питерс как персональный магнит


А спускаясь по темной лестнице, зажигайте ... хотя бы спичку.

А.С. Грин. Зеленая лампа

В этом параграфе описываются модификации градиентного метода и метода Ньютона, а также ряд методов, основанных на других идеях.

5.1. Метод тяжелого шарика.

Градиентный метод решения задачи безусловной минимизации

f(x) ® min,(1)

где f: Rm ® R, можно интерпретировать в терминах обыкновенных дифференциальных уравнений следующим образом. Рассмотрим дифференциальное уравнение

px·+ f ў(x) = 0 (2)

(здесь точка над x обозначает производную по независимой переменной t, а f ў(x) как обычно обозначает градиент отображения f: Rm ® R; предполагается, что p > 0). Простейший разностный аналог уравнения (2), а именно, явная схема Эйлера

pxn+1 - xn
h

 + f ў(xn) = 0

и есть градиентный метод для задачи (1):


 xn+1 = xn -

h
p

f ў(xn). 

(3)

Рассмотрим теперь вместо уравнения (2) уравнение

mx··+ px·+ f ў(x) = 0,

описывающее движение шарика массы m в потенциальном поле f ў при наличии силы трения. Потери энергии на трение вынудят шарик спуститься в точку минимума потенциала f, а силы инерции не дадут ему осциллировать так, как это изображено на рис. 8. Это позволяет надеяться, что изменение уравнения (2) введением в него инерционного члена mx··улучшит сходимость градиентного метода (3). Конечно-разностный аналог уравнения, описыавющего движение шарика — это, например, уравнение

mxn+1 - 2xn + xn-1
h2
 + p xn - xn-1
h

 + f ў(xn) = 0.

После простых преобразований и очевидных обозначений мы получаем

xn+1 = xn - af ў(xn) + b(xn - xn-1).(4)

Итерационная формула (4) задает метод тяжелого шарика решения задачи безусловной оптимизации (см. рис. 14; ср. с рис. 8).

Метод тяжелого шарика
Рис. 14.

Можно доказать, что в условиях теоремы 3.7 метод тяжелого шарика при a = 2/(ЦL + Цl)2 и b = (ЦL - Цl)/(ЦL +Цl)2 сходится со скоростью геометрической прогрессии со знаменателем q = (ЦL -Цl)/(ЦL + Цl).

Если теперь сравнить знаменатели qгм = (L - l)/(L + l) и qмтш = (ЦL - Цl)/(ЦL + Цl), характеризующие скорости сходимости градиентного метода и метода тяжелого шарика, соответственно, то для плохо обусловленных функций, т. е. для функций с m = L /l >> 1, очевидно, qгм » 1- 2/m, а qмтш » 1 - 2/Цm. Поэтому для уменьшения погрешности в e » 2.718 раз градиентный метод с постоянным оптимальным шагом требует -[ln(1 - 2/m)]-1 » m)/2 итераций, а метод тяжелого шарика -ln(1 - 2/Цm)]-1 » Цm/2. Для больших m это весьма значительный выигрыш, поскольку объем вычислений в методе тяжелого шарика почти не отличается от объема вычислений в градиентном методе.

5.2. Методы переменной метрики.

Введем в пространстве Rm новое скалярное произведение б·,·с формулой бx, yс = (Sx, y), где S самосопряженный положительно определенный оператор на Rm. Оно естественно порождает новую норму |||x||| = бx, yс1/2 и метрику на Rm. Операция взятия градиента дифференцируемой функции, как легко видеть, не инвариантна относительно скалярного произведения (метрики): градиент бf сў функции f относительно нового скалярного произведения связан со старым градиентом f ў соотношением

бf сў(x) = S-1f ў(x).

З а д а ч а  5.1. Докажите.

Соответственно градиентный метод в новой метрике имеет вид

xn+1 = xn - aS-1f ў(xn).(5)

Естественно желание подбирать оператор S с целью ускорения сходимости. Например, если f квадратична: f(x) = (Ax, x)/2 + (b, x) + c, то метод (5) имеет вид

xn+1 = xn - aS-1(Axn + b) = xn -a(S-1Axn + S-1

и является градиентным методом в старой метрике для функции f1 = (S-1Ax, x)/2 + (S-1b, x) + c. В соответствии с п. 3.8 при оптимальном выборе шага a его скорость сходимости линейна со знаменателем q = (L - l)/(L + l), где L и l максимальное и минимальное собственные значения оператора S-1A. Поэтому желательно сделать их разброс минимальным.

З а д а ч а  5.2. Покажите, что для квадратичной функции при S = A метод (5) сходится за один шаг.

В методе (5) оператор S можно менять на каждом шаге с той же целью ускорения сходимости. Такие методы называют иногда методами переменных направлений или переменной метрики.

В общем случае неквадратичной функции f в (5) оптимальным выбором в качестве S будет выбор Sn = f ўў(xn). Тогда (5) превращается просто в метод Ньютона со всеми присущими ему недостатками и достоинствами. С целью уменьшения объема вычислительной работы часто поступают следующим образом. Метод (5) записывают в виде

xn+1 = xn - Gnf ў(xn),(6)

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


Gn - [f ўў(xn)]-1 ® Q при n ® Ґ.

(7)

Методы, получающиеся таким способом часто называют квазиньютоновскими.

З а д а ч а  5.3. Докажите, что если f О C2, x* — невырожденная точка минимума функции f, а последовательность операторов Gn удовлетворяет условию (7), то метод (6) локально сверхлинейно сходится.

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

xn+1 = xn - anGnf ў(xn),(8)

имеет место соотношение


ansn = A-1jn, т. е. ansn = [f ўў(xn)]-1jn,

(9)

где

sn = -Gnf ў(xn), jn = f ў(xn+1) - f ў(xn).

З а д а ч а  5.4. Докажите.

Поэтому представляется естественным, чтобы операторы Gn+1, которые должны аппроксимировать [f ўў(xn)]-1, должны удовлетворять аналогу условия (9), а именно, т. н. квазиньютоновскому условию:

ansn = Gn+1jn.

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

Примером метода, построенного на этом пути, может служить один из наиболее эффективных среди квазиньютоновских методов метод Бройдена — Флетчера — Шенно, задаваемый формулой (8), в которой

Gn+1 = Gn + [gnsn(sn)* - sn(jn)*Gn - Gnjn(sn)*]/(sn, jn),

где

gn = an + (Gnjn, jn)/(sn, jn),

а "*" означает операцию транспонирования, в частности, (sn)* — вектор-строка. Шаг an в этом методе часто выбирают как в методе наискорейшего спуска в направлении sn:


an = argminaО[0, Ґ) f(xn + asn).

Для квадратичных функций метод Бройдена — Флетчера — Шенно конечен (выходит на точное решение не более, чем за m шагов). Для неквадратичных же функций он при достаточно общих условиях локально сверхлинейно сходится.

Число различных вариантов и модификаций квазиньютоновских методов весьма велико.

5.3. Методы сопряженных направлений.

Этот класс методов основан на следующем понятии. Векторы x и y в Rm называются сопряженными (относительно положительно определенного самосопряженного оператора A), если (Ax, y) = 0. Сопряженность векторов x и y очевидно означает их ортогональность относительно скалярного произведения бx, yс = (Ax, y); поэтому их называют еще A-ортогональными.

З а д а ч а  5.5. Докажите, что любой набор из m сопряженных векторов — базис в Rm.

Допустим, нам известен набор s1, ..., sm сопряженных векторов оператора A. Для квадратичной функции f(x) = (Ax, x)/2 + (b, x) + c построим релаксационный метод, выбирая sn в качестве направления спуска на n-ом шаге:

xn+1 = xn + ansn,(10)

а в качестве шага an — число

an = argminaОR f(xn + asn).(11)

З а д а ч а  5.6. Покажите, что an = -(f ў(xn), sn)/(f ў(sn), sn) = -(Axn, sn)/(Asn, sn).

Оказывается, метод (10)) при таком выборе sn и an сходится за конечное (Ј m) число шагов.

З а д а ч а  5.7. Докажите.

Отыскание набора сопряженных векторов — отдельная задача, которую решают разными способами. Например, можно взять произвольный базис в Rm и применить к нему известный процесс ортогонализации Грамма — Шмидта относительно скалярного произведения б·,·с.

В популярном методе Флетчера — Ривза минимизации квадратичных функций f(x) = (Ax, y)/2 + (b, x) + c векторы sn вычисляются по рекуррентным формулам

s0 = -f ў(x0),

sn+1 = -f ў(xn) + bn+1sn,

где числа bn+1 выбираются из условия A-ортогональности sn и sn+1:

0 = (Asn, sn+1) = -(Asn, f ў(xn)) + bn+1(Asn, sn),

откуда

bn+1 = (f ў(sn), f ў(xn))/(f ў(sn), sn).

(Напомним, что здесь f ў(x) = Ax (см. задачу 2.2)). Затем по формулам (10)(11) находят xn+1 и переходят к следующему шагу.

Для неквадратичных функций метод сопряженных направлений (называемый в этом случае обычно методом сопряженных градиентов), разумеется, уже не является конечным. Тем не менее, как итерационный он весьма эффективен. Например, метод Флетчера — Ривза может быть модифицирован для неквадратичных функций. Для этого нужно изменить только способ вычисления bn+1:

bn+1 = (f ў(xn) - f ў(xn-1), f ў(xn))/|| f ў(sn-1)||2.

Однако, поскольку для неквадратичных функций метод уже не конечен, на каждом кратном m шаге процедура построения сопряженных направлений повторяется заново с заменой s0 на sm = -f ў(xn).

Можно доказать, что метод Флетчера — Ривза для сильно выпуклых функций f О C3 локально сходится к невырожденному минимуму с порядком  mЦ2:

||xn+1 - x*|| Ј C||xn -x*||mЦ2,

или, за m шагов,

||xn+m - x*|| Ј C1||xn - x*||2.

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

5.4. Методы нулевого порядка.

В случаях, когда вычисление производных стóит дорого, применяют методы, не требующие вычисления производныхметоды нулевого порядка. Общая их схема такова. Итерации строятся в виде

xn+1 = xn - ansn,(12)

где an — длина шага, а направление sn шага обычно выбирается в виде


sn = 

1
h
n
е
i = 1

[f(xn + hei) - f(xn)]ei; 

(13)

здесь ei О Rm, i = 1, ..., k — некоторый, обычно фиксированный набор.

Часто в качестве ei выбираются векторы базиса в Rm. Очевидно, в этом случае при h ® Q выражение, стоящее в правой части (13) аппроксимирует f ў(xn). Тогда метод (12)(13) принимает вид


xn+1 = xn - 

an
h
n
е
i = 1
[f(xn + hei) - f(xn)]ei;

и называется разностным вариантом градиентного метода.

Если в (13) k = 1 и векторы базиса ei "циклически" меняются с изменением n, то получается метод покоординатного спуска:


xn+1 = xn - 

an
h

[f(xn + he{n/mm+1) -f(xn)]e{n/mm+1. 

где {a} — дробная часть числа a. На каждом шаге в нем изменяется только одна координата вектора xn (см. рис. 15).

Метод покоординатного спуска
Рис. 15.

Имеется масса модификаций методов нулевого порядка. В них по разному выбираются векторы ei, шаги ai, параметр h и  т. д.

5.5. Одномерная оптимизация.

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

f(x) ® min,   x О [a, b],(14)

в которой относительно функции f: [a, b]® R предполагается только, что на этом отрезке она имеет единственный миинимум x*, причем сдева от x* функция монотонно убывает, а справа — возрастает (такие функции называются унимодальными). Непрерывность f не предполагается.

З а д а ч а  5.8. Докажите, что всякая выпуклая функция унимодальна.

Основное свойство унимодальных функций, играющее для нас важную роль таково. Если мы знаем значения f в точках x1 и x2 отрезка [a, b], x1 < x2, то мы можем определенно сказать на каком из отрезков [a, x2], [x1, x2], [x1, b] лежит минимум x* функции f (см. рис. 16).

К локализации минимума унимодальной функции
Рис. 16.

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

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

5.6. Метод дихотомии.

Этот метод является аналогом известного метода дихотомии (половинного деления) отыскания решения уравнения F(x) = 0 на отрезке и выглядит так. Обозначим через Ln = [an, bn] отрезок, получающийся на n-ом шаге (естественно, L0 = [a, b]). Переход от Ln к Ln+1 сводится к размещению на Ln двух пробных точек xn и yn так, чтобы расстояние между ними было равно e и они были симметричны относительно центра отрезка:


xn = 

an + bn - e
2

,   yn = 

an + bn + e
2

и выбору в качестве Ln+1 отрезка [an, yn], если f(xn) Ј f(yn) и [xn, bn] в противном случае. Таким образом, каждый шаг требует двух вычислений функции и за n вычислений (n четно) мы получаем отрезок неопределенности длины


ln = 

l
2n/2
 +  ж
и
1 – l
2n/2
ц
ш
e,

где l = b - a — длина начального отрезка.

З а д а ч а  5.9. Докажите.

Например, если l = 1 и e << 1, то для получения отрезка длины 10-6 нужно порядка 40 вычислений функции.

5.7. Методы золотого сечения и Фибоначчи.

Эти методы основываются на одном и том же алгоритме перехода к следующей пробной точке и отличаются лишь способом выбора первой точки. В результате каждого шага (добавления одной пробной точки) мы получаем отрезок Ln = [an, bn] с одной "старой" пробной точкой xn О (an, bn). "Новая" пробная точка yn определяется как точка симметричная xn относительно центра отрезка Ln:


yn = 

an + bn
2
 + ж
и
an + bn
2

 - xn

ц
ш

 = an + bn - xn

Затем, в зависимости от того, какое из неравенств f(xn) Ј f(yn) или f(xn) > f(yn) выполняется, в качестве Ln+1 выбирается или отрезок [an, yn], или [xn, bn] (здесь мы считаем для определенности, что xn < yn), а в качестве старой пробной точки xn+1 выбирается либо xn, либо yn.

З а д а ч а  5.10*. Докажите, что ln-1 = ln + ln+1, где ln — длина отрезка Ln.

В методе золотого сечения x0 делит отрезок [a, b] по правилу золотого сечения:

b - a
x0 - a
 = t,

где t — положительный корень уравнения

t2 = t + 1

(золотое сечение): t = (1 + Ц5)/2 » 1.618). Таким образом, l0 = l1t. Но тогда (см. задачу 5.10) l2 = l0 - l1 = l1(t- 1) = l1/t, т. е. l1 = l2t.

З а д а ч а  5.11. Докажите, что ln+1 = ln/t.

Поэтому после n + 1 эксперимента


ln = 

l
tn
.

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

В методе Фибоначчи сначала задается число n вычислений функции и затем первая пробная точка подбирается таким образом, чтобы последняя пара экспериментов давала отрезок неопределенности наименьшей длины, т. е. чтобы последние две точки составляли симметричную относительно центра отрезка пару, расстояние между которыми равно e. Таким образом,

ln-1 = 2ln - e.

Учитывая, что (см. задачу 5.10)

ln-2 = ln-1 + ln,

получаем ln-2 = 3ln - e, ln-3 = 5ln- 2e, ln-4 = 8ln - 3e, ln-5 = 13ln - 5e и  т. д. Индукцией по i легко доказывается, что

ln-i = Fi+1ln - Fi-1e,(15)

где Fi — так называемые числа Фибоначчи, определяемые соотношениями F0 = F1 = 1, Fi = Fn-1 + Fn-2, i = 2,3, ..., введенные в XIII веке Леонардом Пизанским (Фибоначчи) совсем по другому поводу1). Из (15) мы получаем, во-первых, длину ln отрезка неопределенности, получающегося в результате вычисления f в n + 1 точках:

l = l0 = Fn+1ln - Fn-1e,

откуда


ln = 

l
Fn+1
 +  Fn-1
Fn+1
e,

и, во-вторых, длину отрезка l1, задающего положение первой (стартовой) точки:


l1 = Fnln - Fn-2e = 

Fn
Fn+1
lFnFn-1 - Fn-2Fn+1
Fn+1
e.
(16)

З а д а ч а  5.12. Докажите, что (Fi)2 = Fi-1Fi+1 + (-1)i. Используя это тождество, получите из (16) более простое соотношение, определяющее ln:


ln = 

Fn
Fn+1
l(-1)n
Fn+1
e.
(17)

Таким образом, чтобы запустить метод Фибоначчи необходимо заранее знать количество точек, в которых будет вычисляться функция f. В то же время метод Фибоначчи является наилучшим алгоритмом в таком классе методов2), в частности, по сравнению с методом золотого сечения, за одно и то же число шагов он дает отрезок неопределенности в t2/Ц5 » 1.17 раз меньший. На практике чаще используют первый из них, поскольку выигрыш в методе Фибоначчи не велик, а необходимость заранее знать число n вычислений функции зачастую обременительна.

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

5.8. Проблема нелокальности.

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

З а д а ч а  5.13. Докажите это утверждение.

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

К проблеме поиска глобального минимума
Рис. 17.

для которых, по-видимому, невозможно предложить ничего другого, кроме перебора значений на достаточно мелкой сетке с последующим уточнением решения, например, методом Ньютона. Представьте себе, что еще и m = dim Rm достаточно велико.

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

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

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

Еще один метод — так называемый овражный метод Гельфанда — основан на следующей идее: если функция имеет овражную структуру, то надо спуститься на дно оврага, а затем идти вдоль дна оврага. Шаги овражного метода двух типов: шаги спуска и овражные шаги. Шаги спуска начинаются из двух произвольных достаточно близких точек x0 и x1. Каким-нибудь методом, чаще всего, градиентным, спускаются на дно оврага, и получают точки X0 и X1. Затем делают овражный шаг из X0 в направлении X1 - X0 и получают точку x2:

x2 = X0 + a(X1 - X0).

Из точки x2 спускаются опять на дно в точку X2 и делают овражный шаг из X1 в направлении X2 - X1 и  т. д. (см. рис. 18).

Овражный метод Гельфанда
Рис. 18.

Овражный метод позволяет "просмотреть" окрестности дна оврага. Полученные точки можно затем уточнить с помощью локальных методов высокого порядка. Разумеется, применение овражного метода порождает массу вопросов. Например, одним из важных моментов является выбор длины овражного шага. Если шаг велик, то метод "проскакивает" локальные минимумы, если же шаг мал, то регулярного движения вдоль дна оврага не получается, особенно, если дно оврага имеет "уплощения".

В заключение параграфа опишем еще один метод поиска глобального минимума — метод спуск – подъем – перевал. Этот метод состоит из трех циклически повторяющихся этапов. Этап спуска — это обычный спуск из некоторой начальной точки x0 в ближайший локальный минимум, обычно методом сопряженных градиентов. Критерием попадания в малую окрестность точки локального минимума является "сходимость" итераций. В результате первого этапа мы получаем точку xN = y0, которая служит начальным приближением следующего этапа — этапа подъема. На этом этапе шаг совершается в направлении "наимедленнейшего" возрастания функции — направлении собственного вектора оператора f ўў(yn), отвечающего наименьшему собственному значению. Это направление находят с помощью специальной процедуры. Таким образом строят последовательность yn, "поднимающуюся по ложбине" до ближайшего "перевала" — седловой точки функции f. О выходе на "перевал" судят по смене знака некоторой квадратичной формы. После этого описанный трехэтапный цикл повторяется, начиная с последней ("перешедшей перевал") точки yM этапа подъема.


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

2) Думал ли об этом Фибоначчи?!


File based on translation from TEX by TTH, version 3.05.
Created 8 Jube 2002, 8:06.