Участник:Василий Ломакин/Решение переопределенной СЛАУ

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

< Участник:Василий Ломакин(Различия между версиями)
Перейти к: навигация, поиск
м
 
(5 промежуточных версий не показаны.)
Строка 1: Строка 1:
== Постановка задачи ==
== Постановка задачи ==
-
Рассмотрим прямоугольную матрицу размером <tex></tex>:
+
Рассмотрим прямоугольную матрицу размером <tex>m \times n</tex>:
-
<p align="center"><tex>f{A}= \left( \begin{array}{cccc} a_{11} & a_{12} & \ldots & a_{1n}\\ a_{21} & a_{22} & \ldots & a_{2n}\\ \ldots & \ldots & \ldots & \ldots \\ a_{m1} & a_{m2} & \ldots & a_{mn}\\ \end{array}\right).</tex></p>
+
<p align="center"><tex>
 +
A= \left( \begin{array}{cccc} a_{11} & a_{12} & \ldots & a_{1n}\\ a_{21} & a_{22} & \ldots & a_{2n}\\ \ldots & \ldots & \ldots & \ldots \\ a_{m1} & a_{m2} & \ldots & a_{mn}\\ \end{array}\right).
 +
</tex></p>
 +
 
Пусть в матрице число строк превышает число столбцов (<tex>m > n</tex>), причём все строки линейно независимы.
Пусть в матрице число строк превышает число столбцов (<tex>m > n</tex>), причём все строки линейно независимы.
Систему уравнений вида
Систему уравнений вида
{{ eqno | 1 }}
{{ eqno | 1 }}
-
<p align="center"><tex>Ax=f</tex>,</p>
+
<p align="center"><tex>Au=f</tex>,</p>
-
где А - описанная выше, <tex>{u}={\{u_1, \ldots , u_n \}}^T</tex> — вектор-столбец решения, <tex>{f}={\{f_1, \ldots , f_n \}}^T</tex> — вектор-столбец правой части, назовём переопределённой. Как можно видеть, в такой системе число уравнений превышает число неизвестных, и для неё не существует "классического" решения, например методом Гаусса.
+
где <tex>A</tex> - описанная выше, <tex>{u}={\{u_1, \ldots , u_n \}}^T</tex> — вектор-столбец решения, <tex>{f}={\{f_1, \ldots , f_n \}}^T</tex> — вектор-столбец правой части, назовём '''переопределённой'''. Как можно видеть, в такой системе число уравнений превышает число неизвестных, и для неё не существует "классического" решения.
 +
 
== Изложение метода ==
== Изложение метода ==
 +
Приведем простой пример получения переопределённой системы линейных уравнений. Такого рода задачи часто встречаются, например, при обработке результатов экспериментов. Пусть <tex>f</tex> — линейная (или близкая к линейной) функция аргумента <tex>x:\ f(x) = u_1x + u_0</tex>. В точках <tex>x_k</tex> известны значения функции <tex>f(x_k)</tex>. Тогда <tex>u_0,\ u_1</tex> — коэффициенты, которые необходимо подобрать так, чтобы выполнялись условия <tex>u_1x_k + u_0 = f_k,\ k = 0,1,2,3,4,\ f_k = f(x_k)</tex>. Получим систему пяти уравнений относительно двух неизвестных. Это — переопределённая система. Она не имеет классического решения, так как в общем случае не существует прямой, проходящей через все 5 точек (это возможно только тогда, когда какие - либо три уравнения полученной системы линейными преобразованиями сводятся к двум другим — система линейно зависима). Необходимо провести аппроксимирующую кривую, которая не проходит через экспериментальные точки, но в то же время отражает экспериментальную зависимость, сглаживает возможные выбросы за счёт погрешности эксперимента.
 +
 +
Рассмотрим более общий случай. Пусть коэффициенты <tex>{u_0,\ u_1}</tex> необходимо определить по результатам <tex>n + 1</tex> измерения. Для каждого уравнения рассмотрим невязку <tex>r_k = u_1x_k + u_0 - f_k</tex> - разность левой и правой части. Невязка может принимать положительные и отрицательные значения. Чтобы не учитывать знаки, применим возведение в квадрат. Введем функцию, равную сумме квадратов невязок:
 +
{{ eqno | 2 }}
 +
<p align="center"><tex>\Phi (u_1,u_0) = \sum\limits_{k = 0}^n {r_k^2} = \sum\limits_{k = 0}^n {(u_1 x_k + u_0 - f_k)^2}</tex></p>
 +
Примем за обобщённое решение переопределённой СЛАУ такие <tex>{u_0,\ u_1}</tex>, для которых <tex>\Phi(u_0,\ u_1)</tex> принимает наименьшие значение. Для определения обобщённого решения из условия минимума суммы квадратов невязки получаем систему двух уравнений, уже имеющую классическое решение:
 +
<p align="center"><tex>\frac{\partial \Phi }{\partial u_0} = 0,\ \frac{\partial \Phi }{\partial u_1} = 0.</tex></p>
 +
 +
Рассмотрим теперь '''общий случай'''. Определим невязку <tex>r_k</tex> в виде
 +
<p align="center"><tex>r_k = \sum\limits_{j = 0}^p {u_j\varphi_j (x_k)} - f(x_k),\ k = 1, \ldots, n,</tex></p>
 +
где <tex>\varphi_j (x)</tex> — некоторые функции, образующие базис, например, тригонометрические: <tex>\varphi_j (x) = \sin (jx)</tex> . Выражение <tex>\sum\limits_{j = 0}^p {u_j \varphi_j (x)}</tex> называется обобщённым полиномом. В приведенном выше примере в качестве базисных функций были выбраны степенные функции <tex>\varphi_j (x) = x^j</tex> . Обобщённый полином превратился в алгебраический.
 +
 +
В случае выбора произвольной системы базисных функций переопределенная СЛАУ и функционал <tex>\Phi(u_0, \dots, u_p)</tex> будут иметь вид
 +
<p align="center"><tex>\begin{gather*} u_0 \varphi_0 (x_0) + \ldots + u_p \varphi_p(x_0) = f_0, \\ \ldots \\ u_0 \varphi_0 (x_n) + \ldots + u_p\varphi_p (x_n) = f_n,\\ \Phi (u_0,\ldots,u_n) = \sum\limits_{i = 0}^n (\sum\limits_{j = 0}^p u_j \varphi_j(x_i) - f_i)^2 \end{gather*}</tex></p>
 +
 +
Отыщем обобщенное решение методом наименьших квадратов: приравняем все частные производные по компонентам обобщенного решения к нулю (условия минимума):
 +
<p align="center"><tex>\frac{\partial \Phi }{\partial u_k} = 2 \sum\limits_{i = 0}^n {\varphi_k(x_i)\left({\sum\limits_{j = 0}^p {u_j\varphi_j (x_i) - f_i}} \right)}=0</tex></p>
 +
и изменяя порядок суммирования, получаем СЛАУ:
 +
<p align="center"><tex>\sum\limits_{j = 0}^p {\left({\sum\limits_{i = 0}^n{\varphi_j (x_i)\varphi_k (x_i)}}\right)u_j = \sum\limits_{i = 0}^n {f_i\varphi_k (x_i)}},\ k = 0, \ldots, p,</tex></p>
 +
или, в виде системы уравнений:
 +
<p align="center"><tex>\begin{gather*} (\varphi_0, \varphi_0) u_0 + (\varphi_0, \varphi_1)u_1 + \ldots + (\varphi_0, \varphi_p)u_p = (\varphi_0, f), \\ (\varphi_1, \varphi_0) u_0 + (\varphi_1, \varphi_1)u_1 + \ldots + (\varphi_1, \varphi_p)u_p = (\varphi_1, f), \\ \ldots \\ (\varphi_p, \varphi_0) u_0 + (\varphi_p, \varphi_1)u_1 + \ldots + (\varphi_p, \varphi_p)u_p = (\varphi_p, f), \end{gather*}</tex></p>
 +
 +
Система метода наименьших квадратов имеет вид <tex> \mathbf{Du} = \mathbf{f}</tex> с матрицей <tex>\mathbf{D}</tex>, элементами которой являются скалярные произведения <tex>(\varphi_i, \varphi_j) = \sum\limits_{i = 0}^n \varphi_j (x) \varphi_k (x_i)</tex>:
 +
 +
<p align="center"><tex>
 +
\begin{vmatrix}
 +
(\varphi_0,\ \varphi_0) & (\varphi_0,\ \varphi_1) & \ldots & (\varphi_0,\ \varphi_p)\\
 +
(\varphi_1,\ \varphi_0) & (\varphi_1,\ \varphi_1) & \ldots & (\varphi_1,\ \varphi_p)\\
 +
\ldots & \ldots & \ldots & \ldots\\
 +
(\varphi_p,\ \varphi_0) & (\varphi_p,\ \varphi_1) & \ldots & (\varphi_p,\ \varphi_p)\\
 +
\end{vmatrix}.
 +
</tex></p>
 +
Это — матрица Грама. В правой части системы стоят проекции свободного члена исходной задачи на подпространство базисных функций <tex>(\varphi,f) = \sum\limits_{i = 0}^n {\varphi_j(x_i)f_i}</tex>. Матрица симметричная и положительно определенная, таким образом, решение исследуемой СЛАУ существует и единственно. Находится, например, с помощью итерационного метода Гаусса.
-
== Анализ метода и оценка ошибок ==
+
== Советы программисту ==
-
== Числовой пример ==
+
Отметим основные свойства матрицы Грама, полезные при программировании:
 +
# Матрица симметрична, что позволяет сократить объём вычислений при заполнении матрицы.
 +
# Матрица является положительно определённой, следовательно, при решении системы нормальных уравнений методом исключения Гаусса можно отказаться от процедуры выбора главного элемента.
 +
# Определитель матрицы будет отличен от нуля, если в качестве базиса выбраны линейно независимые функции, при этом система имеет единственное решение.
 +
# При обработке экспериментальных данных, определённых с погрешностью <tex>\varepsilon</tex> в каждой точке, обычно сначала берут небольшое (одну-две) число базисных функций. После вычисления приближённого решения, вычисляют сумму квадратов невязок по формуле, аналогичной {{eqref|2}}. Если получается, что <tex>sqrt{\Phi} > \varepsilon</tex>, то необходимо расширить базис добавлением новых функций. Расширение базиса необходимо производить до тех пор, пока не выполнится условие <tex>\sqrt{\Phi} \approx \varepsilon</tex>.
 +
# Выбор конкретных базисных функций зависит от свойств функции <tex>f</tex>, таких как периодичность, экспоненциальный или логарифмический характер, свойства симметрии, наличие асимптотики и т.д.
== Список литературы ==
== Список литературы ==
-
* ''Н.Н.Калиткин.'' Численные методы М.: Наука, 1978.
+
* ''А.Е.Мудров'' Численные методы для ПЭВМ. Томск: Раско, 1991.
* ''А.А.Самарский, А.В.Гулин.''&nbsp; Численные методы М.: Наука, 1989.
* ''А.А.Самарский, А.В.Гулин.''&nbsp; Численные методы М.: Наука, 1989.
== См. также ==
== См. также ==
* [[Практикум ММП ВМК, 4й курс, осень 2008|Практикум ММП ВМК, 4й курс, осень 2008]]
* [[Практикум ММП ВМК, 4й курс, осень 2008|Практикум ММП ВМК, 4й курс, осень 2008]]
 +
* [[Метод наименьших квадратов|Метод наименьших квадратов]]

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

Содержание

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

Рассмотрим прямоугольную матрицу размером m \times n:

A= \left( \begin{array}{cccc} a_{11} & a_{12} & \ldots & a_{1n}\\ a_{21} & a_{22} & \ldots & a_{2n}\\ \ldots & \ldots & \ldots & \ldots \\ a_{m1} & a_{m2} & \ldots & a_{mn}\\ \end{array}\right).

Пусть в матрице число строк превышает число столбцов (m > n), причём все строки линейно независимы. Систему уравнений вида

( 1 )

Au=f,

где A - описанная выше, {u}={\{u_1, \ldots , u_n \}}^T — вектор-столбец решения, {f}={\{f_1, \ldots , f_n \}}^T — вектор-столбец правой части, назовём переопределённой. Как можно видеть, в такой системе число уравнений превышает число неизвестных, и для неё не существует "классического" решения.

Изложение метода

Приведем простой пример получения переопределённой системы линейных уравнений. Такого рода задачи часто встречаются, например, при обработке результатов экспериментов. Пусть f — линейная (или близкая к линейной) функция аргумента x:\ f(x) = u_1x + u_0. В точках x_k известны значения функции f(x_k). Тогда u_0,\ u_1 — коэффициенты, которые необходимо подобрать так, чтобы выполнялись условия u_1x_k + u_0 = f_k,\ k = 0,1,2,3,4,\ f_k = f(x_k). Получим систему пяти уравнений относительно двух неизвестных. Это — переопределённая система. Она не имеет классического решения, так как в общем случае не существует прямой, проходящей через все 5 точек (это возможно только тогда, когда какие - либо три уравнения полученной системы линейными преобразованиями сводятся к двум другим — система линейно зависима). Необходимо провести аппроксимирующую кривую, которая не проходит через экспериментальные точки, но в то же время отражает экспериментальную зависимость, сглаживает возможные выбросы за счёт погрешности эксперимента.

Рассмотрим более общий случай. Пусть коэффициенты {u_0,\ u_1} необходимо определить по результатам n + 1 измерения. Для каждого уравнения рассмотрим невязку r_k = u_1x_k + u_0 - f_k - разность левой и правой части. Невязка может принимать положительные и отрицательные значения. Чтобы не учитывать знаки, применим возведение в квадрат. Введем функцию, равную сумме квадратов невязок:

( 2 )

\Phi (u_1,u_0) = \sum\limits_{k = 0}^n {r_k^2} = \sum\limits_{k = 0}^n {(u_1 x_k + u_0 - f_k)^2}

Примем за обобщённое решение переопределённой СЛАУ такие {u_0,\ u_1}, для которых \Phi(u_0,\ u_1) принимает наименьшие значение. Для определения обобщённого решения из условия минимума суммы квадратов невязки получаем систему двух уравнений, уже имеющую классическое решение:

\frac{\partial \Phi }{\partial u_0} = 0,\ \frac{\partial \Phi }{\partial u_1} = 0.

Рассмотрим теперь общий случай. Определим невязку r_k в виде

r_k = \sum\limits_{j = 0}^p {u_j\varphi_j (x_k)} - f(x_k),\ k = 1, \ldots, n,

где \varphi_j (x) — некоторые функции, образующие базис, например, тригонометрические: \varphi_j (x) = \sin (jx) . Выражение \sum\limits_{j = 0}^p {u_j \varphi_j (x)} называется обобщённым полиномом. В приведенном выше примере в качестве базисных функций были выбраны степенные функции \varphi_j (x) = x^j . Обобщённый полином превратился в алгебраический.

В случае выбора произвольной системы базисных функций переопределенная СЛАУ и функционал \Phi(u_0, \dots, u_p) будут иметь вид

\begin{gather*} u_0 \varphi_0 (x_0) + \ldots + u_p \varphi_p(x_0) = f_0, \\ \ldots \\ u_0 \varphi_0 (x_n) + \ldots + u_p\varphi_p (x_n) = f_n,\\ \Phi (u_0,\ldots,u_n) = \sum\limits_{i = 0}^n (\sum\limits_{j = 0}^p u_j \varphi_j(x_i) - f_i)^2 \end{gather*}

Отыщем обобщенное решение методом наименьших квадратов: приравняем все частные производные по компонентам обобщенного решения к нулю (условия минимума):

\frac{\partial \Phi }{\partial u_k} = 2 \sum\limits_{i = 0}^n {\varphi_k(x_i)\left({\sum\limits_{j = 0}^p {u_j\varphi_j (x_i) - f_i}} \right)}=0

и изменяя порядок суммирования, получаем СЛАУ:

\sum\limits_{j = 0}^p {\left({\sum\limits_{i = 0}^n{\varphi_j (x_i)\varphi_k (x_i)}}\right)u_j = \sum\limits_{i = 0}^n {f_i\varphi_k (x_i)}},\ k = 0, \ldots, p,

или, в виде системы уравнений:

\begin{gather*} (\varphi_0, \varphi_0) u_0 + (\varphi_0, \varphi_1)u_1 + \ldots + (\varphi_0, \varphi_p)u_p = (\varphi_0, f), \\ (\varphi_1, \varphi_0) u_0 + (\varphi_1, \varphi_1)u_1 + \ldots + (\varphi_1, \varphi_p)u_p = (\varphi_1, f), \\ \ldots \\ (\varphi_p, \varphi_0) u_0 + (\varphi_p, \varphi_1)u_1 + \ldots + (\varphi_p, \varphi_p)u_p = (\varphi_p, f), \end{gather*}

Система метода наименьших квадратов имеет вид \mathbf{Du} = \mathbf{f} с матрицей \mathbf{D}, элементами которой являются скалярные произведения (\varphi_i, \varphi_j) = \sum\limits_{i = 0}^n \varphi_j (x) \varphi_k (x_i):

\begin{vmatrix} (\varphi_0,\ \varphi_0) & (\varphi_0,\ \varphi_1) & \ldots & (\varphi_0,\ \varphi_p)\\ (\varphi_1,\ \varphi_0) & (\varphi_1,\ \varphi_1) & \ldots & (\varphi_1,\ \varphi_p)\\ \ldots & \ldots & \ldots & \ldots\\ (\varphi_p,\ \varphi_0) & (\varphi_p,\ \varphi_1) & \ldots & (\varphi_p,\ \varphi_p)\\ \end{vmatrix}.

Это — матрица Грама. В правой части системы стоят проекции свободного члена исходной задачи на подпространство базисных функций (\varphi,f) = \sum\limits_{i = 0}^n {\varphi_j(x_i)f_i}. Матрица симметричная и положительно определенная, таким образом, решение исследуемой СЛАУ существует и единственно. Находится, например, с помощью итерационного метода Гаусса.

Советы программисту

Отметим основные свойства матрицы Грама, полезные при программировании:

  1. Матрица симметрична, что позволяет сократить объём вычислений при заполнении матрицы.
  2. Матрица является положительно определённой, следовательно, при решении системы нормальных уравнений методом исключения Гаусса можно отказаться от процедуры выбора главного элемента.
  3. Определитель матрицы будет отличен от нуля, если в качестве базиса выбраны линейно независимые функции, при этом система имеет единственное решение.
  4. При обработке экспериментальных данных, определённых с погрешностью \varepsilon в каждой точке, обычно сначала берут небольшое (одну-две) число базисных функций. После вычисления приближённого решения, вычисляют сумму квадратов невязок по формуле, аналогичной (2). Если получается, что sqrt{\Phi} > \varepsilon, то необходимо расширить базис добавлением новых функций. Расширение базиса необходимо производить до тех пор, пока не выполнится условие \sqrt{\Phi} \approx \varepsilon.
  5. Выбор конкретных базисных функций зависит от свойств функции f, таких как периодичность, экспоненциальный или логарифмический характер, свойства симметрии, наличие асимптотики и т.д.

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

  • А.Е.Мудров Численные методы для ПЭВМ. Томск: Раско, 1991.
  • А.А.Самарский, А.В.Гулин.  Численные методы М.: Наука, 1989.

См. также

Источник — «http://machinelearning.ru/wiki/index.php?title=%D0%A3%D1%87%D0%B0%D1%81%D1%82%D0%BD%D0%B8%D0%BA:%D0%92%D0%B0%D1%81%D0%B8%D0%BB%D0%B8%D0%B9_%D0%9B%D0%BE%D0%BC%D0%B0%D0%BA%D0%B8%D0%BD/%D0%A0%D0%B5%D1%88%D0%B5%D0%BD%D0%B8%D0%B5_%D0%BF%D0%B5%D1%80%D0%B5%D0%BE%D0%BF%D1%80%D0%B5%D0%B4%D0%B5%D0%BB%D0%B5%D0%BD%D0%BD%D0%BE%D0%B9_%D0%A1%D0%9B%D0%90%D0%A3»
Личные инструменты