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

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

(Различия между версиями)
Перейти к: навигация, поиск
(Вычисление матрицы Гессе)
 
(15 промежуточных версий не показаны.)
Строка 5: Строка 5:
==== Вычисление матрицы Якоби ====
==== Вычисление матрицы Якоби ====
-
Пусть задана система <tex>m</tex> функций <tex>y_1(x_1, x_2, \dots x_n) \dots y_m(x_1, x_2, \dots x_n)</tex> от <tex>m</tex> переменных. '''Матрицей Якоби''' данной системы функций называется
+
Пусть задана система <tex>n</tex> функций <tex>y_1(x_1, x_2, \dots x_m) \dots y_n(x_1, x_2, \dots x_m)</tex> от <tex>m</tex> переменных. '''Матрицей Якоби''' или '''якобианом''' данной системы функций называется
матрица, составленная из частных производных этих функций по всем переменным.
матрица, составленная из частных производных этих функций по всем переменным.
-
<p align = "center">
 
<tex>
<tex>
-
J=\begin{bmatrix} \frac{\partial y_1}{\partial x_1} & \cdots & \frac{\partial y_1}{\partial x_n} \\ \vdots & \ddots & \vdots \\ \frac{\partial y_m}{\partial x_1} & \cdots & \frac{\partial y_m}{\partial x_n} \end{bmatrix}.
+
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}.
</tex>
</tex>
-
</p>
 
-
Если в некоторой точке <tex>x_1 \dots x_m </tex> очень сложно или невозможно вычислить частные производные, <tex>\frac{\partial y_i}{\partial x_j} </tex>, то для вsчисления матрицы Якоби применяются методы численного дифференцирования.
+
Если в некоторой точке <tex>x_1 \dots x_m </tex> очень сложно или невозможно вычислить частные производные, <tex>\frac{\partial y_i}{\partial x_j} </tex>, то для вычисления матрицы Якоби применяются методы численного дифференцирования.
 +
==== Вычисление матрицы Гессе ====
 +
'''Матрицей Гессе''' функции <tex>m</tex> переменных <tex>y(x_1 \dots x_m) </tex> называется матрица, составленная из вторых производных функции <tex>y(x_1 \dots x_m) </tex> по всем переменным
 +
<tex>
 +
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}
 +
</tex>
-
== Изложение метода ==
+
Если в некоторой точке <tex>x_1 \dots x_m </tex> очень сложно или невозможно вычислить частные производные, <tex>\frac{\partial^2 y}{\partial x_i\ \partial x_j} </tex>, то для вычисления матрицы Гессе применяются методы численного дифференцирования.
 +
 
 +
== Методы вычисления матрицы Якоби ==
 +
=== Прямое вычисление частных производных ===
 +
Для вычисления матрицы Якоби в заданной необходимо найти частные производные всех функций системы по всем переменным. Для вычисления
 +
производной <tex>\frac{\partial y_i}{\partial x_j} </tex> применяются методы [[ Вычисление первой производной|вычисления первой производной]].
 +
 
 +
Формула для элемента якобиана при использовании правой разностной производной: <br/>
 +
<tex> J_{ij} = \frac {y_i(x_1 \dots x_j + h_{ij} \dots x_m) - y_i(x_1 \dots x_m)} {h_{ij} } </tex>
 +
 
 +
Формула для элемента якобиана при использовании центральной разностной производной:
 +
<br/>
 +
<tex> 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}} </tex>
 +
 
 +
<br/>
 +
 
 +
 
 +
Вычисление якобиана с использованием правой разностной производной требует вычислять значения функций в <tex>m * (n + 1) </tex>
 +
точках. Если использовать центральную производную, то нужно находить значения функций <tex> y_1 \dots y_n </tex> в <tex>2 * m * n </tex> точках. С другой, стороны погрешность правой производной имеет порядок <tex>O(h)</tex> а центральной - <tex>O(h^2)</tex>.
 +
В большинстве случаев вычисление значения функции <tex>y_i</tex> - это затратная по времени операция, поэтому используется правая разностная производная.
 +
 
 +
====Оценка погрешности метода====
 +
Основная проблема при вычислении каждого элемента матрицы Якоби - как правильно выбрать шаг метода <tex> h_{ij} </tex>. Шаг выбирается независимо для каждого элемента матрицы.
 +
 
 +
Проанализируем зависимость погрешности метода от величины шага в случае использования правой разностной производной.
 +
Для сокращения записи введём обозначения <tex>f(\xi_j) = y_i(x_1 \dots \xi_j \dots x_m)</tex>.
 +
Остаточный член в соотношении <tex> f(x)' \approx \frac{f(x+ h) - f1}{h} </tex> имеет вид
 +
<tex>r_1 = - f''(\xi)h/2 </tex>.
 +
Если <tex> |f''(\xi)| \le M_2 </tex>, то <tex>|r_1| \le M_2*h/2 </tex>
 +
Если значения <tex>f1</tex> и <tex>f2</tex> заданы с погрешностями <tex>\vareps_1, \vareps_2 \le E</tex> , то погрешность будет
 +
содержать ещё одно слагаемое <tex> r_2 = - \frac{\vareps_1 - \vareps_2}{h} \le 2E/h </tex>. Таким образом, оценка суммарной погрешности имеет вид <tex>|r| \le |r_1| + |r_2| \le f(h) = M_2 h/2 + 2E/h </tex>. Эта оценка достигает минимума при
 +
<tex>h_0 = 2 \sqrt {E/M_2}</tex>. При этом <tex> f(h_0) = 2 \sqrt{M_2 E} </tex>.
 +
Оценка погрешности имеет один глобальный миниум. Поэтому выбор очень маленького шага не привидёт к росту точности. При величине шага, близкой к <tex>h_0</tex> погрешность имеет порядок <tex>O(\sqrt E) </tex>.
 +
 
 +
===Метод Бройдена===
 +
Чаще всего вычисление якобиана является одной из подзадач в различных методах оптимизации и решения систем нелинейных уравнений.
 +
При решении систем нелинейных уравнений методом Ньютона требуется вычислять якобиан на каждой итерации.
 +
Вычисление якобиана требует вычисления функций в <tex>O(m * n)</tex> точках. Это сложная и затратная по времени операция. Суть
 +
метода Бройдена состоит в том, чтобы вычислить якобиан аналитически или с помощью метода конечных разностей на первой итерации,
 +
а после этого на каждой итерации обновлять якобиан, не вычисляя значения функций <tex>y_1 \dots y_n</tex> и их производных.
 +
 
 +
Пусть задана система нелинейных уравнений <tex>F(x) = 0</tex>, где <tex> F:R^m \to R^n </tex>.
 +
Тогда якобиан на <tex>n</tex>-ой итерации выражается по формуле <br/>
 +
<tex>
 +
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
 +
</tex>
 +
 
 +
После этого следующее приближение вычисляется по формуле
 +
 
 +
<tex>x_{n + 1} = x_n - J_n^{-1} F(x_n) </tex>
 +
 
 +
== Методы вычисления матрицы Гессе ==
 +
Как и матрица Якоби, матрица Гессе может быть вычислена с помощью разностной аппроксимации производных.
 +
<tex> \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}</tex>, где
 +
<tex>x </tex> - вектор переменных, а <tex>e_i</tex> и <tex>e_j</tex> - единичные вектора.
 +
Эта формула требует вычисления значений функции <tex> f </tex> в <tex> \frac {n (n + 1)} {2} </tex> точках.
 +
Погрешность формулы имеет порядок <tex>O(h)</tex>.
 +
 
 +
 
 +
== Численный эксперимент==
 +
 
 +
В качестве примера рассчитаем с помощью вышеизложенного метода матрицу Гессе функции <tex>x^2 * \sqrt y </tex> в точке (1, 1)
 +
 
 +
{| class="standard"
 +
|-
 +
! Величина шага||Численный результат || Истинное значение || Погрешность
 +
|-
 +
|0.01
 +
||
 +
<tex>
 +
\begin{bmatrix} 1.999999999999780 & 1.002499984530392 \\ \\ 1.002499984530392 & -0.250007812910846 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 2 & 1 \\ \\ 1 & -0.25 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix}
 +
0.000000000000220 & 0.002499984530392 \\ \\ 0.002499984530392 & 0.500007812910846 \end{bmatrix}
 +
</tex>
 +
|-
 +
| 0.001||
 +
<tex>
 +
\begin{bmatrix} 1.999999999724444 & 1.000249999938419 \\1.000249999938419 & -0.250000078083623 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 2 & 1 \\ \\ 1 & -0.25 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix}
 +
0.000000000275556 & 0.000249999938418 \\ \\ 0.000249999938418 & 0.500000078083623 \end{bmatrix}
 +
</tex>
 +
|-
 +
| 0.0001||
 +
<tex>
 +
\begin{bmatrix} 1.999999998947288 & 1.000024996145044 \\ \\ 1.000024996145044 & -0.250000009582863 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 2 & 1 \\ \\ 1 & -0.25 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 0.000000001052712 & 0.000024996145044 \\ \\ 0.000024996145044 & 0.500000009582863 \end{bmatrix}
 +
</tex>
 +
|-
 +
|0.00001 ||
 +
<tex>
 +
\begin{bmatrix} 2.000002385926791 & 1.000002303186420 \\ \\ 1.000002303186420 & -0.250000020685093 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 2 & 1 \\ \\ 1 & -0.25 \end{bmatrix}
 +
</tex>
 +
||
 +
<tex>
 +
\begin{bmatrix} 0.000002385926791 & 0.000002303186420 \\ \\ 0.000002303186420 & 0.500000020685093 \end{bmatrix}
 +
</tex>
 +
 
 +
|}
-
== Числовой пример ==
 
== Рекомендации программисту ==
== Рекомендации программисту ==
 +
Вычисление матрицы Якоби требует большого числа операций. Поэтому при выборе метода её вычисления следует найти оптимальный баланс
 +
между погрешностью метода и его скоростью.
 +
 +
[[Media:Slimper_Jacobian.zip| Программа для вычисления якобиана, исходный код [1кб]]]
 +
 +
[[Media:Slimper_Hessian.zip| Программа для вычисления гессиана, исходный код [1кб] ]]
 +
 +
[[Media:Slimper_Ublas.zip| Библиотека алгоритмов линейной алгебры [180кб] ]]
 +
 +
== Заключение ==
== Заключение ==
== Список литературы ==
== Список литературы ==
-
* А. А. Самарский, А. В. Гулин. Численные методы. Москва «Наука», 1989.
+
* Самарский А. А. , Гулин А. В. Численные методы. Москва «Наука», 1989.
-
* Н. С. Бахвалов, Н. П. Жидков, Г. М. Кобельков. Численные методы. Лаборатория Базовых Знаний, 2003.
+
* Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. Численные методы. Лаборатория Базовых Знаний, 2003.
-
* Магнус Я. Р., Нейдеккер Х. Матричное дифференциальное исчисление с приложениями к статистике и эконометрике: Пер. с англ./ Под ред. С. А. Айвазяна.&nbsp;— М.:ФИЗМАТЛИТ, 2002.&nbsp;— 496&nbsp;с.
+
* Магнус Я. Р., Нейдеккер Х. Матричное дифференциальное исчисление с приложениями к статистике и эконометрике: Пер. с англ./ Под ред. С. А. Айвазяна.&nbsp;— М.:ФИЗМАТЛИТ, 2002.&nbsp;— 496&nbsp;с.
-
 
+
* Хайкин С. Нейронные сети, полный курс. 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
 +
* [http://www.inference.phy.cam.ac.uk/mackay/Bayes_FAQ.html#hessian MacKay D. How do you evaluate the Hessian in the first place?]
{{stub}}
{{stub}}
[[Категория:Численное дифференцирование]]
[[Категория:Численное дифференцирование]]
[[Категория:Учебные задачи]]
[[Категория:Учебные задачи]]
 +
[[Категория:Регрессионный анализ]]

Текущая версия

Содержание

Введение

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

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

Пусть задана система 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?
Личные инструменты