Алгоритм Бройдена–Флетчера–Гольдфарба–Шенно

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

В числовой оптимизации алгоритм Бройдена-Флетчера-Гольдфарба-Шанно ( BFGS ) представляет собой итерационный метод решения задач нелинейной оптимизации без ограничений . ^[1] Как и родственный метод Дэвидона-Флетчера-Пауэлла , BFGS определяет направление спуска , предварительно определяя градиент с помощью информации о кривизне . Это достигается путем постепенного улучшения аппроксимации матрицы Гессе функции потерь , полученной только из оценок градиента (или приблизительных оценок градиента) с помощью обобщенного метода секущих . ^[2]

Поскольку обновления матрицы кривизны BFGS не требуют инверсии матрицы , ее вычислительная сложность составляет всего 1,0 по сравнению с методом Ньютона . Также широко используется L-BFGS , версия BFGS с ограниченной памятью, которая особенно подходит для задач с очень большим количеством переменных (например, >1000). Вариант BFGS-B обрабатывает простые ограничения блока. ^[3] ${\displaystyle {\mathcal {O}}(n^{2})}$ ${\displaystyle {\mathcal {O}}(n^{3})}$

Алгоритм назван в честь Чарльза Джорджа Бройдена , Роджера Флетчера , Дональда Голдфарба и Дэвида Шенно . ^{[4] [5] [6] [7]}

Обоснование

Задача оптимизации состоит в минимизации , где – вектор в и – дифференцируемая скалярная функция. Нет никаких ограничений на значения, которые могут принимать. ${\displaystyle f(\mathbf {x} )}$ ${\displaystyle \mathbf {x} }$ ${\displaystyle \mathbb {R} ^{n}}$ ${\displaystyle f}$ ${\displaystyle \mathbf {x} }$

Алгоритм начинается с начальной оценки оптимального значения и продолжается итеративно, чтобы получить лучшую оценку на каждом этапе. ${\displaystyle \mathbf {x} _{0}}$

Направление поиска p _k на этапе k задается решением аналога уравнения Ньютона:

${\displaystyle B_{k}\mathbf {p} _{k}=-\nabla f(\mathbf {x} _{k}),}$

где – аппроксимация матрицы Гессиана в точке , которая итеративно обновляется на каждом этапе, и – градиент функции, вычисляемой в точке x _k . Затем используется поиск линии в направлении p k для поиска следующей точки x _k +1 _{путем минимизации} по скаляру ${\displaystyle B_{k}}$ ${\displaystyle \mathbf {x} _{k}}$ ${\displaystyle \nabla f(\mathbf {x} _{k})}$ ${\displaystyle f(\mathbf {x} _{k}+\gamma \mathbf {p} _{k})}$ ${\displaystyle \gamma >0.}$

Условие квазиньютона, налагаемое на обновление, равно ${\displaystyle B_{k}}$

${\displaystyle B_{k+1}(\mathbf {x} _{k+1}-\mathbf {x} _{k})=\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k}).}$

Пусть и , тогда удовлетворяет ${\displaystyle \mathbf {y} _{k}=\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}$ ${\displaystyle \mathbf {s} _{k}=\mathbf {x} _{k+1}-\mathbf {x} _{k}}$ ${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}\mathbf {s} _{k}=\mathbf {y} _{k}}$,

что является секущим уравнением.

Чтобы быть положительно определенным, условие кривизны должно выполняться , что можно проверить, предварительно умножив уравнение секущего на . Если функция не является сильно выпуклой , то условие должно быть выполнено явно, например, путем нахождения точки x _{k +1} , удовлетворяющей условиям Вульфа , которые влекут за собой условие кривизны, с использованием поиска по прямой. ${\displaystyle \mathbf {s} _{k}^{\top }\mathbf {y} _{k}>0}$ ${\displaystyle B_{k+1}}$ ${\displaystyle \mathbf {s} _{k}^{T}}$

Вместо того, чтобы требовать вычисления полной матрицы Гессиана в точке как , приблизительный гессиан на этапе k обновляется путем добавления двух матриц: ${\displaystyle \mathbf {x} _{k+1}}$ ${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}=B_{k}+U_{k}+V_{k}.}$

Обе и являются симметричными матрицами первого ранга, но их сумма представляет собой матрицу обновления второго ранга. Матрица обновления BFGS и DFP отличается от своей предшественницы матрицей второго ранга. Другой более простой метод ранга один известен как симметричный метод ранга один , который не гарантирует положительную определенность . Чтобы сохранить симметрию и положительную определенность , форму обновления можно выбрать как . Наложение секущего условия, . Выбирая и , получим: ^[8] ${\displaystyle U_{k}}$ ${\displaystyle V_{k}}$ ${\displaystyle B_{k+1}}$ ${\displaystyle B_{k+1}=B_{k}+\alpha \mathbf {u} \mathbf {u} ^{\top }+\beta \mathbf {v} \mathbf {v} ^{\top }}$ ${\displaystyle B_{k+1}\mathbf {s} _{k}=\mathbf {y} _{k}}$ ${\displaystyle \mathbf {u} =\mathbf {y} _{k}}$ ${\displaystyle \mathbf {v} =B_{k}\mathbf {s} _{k}}$

${\displaystyle \alpha ={\frac {1}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}},}$

${\displaystyle \beta =-{\frac {1}{\mathbf {s} _{k}^{T}B_{k}\mathbf {s} _{k}}}.}$

Наконец, мы подставляем и в и получаем уравнение обновления : ${\displaystyle \alpha }$ ${\displaystyle \beta }$ ${\displaystyle B_{k+1}=B_{k}+\alpha \mathbf {u} \mathbf {u} ^{\top }+\beta \mathbf {v} \mathbf {v} ^{\top }}$ ${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}=B_{k}+{\frac {\mathbf {y} _{k}\mathbf {y} _{k}^{\mathrm {T} }}{\mathbf {y} _{k}^{\mathrm {T} }\mathbf {s} _{k}}}-{\frac {B_{k}\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} }B_{k}^{\mathrm {T} }}{\mathbf {s} _{k}^{\mathrm {T} }B_{k}\mathbf {s} _{k}}}.}$

Алгоритм

Рассмотрим следующую задачу неограниченной оптимизации.

$${\displaystyle {\begin{aligned}{\underset {\mathbf {x} \in \mathbb {R} ^{n}}{\text{minimize}}}\quad &f(\mathbf {x} ),\end{aligned}}}$$

${\displaystyle f:\mathbb {R} ^{n}\to \mathbb {R} }$

На основе первоначального предположения и первоначального предположения матрицы Гессе следующие шаги повторяются по мере приближения к решению: ${\displaystyle \mathbf {x} _{0}\in \mathbb {R} ^{n}}$ ${\displaystyle B_{0}\in \mathbb {R} ^{n\times n}}$ ${\displaystyle \mathbf {x} _{k}}$

1. Получите направление , решив . ${\displaystyle \mathbf {p} _{k}}$ ${\displaystyle B_{k}\mathbf {p} _{k}=-\nabla f(\mathbf {x} _{k})}$
2. Выполните одномерную оптимизацию ( линейный поиск ), чтобы найти приемлемый размер шага в направлении, найденном на первом шаге. Если производится точный поиск строки, то . На практике обычно бывает достаточно неточного поиска по строке с приемлемым выполнением условий Вульфа . ${\displaystyle \alpha _{k}}$ ${\displaystyle \alpha _{k}=\arg \min f(\mathbf {x} _{k}+\alpha \mathbf {p} _{k})}$ ${\displaystyle \alpha _{k}}$
3. Установите и обновите . ${\displaystyle \mathbf {s} _{k}=\alpha _{k}\mathbf {p} _{k}}$ ${\displaystyle \mathbf {x} _{k+1}=\mathbf {x} _{k}+\mathbf {s} _{k}}$
4. ${\displaystyle \mathbf {y} _{k}={\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}}$.
5. ${\displaystyle B_{k+1}=B_{k}+{\frac {\mathbf {y} _{k}\mathbf {y} _{k}^{\mathrm {T} }}{\mathbf {y} _{k}^{\mathrm {T} }\mathbf {s} _{k}}}-{\frac {B_{k}\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} }B_{k}^{\mathrm {T} }}{\mathbf {s} _{k}^{\mathrm {T} }B_{k}\mathbf {s} _{k}}}}$.

Сходимость можно определить, наблюдая за нормой градиента; учитывая некоторые , можно остановить алгоритм, если инициализируется с помощью , первый шаг будет эквивалентен градиентному спуску , но дальнейшие шаги все более и более уточняются с помощью , приближения к гессиану. ${\displaystyle \epsilon >0}$ ${\displaystyle ||\nabla f(\mathbf {x} _{k})||\leq \epsilon .}$ ${\displaystyle B_{0}}$ ${\displaystyle B_{0}=I}$ ${\displaystyle B_{k}}$

Первый шаг алгоритма выполняется с использованием обратной матрицы , которую можно эффективно получить, применив к шагу 5 алгоритма формулу Шермана–Моррисона , что дает ${\displaystyle B_{k}}$

${\displaystyle B_{k+1}^{-1}=\left(I-{\frac {\mathbf {s} _{k}\mathbf {y} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}\right)B_{k}^{-1}\left(I-{\frac {\mathbf {y} _{k}\mathbf {s} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}\right)+{\frac {\mathbf {s} _{k}\mathbf {s} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}.}$

Это можно эффективно вычислить без временных матриц, учитывая, что это симметрично, а что и являются скалярами, используя такое расширение, как ${\displaystyle B_{k}^{-1}}$ ${\displaystyle \mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}\mathbf {y} _{k}}$ ${\displaystyle \mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}$

${\displaystyle B_{k+1}^{-1}=B_{k}^{-1}+{\frac {(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}+\mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}\mathbf {y} _{k})(\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} })}{(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k})^{2}}}-{\frac {B_{k}^{-1}\mathbf {y} _{k}\mathbf {s} _{k}^{\mathrm {T} }+\mathbf {s} _{k}\mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}}{\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}}.}$

Следовательно, чтобы избежать инверсии матрицы, можно аппроксимировать обратную гессиану вместо самого гессиана: ^[9] ${\displaystyle H_{k}{\overset {\operatorname {def} }{=}}B_{k}^{-1}.}$

На основе первоначального предположения и приближенной обращенной матрицы Гессе следующие шаги повторяются по мере приближения к решению: ${\displaystyle \mathbf {x} _{0}}$ ${\displaystyle H_{0}}$ ${\displaystyle \mathbf {x} _{k}}$

1. Получите направление , решив . ${\displaystyle \mathbf {p} _{k}}$ ${\displaystyle \mathbf {p} _{k}=-H_{k}\nabla f(\mathbf {x} _{k})}$
2. Выполните одномерную оптимизацию ( линейный поиск ), чтобы найти приемлемый размер шага в направлении, найденном на первом шаге. Если производится точный поиск строки, то . На практике обычно бывает достаточно неточного поиска по строке с приемлемым выполнением условий Вульфа . ${\displaystyle \alpha _{k}}$ ${\displaystyle \alpha _{k}=\arg \min f(\mathbf {x} _{k}+\alpha \mathbf {p} _{k})}$ ${\displaystyle \alpha _{k}}$
3. Установите и обновите . ${\displaystyle \mathbf {s} _{k}=\alpha _{k}\mathbf {p} _{k}}$ ${\displaystyle \mathbf {x} _{k+1}=\mathbf {x} _{k}+\mathbf {s} _{k}}$
4. ${\displaystyle \mathbf {y} _{k}={\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}}$.
5. ${\displaystyle H_{k+1}=H_{k}+{\frac {(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}+\mathbf {y} _{k}^{\mathrm {T} }H_{k}\mathbf {y} _{k})(\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} })}{(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k})^{2}}}-{\frac {H_{k}\mathbf {y} _{k}\mathbf {s} _{k}^{\mathrm {T} }+\mathbf {s} _{k}\mathbf {y} _{k}^{\mathrm {T} }H_{k}}{\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}}}$.

В задачах статистической оценки ( таких как максимальное правдоподобие или байесовский вывод) достоверные интервалы или доверительные интервалы для решения могут быть оценены на основе обратной окончательной ^{матрицы Гессе .} Однако эти величины технически определяются истинной матрицей Гессе, и приближение BFGS может не сходиться к истинной матрице Гессе. ^[10]

Дальнейшие разработки

Формула обновления BFGS в значительной степени полагается на то, что кривизна строго положительна и не равна нулю. Это условие выполняется, когда мы выполняем поиск по прямой с условиями Вульфа на выпуклой цели. Однако некоторые реальные приложения (например, методы последовательного квадратичного программирования) обычно создают отрицательную или почти нулевую кривизну. Это может произойти при оптимизации невыпуклой цели или при использовании подхода доверительной области вместо линейного поиска. Также возможно получение ложных значений из-за шума в цели. ${\displaystyle \mathbf {s} _{k}^{\top }\mathbf {y} _{k}}$

В таких случаях можно использовать одно из так называемых демпфированных обновлений BFGS (см. ^[11] ), которые модифицируют и/или получают более надежное обновление. ${\displaystyle \mathbf {s} _{k}}$ ${\displaystyle \mathbf {y} _{k}}$

Известные реализации

Известные реализации с открытым исходным кодом:

• ALGLIB реализует BFGS и его версию с ограниченной памятью на C++ и C#.
• GNU Octave использует в своих функциях разновидность BFGS fsolveс расширениями доверенных регионов .
• GSL реализует BFGS как gsl_multimin_fdfminimizer_vector_bfgs2. ^[12]
• В R алгоритм BFGS (и версия L-BFGS-B, допускающая ограничения блоков) реализован как опция базовой функции optim(). ^[13]
• В SciPy функция scipy.optimize.fmin_bfgs реализует BFGS. ^[14] Также возможно запустить BFGS с использованием любого из алгоритмов L-BFGS , установив для параметра L очень большое число.
• В Julia пакет Optim.jl реализует BFGS и L-BFGS в качестве опции решателя для функцииоптимизации() (среди других опций). ^[15]

Известные собственные реализации включают:

• Программное обеспечение для крупномасштабной нелинейной оптимизации Artelys Knitro реализует, среди прочего, алгоритмы BFGS и L-BFGS.
• В MATLAB Optimization Toolbox функция fminunc использует BFGS с поиском кубической линии, когда размер проблемы установлен на «средний масштаб».
• Mathematica включает BFGS.

Смотрите также

• Алгоритм БХХХ
• Формула Дэвидона – Флетчера – Пауэлла
• Градиентный спуск
• Л-БФГС
• Алгоритм Левенберга – Марквардта
• Метод Нелдера-Мида
• Поиск по шаблону (оптимизация)
• Квазиньютоновские методы
• Симметричный первого ранга

Рекомендации

1. Флетчер, Роджер (1987), Практические методы оптимизации (2-е изд.), Нью-Йорк: John Wiley & Sons , ISBN 978-0-471-91547-8
2. Деннис, Дж. Э. младший ; Шнабель, Роберт Б. (1983), «Методы секущих для неограниченной минимизации», Численные методы для неограниченной оптимизации и нелинейных уравнений , Энглвуд Клиффс, Нью-Джерси: Прентис-Холл, стр. 194–215, ISBN 0-13-627216-9
3. Берд, Ричард Х.; Лу, Пэйхуан; Носедаль, Хорхе; Чжу, Цию (1995), «Алгоритм с ограниченной памятью для оптимизации с ограниченными ограничениями», SIAM Journal on Scientific Computing , 16 (5): 1190–1208, CiteSeerX 10.1.1.645.5814 , doi : 10.1137/0916069 
4. Бройден, К.Г. (1970), «Сходимость класса двухранговых алгоритмов минимизации», Журнал Института математики и ее приложений , 6 : 76–90, doi : 10.1093/imamat/6.1.76
5. Флетчер, Р. (1970), «Новый подход к алгоритмам с переменной метрикой», Computer Journal , 13 (3): 317–322, doi : 10.1093/comjnl/13.3.317
6. Гольдфарб, Д. (1970), «Семейство обновлений переменных метрик, полученных вариационными средствами», Mathematics of Computing , 24 (109): 23–26, doi : 10.1090/S0025-5718-1970-0258249-6
7. Шанно, Дэвид Ф. (июль 1970 г.), «Обусловление квазиньютоновских методов для минимизации функций», Mathematics of Computation , 24 (111): 647–656, doi : 10.1090/S0025-5718-1970-0274029-X , МР  0274029
8. Флетчер, Роджер (1987), Практические методы оптимизации (2-е изд.), Нью-Йорк: John Wiley & Sons , ISBN 978-0-471-91547-8
9. Носедаль, Хорхе; Райт, Стивен Дж. (2006), Численная оптимизация (2-е изд.), Берлин, Нью-Йорк: Springer-Verlag , ISBN 978-0-387-30303-1
10. Ге, Рен-пу; Пауэлл, MJD (1983). «Сходимость матриц переменных показателей при безусловной оптимизации». Математическое программирование . 27 (2). 123. дои : 10.1007/BF02591941. S2CID  8113073.
11. Хорхе Носедал; Стивен Дж. Райт (2006), Численная оптимизация
12. «Научная библиотека GNU — документация GSL 2.6» . www.gnu.org . Проверено 22 ноября 2020 г.
13. «R: Универсальная оптимизация». stat.ethz.ch. _ Проверено 22 ноября 2020 г.
14. «scipy.optimize.fmin_bfgs — Справочное руководство SciPy v1.5.4» . docs.scipy.org . Проверено 22 ноября 2020 г.
15. «Настраиваемые параметры Optim.jl» . Джулианлсолверс .

дальнейшее чтение

• Авриэль, Мордекай (2003), Нелинейное программирование: анализ и методы , Dover Publishing, ISBN 978-0-486-43227-4
• Боннан, Ж. Фредерик; Гилберт, Дж. Чарльз; Лемарешаль, Клод ; Сагастисабал, Клаудия А. (2006), «Ньютоновские методы», Численная оптимизация: теоретические и практические аспекты (второе изд.), Берлин: Springer, стр. 51–66, ISBN 3-540-35445-Х
• Флетчер, Роджер (1987), Практические методы оптимизации (2-е изд.), Нью-Йорк: John Wiley & Sons , ISBN 978-0-471-91547-8
• Люенбергер, Дэвид Г .; Йе, Инью (2008), Линейное и нелинейное программирование , Международная серия по исследованию операций и науке управления, том. 116 (Третье изд.), Нью-Йорк: Springer, стр. xiv+546, ISBN. 978-0-387-74502-2, МР  2423726
• Келли, Коннектикут (1999), Итеративные методы оптимизации , Филадельфия: Общество промышленной и прикладной математики, стр. 71–86, ISBN. 0-89871-433-8