Метод коллинеарных градиентов

Метод коллинеарных градиентов (МКГ)^[1] — итерационный метод направленного поиска локального экстремума гладкой функции многих переменных ${\displaystyle J(u)\colon \mathbb {R} ^{n}\to \mathbb {R} }$ с движением к экстремуму вдоль вектора ${\displaystyle d\in \mathbb {R} ^{n}}$ такого, где градиенты ${\displaystyle \nabla J(u)}$ и ${\displaystyle \nabla J(u+d)}$ коллинеарные. Это метод перового порядка (использует только первые производные ${\displaystyle \nabla J}$) с квадратичной скоростью сходимости. Может применяться к функциям высокой размерности ${\displaystyle n}$ с несколькими локальными экстремумами. МКГ можно отнести к семейству методов Truncated Newton method

Для гладкой функции ${\displaystyle J(u)}$ в относительно большой окрестности точки ${\displaystyle u^{k}}$ найдётся точка ${\displaystyle u^{k_{\ast }}}$, где градиенты ${\displaystyle \nabla J^{k}\,{\overset {\textrm {def}}{=}}\,\nabla J(u^{k})}$ и ${\displaystyle \nabla J^{k_{\ast }}\,{\overset {\textrm {def}}{=}}\,\nabla J(u^{k_{\ast }})}$ коллинеарные. Направлением на экстремум ${\displaystyle u_{\ast }}$ из точки ${\displaystyle u^{k}}$ будет направление ${\displaystyle d^{k}=(u^{k_{\ast }}-u^{k})}$. Вектор ${\displaystyle d^{k}}$ указывает на максимум или на минимум в зависимости от положения точки ${\displaystyle u^{k_{\ast }}}$. Она может быть спереди или сзади от ${\displaystyle u^{k}}$ относительно направления на ${\displaystyle u_{\ast }}$ (см. рисунок). Далее будем рассматривать минимизацию.

Очередная итерация МКГ:

(1) ${\displaystyle \quad u^{k+1}=u^{k}+b^{k}d^{k},\quad k\in \left\{0,1\dots \right\},}$

где оптимальное ${\displaystyle b^{k}\in \mathbb {R} }$ находится аналитически из предположения квадратичности одномерной функции ${\displaystyle J(u^{k}+bd^{k})}$:

(2) ${\displaystyle \quad b^{k}=\left(1-{\frac {\langle \nabla J(u^{k_{\ast }},d^{k}\rangle }{\langle \nabla J(u^{k}),d^{k}\rangle }}\right)^{-1},\quad \forall u^{k_{\ast }}.}$

Угловые скобки — это скалярное произведение в евклидовом пространстве ${\displaystyle \mathbb {R} ^{n}}$. Если ${\displaystyle J(u)}$ выпуклая функция в окрестности ${\displaystyle u^{k}}$, то для передней точки ${\displaystyle u^{k_{\ast }}}$ получаем число ${\displaystyle b^{k}>0}$, для задней ${\displaystyle b^{k}<0}$. Делаем шаг (1).

Для строго выпуклой квадратичной функции ${\displaystyle J(u)}$ шаг МКГ

 ${\displaystyle b^{k}d^{k}=-H^{-1}\nabla J^{k},}$

т.е. это шаг метода Ньютона (метод второго порядка с квадратичной скоростью сходимости), где ${\displaystyle H}$ — матрица Гессе. Такие шаги обеспечивают МКГ квадратичную скорость сходимости.

В общем случае, если ${\displaystyle J(u)}$ имеет переменную выпуклость и возможны седловые точки, то следует контролировать направление минимизации по углу ${\displaystyle \gamma }$ между векторами ${\displaystyle \nabla J^{k}}$ и ${\displaystyle d^{k}}$. Если ${\displaystyle \cos(\gamma )={\frac {\langle \nabla J^{k},d^{k}\rangle }{||\nabla J(u^{k})||\;||d^{k}||}}\geq 0}$, то ${\displaystyle d^{k}}$ — это направление максимизации и в (1) следует брать ${\displaystyle b^{k}}$ с обратным знаком.

Коллинеарность градиентов оценивается невязкой их ортов, которая имеет вид системы ${\displaystyle n}$ уравнений для поиска корня ${\displaystyle u=u^{k_{\ast }}}$:

(3) ${\displaystyle \quad r^{k}(u)={\frac {\nabla J(u)}{||\nabla J(u)||}}s-{\frac {\nabla J(u^{k})}{||\nabla J(u^{k})||}}=0\in \mathbb {R} ^{n},}$

где знак ${\displaystyle s=\operatorname {sgn} \langle \nabla J(u),\nabla J(u^{k})\rangle }$ позволяет одинаково оценивать коллинеарность градиентов по одну или разные стороны от минимума ${\displaystyle u_{\ast }}$, ${\displaystyle ||r^{k}(u)||\leq {\sqrt {2}}}$.

Система (3) решается итерационно (подитерации ${\displaystyle l\,}$) методом сопряжённых градиентов в предположении, что она линейна в окрестности ${\displaystyle u^{k}}$:

(4) ${\displaystyle \quad u^{k_{l+1}}=u^{k_{l}}+\tau ^{l}p^{l},\quad l=1,2\ldots ,}$

где вектор ${\displaystyle \;p^{l}\;{\overset {\textrm {def}}{=}}\,p(u^{k_{l}})=-r^{l}+{\beta ^{l}p}^{l-1}}$, ${\displaystyle \;r^{l}\,{\overset {\textrm {def}}{=}}\,r(u^{k_{l}})}$, ${\displaystyle \;\beta ^{l}=||r^{l}||^{2}/||r^{l-1}||^{2},\ \beta ^{1,n,2n...}=0}$, ${\displaystyle \;\tau ^{l}=||r^{l}||^{2}/\langle p^{l},H^{l}p^{l}\rangle }$, произведение матрицы Гессе ${\displaystyle H^{l}}$ на ${\displaystyle p^{l}}$ находится численным дифференцированием:

(5) ${\displaystyle \quad H^{l}p^{l}\approx {\frac {r(u^{k_{h}})-r(u^{k_{l}})}{h/||p^{l}||}},}$

где ${\displaystyle u^{k_{h}}=u^{k_{l}}+hp^{l}/||p^{l}||}$, ${\displaystyle h}$ — малое положительное число такое, что ${\displaystyle \langle p^{l},H^{l}p^{l}\rangle \neq 0}$.

Начальное приближение задаётся под 45° ко всем осям координат длинной ${\displaystyle \delta ^{k}}$:

(6) ${\displaystyle \quad u_{i}^{k_{1}}=u_{i}^{k}+{\frac {\delta ^{k}}{\sqrt {n}}}\operatorname {sgn} {\ \nabla _{i}J}^{k},\quad i=1\ldots n.}$

Начальный радиус ${\displaystyle \delta ^{k}}$-окрестности точки ${\displaystyle u^{k}}$ корректируется:

(7) ${\displaystyle \quad \delta ^{k}=\max \left[\min \left(\delta ^{k-1}{\frac {||\nabla J(u^{k})||}{||\nabla J(u^{k-1})||}},\delta ^{0}\right),\delta _{m}\right],\quad k>0.}$

Необходимо ${\displaystyle ||u^{k_{l}}-u^{k}||\geq \delta ^{m},\;l\geq 1}$. Здесь малое положительное число ${\displaystyle \delta _{m}}$ заметно больше машинного эпсилон.

Подитерации ${\displaystyle l}$ завершаются при выполнении хотя бы одного из условий:

1. ${\displaystyle ||r^{l}||\leq c_{1}{\sqrt {2}},\quad 0\leq c_{1}<1}$ — достигнута точность;
2. ${\displaystyle \left|{\frac {||r^{l}||-||r^{l-1}||}{||r^{l}||}}\right|\leq c_{1},\quad l>1}$ — прекратилась сходимость;
3. ${\displaystyle l\leq l_{max}=\operatorname {integer} \left|c_{2}\ln c_{1}\ln n\right|,\quad c_{2}\geq 1}$ — избыточность подитераций.

• Параметры: ${\displaystyle c_{1},c_{2},\delta ^{0},\delta _{m}=10^{-15}\delta ^{0},h=10^{-5}}$.
• Входные данные: ${\displaystyle n,k=0,u^{k},J(u^{k}),\nabla J(u^{k}),\nabla J^{k}}$.

1. ${\displaystyle l=1}$. Если ${\displaystyle k>0}$ задаём ${\displaystyle \delta ^{k}}$ из (7).
2. Находим ${\displaystyle u^{k_{l}}}$ из (6).
3. Вычисляем ${\displaystyle \nabla J^{l},||\nabla J^{l}||}$ и находим ${\displaystyle r^{l}}$ из (3) при ${\displaystyle u=u^{k_{l}}}$.
4. Если ${\displaystyle ||r^{l}||\leq c_{1}{\sqrt {2}}\,}$ или ${\displaystyle l\geq l_{max}}$, или ${\displaystyle ||u^{k_{l}}-u^{k}||<\delta _{m}}$, или {${\displaystyle \,l>1}$ и ${\displaystyle \left|{\frac {||r^{l}||-||r^{l-1}||}{||r^{l}||}}\right|\leq c_{1}}$}, то принимаем ${\displaystyle u^{k_{\ast }}=u^{k_{l}}}$, возвращаем ${\displaystyle \nabla J\left(u^{k_{\ast }}\right)}$, ${\displaystyle d^{k}={(u^{k_{\ast }}-u}^{k})}$, стоп.
5. Если ${\displaystyle l\neq 1,n,2n,3n\ldots }$, задаём ${\displaystyle \beta ^{l}=||r^{l}||^{2}/||r^{l-1}||^{2}}$, иначе ${\displaystyle \beta ^{l}=0}$.
6. Вычисляем ${\displaystyle p^{l}=-r^{l}+\beta ^{l}p^{l-1}}$.
7. Находим шаговый множитель ${\displaystyle \tau ^{l}}$ для подитераций:
   1. запоминаем ${\displaystyle u^{k_{l}}}$, ${\displaystyle \nabla J^{l}}$, ${\displaystyle ||\nabla J^{l}||}$, ${\displaystyle r^{l}}$, ${\displaystyle ||r^{l}||}$;
   2. задаём ${\displaystyle u^{k_{h}}=u^{k_{l}}+hp^{l}/||p^{l}||}$, вычисляем ${\displaystyle \nabla J(u^{k_{h}})}$, ${\displaystyle r\left(u^{k_{h}}\right)}$ и находим ${\displaystyle H^{l}p^{l}}$ из (5), присваиваем ${\displaystyle w\leftarrow \langle p^{l},H^{l}p^{l}\rangle }$;
   3. если ${\displaystyle w=0}$, тогда ${\displaystyle h\leftarrow 10h}$, возвращаемся к шагу 7.2;
   4. восстанавливаем ${\displaystyle u^{k_{l}}}$, ${\displaystyle \nabla J^{l}}$, ${\displaystyle ||\nabla J^{l}||}$, ${\displaystyle r^{l}}$, ${\displaystyle ||r^{l}||}$;
   5. находим ${\displaystyle \tau ^{l}=||r^{l}||^{2}/w}$.
8. Делаем подитерацию ${\displaystyle u^{k_{l+1}}}$ из (4).
9. ${\displaystyle l\leftarrow l+1}$, переходим к шагу 3.

Параметр ${\displaystyle c_{2}=3\div 5}$. Для функций без седловых точек рекомендуется ${\displaystyle c_{1}\approx 10^{-8}}$, ${\displaystyle \delta \approx {10}^{-5}}$. Для «обхода» седловых точек рекомендуется ${\displaystyle c_{1}\approx 0.1}$, ${\displaystyle \delta \approx 0.1}$.

Описанный алгоритм позволяет приблизительно найти коллинеарные градиенты из системы уравнений (3). Полученное направление ${\displaystyle b^{k}d^{k}}$ для алгоритма МКГ (1) будет приблизительным направлением Ньютона (truncated Newton method).

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

Тестовая функция «повёрнутый эллипсоид»

Строго выпуклая квадратичная функция:

${\displaystyle J(u)=\sum _{i=1}^{n}\left(\sum _{j=1}^{i}u_{j}\right)^{2},\quad u_{\ast }=(0...0).}$

На рисунке для ${\displaystyle {\color {red}n=2}}$ заданы три чёрные стартовые точки ${\displaystyle u^{0}}$. Серые точки — подитерации ${\displaystyle u^{0_{l}}}$ с ${\displaystyle \delta ^{0}=0.5}$ (показано пунктиром, завышено для демонстрации). Параметры ${\displaystyle c_{1}=10^{-8}}$, ${\displaystyle c_{2}=4}$. Для всех ${\displaystyle u^{0}}$ потребовалась одна итерация и подитераций ${\displaystyle l}$ не более двух.

При ${\displaystyle {\color {red}n=1000}}$ (параметр ${\displaystyle \delta ^{0}={10}^{-5}}$) с начальной точкой ${\displaystyle u^{0}=(-1...1)}$ МКГ достиг ${\displaystyle u_{\ast }}$ с точностью 1 % за 3 итерации и 754 вычисления ${\displaystyle J}$ и ${\displaystyle \nabla J}$. Другие методы первого порядка: Квазиньютоновский BFGS (работа с матрицами) потребовал 66 итераций и 788 вычислений; сопряжённых градиентов (Fletcher-Reeves) — 274 итерации и 2236 вычислений; конечно-разностный метод Ньютона — 1 итерация и 1001 вычислений. Метод Ньютона второго порядка — 1 итерация.

С ростом размерности ${\displaystyle \color {red}n}$, вычислительные погрешности при реализации условия коллинеарности (3), могут заметно возрастать. Поэтому МКГ, по сравнению с методом Ньютона, в рассматриваемом примере потребовал более одной итерации.

Тестовая функция Химмельблау

${\displaystyle J(u)=(u_{1}^{2}+u_{2}-11)^{2}+(u_{1}+u_{2}^{2}-7)^{2}.}$

Параметры ${\displaystyle c_{1}=0.1}$, ${\displaystyle c_{2}=4}$, ${\displaystyle \delta ^{0}=0.05}$.

МКГ является очень экономичным по количеству вычислений ${\displaystyle J}$ и ${\displaystyle \nabla J}$. Благодаря формуле (2), он не требует затратных вычислений шагового множителя ${\displaystyle b^{k}}$ посредством линейного поиска (например, методом золотого сечения и т.п.).