Вычисление матриц Якоби и Гессе

Материал из MachineLearning.

Перейти к: навигация, поиск

Содержание

Введение

Постановка математической задачи

Вычисление матрицы Якоби

Пусть задана система n функций y_1(x_1, x_2, \dots x_m) \dots y_n(x_1, x_2, \dots x_m) от m переменных. Матрицей Якоби или якобианом данной системы функций называется матрица, составленная из частных производных этих функций по всем переменным.


J=\begin{bmatrix} \frac{\partial y_1}{\partial x_1} & \cdots & \frac{\partial y_1}{\partial x_m} \\ \vdots & \ddots & \vdots \\ \frac{\partial y_n}{\partial x_1} & \cdots & \frac{\partial y_n}{\partial x_m} \end{bmatrix}.

Если в некоторой точке x_1 \dots x_m очень сложно или невозможно вычислить частные производные, \frac{\partial y_i}{\partial x_j} , то для вычисления матрицы Якоби применяются методы численного дифференцирования.

Вычисление матрицы Гессе

Матрицей Гессе функции m переменных y(x_1 \dots x_m) называется матрица, составленная из вторых производных функции y(x_1 \dots x_m) по всем переменным


H(f) = \begin{bmatrix}
\frac{\partial^2 y}{\partial x_1^2} & \frac{\partial^2 y}{\partial x_1\,\partial x_2} & \cdots & \frac{\partial^2 y}{\partial x_1\,\partial x_m} \\  \\
\frac{\partial^2 y}{\partial x_2\,\partial x_1} & \frac{\partial^2 y}{\partial x_2^2} & \cdots & \frac{\partial^2 y}{\partial x_2\,\partial x_m} \\  \\
\vdots & \vdots & \ddots & \vdots \\  \\
\frac{\partial^2 y}{\partial x_m\,\partial x_1} & \frac{\partial^2 y}{\partial x_m\,\partial x_2} & \cdots & \frac{\partial^2 y}{\partial x_m^2}
\end{bmatrix}

Если в некоторой точке x_1 \dots x_m очень сложно или невозможно вычислить частные производные, \frac{\partial^2 y}{\partial x_i\ \partial x_j} , то для вычисления матрицы Гессе применяются методы численного дифференцирования.

Методы вычисления матрицы Якоби

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

Для вычисления матрицы Якоби в заданной необходимо найти частные производные всех функций системы по всем переменным. Для вычисления производной \frac{\partial y_i}{\partial x_j} применяются методы вычисления первой производной.

Формула для элемента якобиана при использовании правой разностной производной:
 J_{ij} = \frac {y_i(x_1 \dots x_j + h_{ij} \dots x_m) - y_i(x_1 \dots x_m)} {h_{ij} }

Формула для элемента якобиана при использовании центральной разностной производной:
 J_{ij} = \frac {y_i(x_1 \dots x_j + h_{ij} \dots x_m) - y_i(x_1 \dots x_j - h_{ij} \dots x_m)} {2*h_{ij}}



Вычисление якобиана с использованием правой разностной производной требует вычислять значения функций в m * (n + 1) точках. Если использовать центральную производную, то нужно находить значения функций  y_1 \dots y_n в 2 * m * n   точках. С другой, стороны погрешность правой производной имеет порядок O(h) а центральной - O(h^2). В большинстве случаев вычисление значения функции y_i - это затратная по времени операция, поэтому используется правая разностная производная.

Оценка погрешности метода

Основная проблема при вычислении каждого элемента матрицы Якоби - как правильно выбрать шаг метода  h_{ij} . Шаг выбирается независимо для каждого элемента матрицы.

Проанализируем зависимость погрешности метода от величины шага в случае использования правой разностной производной. Для сокращения записи введём обозначения f(\xi_j) = y_i(x_1 \dots \xi_j \dots x_m). Остаточный член в соотношении  f(x)' \approx \frac{f(x+ h) - f1}{h} имеет вид r_1 = - f''(\xi)h/2 . Если  |f''(\xi)| \le M_2 , то |r_1| \le M_2*h/2 Если значения f1 и f2 заданы с погрешностями \vareps_1, \vareps_2 \le E , то погрешность будет содержать ещё одно слагаемое  r_2 = - \frac{\vareps_1 - \vareps_2}{h} \le 2E/h . Таким образом, оценка суммарной погрешности имеет вид |r| \le |r_1| + |r_2| \le f(h) = M_2 h/2 + 2E/h . Эта оценка достигает минимума при h_0 = 2 \sqrt {E/M_2}. При этом  f(h_0) = 2 \sqrt{M_2 E} . Оценка погрешности имеет один глобальный миниум. Поэтому выбор очень маленького шага не привидёт к росту точности. При величине шага, близкой к h_0 погрешность имеет порядок O(\sqrt E) .

Метод Бройдена

Чаще всего вычисление якобиана является одной из подзадач в различных методах оптимизации и решения систем нелинейных уравнений. При решении систем нелинейных уравнений методом Ньютона требуется вычислять якобиан на каждой итерации. Вычисление якобиана требует вычисления функций в O(m * n) точках. Это сложная и затратная по времени операция. Суть метода Бройдена состоит в том, чтобы вычислить якобиан аналитически или с помощью метода конечных разностей на первой итерации, а после этого на каждой итерации обновлять якобиан, не вычисляя значения функций y_1 \dots y_n и их производных.

Пусть задана система нелинейных уравнений F(x) = 0, где  F:R^m \to R^n . Тогда якобиан на n-ой итерации выражается по формуле

J_n = J_{n - 1} + \frac{F(x_n) - F(x_{n-1}) - J_{n - 1} (x_n - x_{n-1})} { ||x_n - x_{n-1}||^2} (x_n - x_{n-1})^T

После этого следующее приближение вычисляется по формуле

x_{n + 1} = x_n - J_n^{-1} F(x_n)

Методы вычисления матрицы Гессе

Как и матрица Якоби, матрица Гессе может быть вычислена с помощью разностной аппроксимации производных.  \frac{\partial^2 f}{\partial x_i\,\partial x_j} = 
\frac{ f(x + (e_i + e_j) h) - f(x + e_i h) - f(x + e_j h) + f(x) } {h^2}, где x - вектор переменных, а e_i и e_j - единичные вектора. Эта формула требует вычисления значений функции  f в  \frac {n (n + 1)} {2} точках. Погрешность формулы имеет порядок O(h).


Численный эксперимент

В качестве примера рассчитаем с помощью вышеизложенного метода матрицу Гессе функции x^2 * \sqrt y в точке (1, 1)

Величина шагаЧисленный результат Истинное значение Погрешность
0.01

 
\begin{bmatrix} 1.999999999999780 &  1.002499984530392 \\ \\  1.002499984530392 & -0.250007812910846 \end{bmatrix}

 
\begin{bmatrix} 2 &  1 \\ \\  1 & -0.25 \end{bmatrix}

 
\begin{bmatrix} 
0.000000000000220 & 0.002499984530392 \\ \\ 0.002499984530392 & 0.500007812910846   \end{bmatrix}

0.001

 
\begin{bmatrix} 1.999999999724444 & 1.000249999938419 \\1.000249999938419 & -0.250000078083623 \end{bmatrix}

 
\begin{bmatrix} 2 &  1 \\ \\  1 & -0.25 \end{bmatrix}

 
\begin{bmatrix}
0.000000000275556 & 0.000249999938418 \\ \\ 0.000249999938418 & 0.500000078083623  \end{bmatrix}

0.0001

 
\begin{bmatrix} 1.999999998947288 & 1.000024996145044 \\ \\ 1.000024996145044 & -0.250000009582863 \end{bmatrix}

 
\begin{bmatrix} 2 &  1 \\ \\  1 & -0.25 \end{bmatrix}

 
\begin{bmatrix} 0.000000001052712 & 0.000024996145044 \\ \\ 0.000024996145044 & 0.500000009582863  \end{bmatrix}

0.00001

 
\begin{bmatrix} 2.000002385926791 & 1.000002303186420 \\ \\ 1.000002303186420 & -0.250000020685093  \end{bmatrix}

 
\begin{bmatrix} 2 &  1 \\ \\  1 & -0.25 \end{bmatrix}

 
\begin{bmatrix} 0.000002385926791 & 0.000002303186420 \\ \\ 0.000002303186420 & 0.500000020685093  \end{bmatrix}


Рекомендации программисту

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

Программа для вычисления якобиана, исходный код [1кб]

Программа для вычисления гессиана, исходный код [1кб]

Библиотека алгоритмов линейной алгебры [180кб]


Заключение

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

  • Самарский А. А. , Гулин А. В. Численные методы. Москва «Наука», 1989.
  • Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. Численные методы. Лаборатория Базовых Знаний, 2003.
  • Магнус Я. Р., Нейдеккер Х. Матричное дифференциальное исчисление с приложениями к статистике и эконометрике: Пер. с англ./ Под ред. С. А. Айвазяна. — М.:ФИЗМАТЛИТ, 2002. — 496 с.
  • Хайкин С. Нейронные сети, полный курс. 2е издание, испр. — М: Вильямс. 2008. — 1103 с ISBN 978-5-8459-0890-2
  • Bishop C. Exact calculation of the Hessian matrix for the multilayer perceptron / Neural Computation. Volume 4, Issue 4 (July 1992). Pp. 494—501. ISSN:0899-7667
  • MacKay D. How do you evaluate the Hessian in the first place?
Личные инструменты