Викия

Математика

Метод наименьших квадратов

1457статей на
этой вики
Добавить новую страницу
Обсуждение3 Поделиться

Обнаружено использование расширения AdBlock.


Викия — это свободный ресурс, который существует и развивается за счёт рекламы. Для блокирующих рекламу пользователей мы предоставляем модифицированную версию сайта.

Викия не будет доступна для последующих модификаций. Если вы желаете продолжать работать со страницей, то, пожалуйста, отключите расширение для блокировки рекламы.

Метод наименьших квадратов(МНК) - математический метод, по нахождению приближающей функции - f(x) по набору данных(точек) \{ (x_1, y_1), (x_2, y_2),\ldots, (x_n, y_n)\}, которая минимизирует сумму квадратов отклонений точек от найденной функции. 
6990 5.gif

Приближающая функция - f(x)

Сущность МНК Править

Требуется по значениям x_1, ...,x_n и y_1, ..., y_nкак можно точнее оценить f(x), точнее означает с минимальными ошибками, с минимальным разбросом(дисперсией), однако зависимость f(x) не должна быть просто построенной по точкам (x_i, y_i).
Из Модели Регрессии следует минимизировать дисперсию ошибок:

D[\boldsymbol{\varepsilon}] = \frac{1}{n}\sum_{i=1}^{n}(y_i - f(\vec x_i, \vec c))^2

где \vec c=\{c_0,..,c_m\} — вектор неизвестных параметров.

Далее будем рассматривать только случай однопараметрического (однофакторного) МНК, когда приближающая функция f(\vec x_i, \vec c) зависит только от одного параметра x, (т.е приближающая функция зависит только от одной переменной), поэтому считаем что дисперсия имеет вид:

D[\boldsymbol{\varepsilon}] = \frac{1}{n}\sum_{i=1}^{n}(y_i - f(x_i, \vec c))^2

Минимум D[\boldsymbol{\varepsilon}] будем искать как:

\frac{\partial D[\boldsymbol{\varepsilon}]}{\partial \vec c} = 0

\frac{\partial D[\boldsymbol{\varepsilon}]}{\partial \vec c} = \frac{\partial D[\boldsymbol{\varepsilon}]}{\partial c_0} \vec e_0 + \ldots + \frac{\partial D[\boldsymbol{\varepsilon}]}{\partial c_m} \vec e_m = 0

где \{\vec e_k\}_{k=0}^m = \{\vec e_0, .., \vec e_m\} — базис из линейно-независимых функций.

Минимизируя дисперсию ошибок D[\boldsymbol{\varepsilon}] по неизвестным параметрам \{c_0,..,c_m\}, базисом может являться система линейно-независимых функций \{\vec e_k\}_{k=0}^m = \{\varphi_k\}_{k=0}^m = \{\varphi_0(x), .., \varphi_m(x)\}, при условии m < n, тогда приближающая функция есть разложение по этому базису: f(x, \vec c) = \sum_{k=0}^m c_k \varphi_k(x), а функция дисперсии ошибки будет иметь вид:

D[\boldsymbol{\varepsilon}] = \frac{1}{n}\sum_{i=1}^{n}(y_i - \sum_{k=0}^m c_k \varphi_k(x_i) )^2

Общий случай вычисления коэффициентов приближающей функции Править


Для этого дифференцируем D[\boldsymbol{\varepsilon}](переобозначение: D[\boldsymbol{\varepsilon}] \equiv D) отдельно по каждому из параметров c_j из \{c_k\}_{k=0}^m, где j пробегает значения от 0 до m:
\frac{\partial D}{\partial c_j } = 0 , где j=0..m
Получаем систему уравнений относительно параметров \{c_k\}_{k=0}^m:

 \begin{cases}
   \frac{\partial D}{\partial c_0 } = 0
   \\
   \frac{\partial D}{\partial c_1 } = 0
   \\
   \ldots
   \\
   \frac{\partial D}{\partial c_m } = 0
 \end{cases}

\frac{\partial D}{\partial c_j }= \frac{\partial}{\partial c_j } \sum_{i=1}^n (y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}))^2 = -2 \sum_{i=1}^n (y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}) ) \sum_{k=0}^m \frac{\partial c_k}{\partial c_j} \varphi_k(x_{i}) = 0

\frac{\partial c_k}{\partial c_j} =  \begin{cases}
1,  j = k \\
0,  j \ne k
 \end{cases}\equiv \delta_{jk}

-2 \sum_{i=1}^n (y_i - \sum_{k=0}^m c_k \varphi_k(x_{i})) \sum_{k=0}^m \delta_{jk} \varphi_k(x_{i}) = 0

Используя:

\sum_{k=0}^m \delta_{jk} \varphi_k(x) = \varphi_j(x)

приходим к:

\sum_{i=1}^n [y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}) ] \varphi_j(x_{i}) = 0,   где j=0 \ldots m


Распишем предыдущую формулу в виде системы m+1 уравнений:


\begin{cases}
    \sum_{i=1}^n [y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}) ] \varphi_0(x_{i}) = 0 \\
    \sum_{i=1}^n [y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}) ] \varphi_1(x_{i}) = 0\\
    \dots\\
    \sum_{i=1}^n [y_i - \sum_{k=0}^m c_k \varphi_k(x_{i}) ] \varphi_m(x_{i}) = 0 \\
\end{cases}

Процесс приведения системы уравнений к матричному виду:
Раскроем скобки и перенесём y_i \varphi_0(x_{i}) вправо:


\begin{cases}
    \sum_{i=1}^n \varphi_0(x_{i}) \sum_{k=0}^m c_k \varphi_k(x_{i}) = \sum_{i=1}^n y_i \varphi_0(x_{i}) \\
    \sum_{i=1}^n \varphi_1(x_{i}) \sum_{k=0}^m c_k \varphi_k(x_{i}) = \sum_{i=1}^n y_i \varphi_1(x_{i}) \\
    \dots\\
    \sum_{i=1}^n \varphi_m(x_{i}) \sum_{k=0}^m c_k \varphi_k(x_{i}) = \sum_{i=1}^n y_i  \varphi_m(x_{i}) \\
\end{cases}

В левой части, раскроем сумму по k и перемножим:


\begin{cases}
    \sum_{i=1}^n [c_0 \varphi_0^2(x_{i}) + c_1 \varphi_0(x_{i}) \varphi_1(x_{i}) + ... + c_m \varphi_0(x_{i}) \varphi_m(x_{i})]  = \sum_{i=1}^n y_i \varphi_0(x_{i})\\
    \sum_{i=1}^n [c_0 \varphi_1(x_{i}) \varphi_0(x_{i}) + c_1 \varphi_1^2(x_{i}) + ... + c_m \varphi_1(x_{i}) \varphi_m(x_{i})]  = \sum_{i=1}^n y_i \varphi_1(x_{i}) \\
    \dots\\
    \sum_{i=1}^n [c_0 \varphi_m(x_{i}) \varphi_0(x_{i}) + c_1 \varphi_m(x_{i}) \varphi_1(x_{i}) + ... + c_m \varphi_m^2(x_{i})]  = \sum_{i=1}^n y_i \varphi_m(x_{i}) \\
\end{cases}

Используя следующие переобозначения:

\langle\varphi_j, \varphi_k\rangle = \sum_{i=1}^n \varphi_j(x_i) \varphi_k(x_i),  

\langle y, \varphi_j\rangle = \sum_{i=1}^n y_i \varphi_j(x_i)

Запишем предыдущую систему уравнений в матричном виде: A \vec c = \vec b,

где 
A = \begin{pmatrix} 
 \langle\varphi_0, \varphi_0\rangle & \langle\varphi_1, \varphi_0\rangle & \cdots & \langle\varphi_m, \varphi_0\rangle \\
 \langle\varphi_0, \varphi_1\rangle & \langle\varphi_1, \varphi_1\rangle & \cdots & \langle\varphi_m, \varphi_1\rangle \\ 
 \vdots & \vdots & \ddots & \vdots \\
 \langle\varphi_0, \varphi_m\rangle & \langle\varphi_1, \varphi_m\rangle & \cdots & \langle\varphi_m, \varphi_m\rangle \\ 
\end{pmatrix}, 
\vec c = \begin{pmatrix} c_0 \\
 c_1 \\ \vdots \\ c_m 
\end{pmatrix}, \vec b = \begin{pmatrix} 
 \langle y, \varphi_0\rangle \\
 \langle y, \varphi_1\rangle \\
 \vdots \\ 
 \langle y, \varphi_m\rangle \\
\end{pmatrix}

Матрица A - называется матрицей Грама.

решаем систему матричным методом и находим вектор искомых коэффициентов приближающей функции: \vec c = \{c_0\;, c_1,\;\ldots,\;c_m\;\}\vec c = \begin{pmatrix} c_0 \\
 c_1 \\ \vdots \\ c_m \end{pmatrix} = A^{-1}b = \begin{pmatrix} 
 \langle\varphi_0, \varphi_0\rangle & \langle\varphi_1, \varphi_0\rangle & \cdots & \langle\varphi_m, \varphi_0\rangle \\
 \langle\varphi_0, \varphi_1\rangle & \langle\varphi_1, \varphi_1\rangle & \cdots & \langle\varphi_m, \varphi_1\rangle \\ 
 \vdots & \vdots & \ddots & \vdots \\
 \langle\varphi_0, \varphi_m\rangle & \langle\varphi_1, \varphi_m\rangle & \cdots & \langle\varphi_m, \varphi_m\rangle \\ 
\end{pmatrix}^{-1} \begin{pmatrix} 
 \langle y, \varphi_0\rangle \\
 \langle y, \varphi_1\rangle \\
 \vdots \\ 
 \langle y, \varphi_m\rangle \\
\end{pmatrix}

В результате, мы нашли приближающую функцию f(x):

f(x) = \begin{pmatrix} c_0 \\
 c_1 \\ \vdots \\ c_m \end{pmatrix} \begin{pmatrix} \varphi_0(x), & \varphi_1(x), & \ldots, & \varphi_m(x) \end{pmatrix}  = c_0 \varphi_0(x) + c_1 \varphi_1(x) + ... + c_m \varphi_m(x)

Основные случаи Править

Степенное разложение(общий случай) Править

Часто в виде системы линейно-независимых функций \{\varphi_k(x)\}_{k=0}^m выбирают набор степенных функций \{1, x, .., x^m\} и приближающую функцию ищут в виде:

f(x, \vec c) = c_0 + c_1 x + c_2 x^2 + ... + c_m x^m = \sum_{k=0}^m c_k x^k

Тогда система уравнений A\vec c = \vec b будет иметь следующий вид:

(далее производится суммирование \sum_{i=1}^n которое просто обозначаем \sum )

где A = \begin{pmatrix} \sum 1 & \sum x_i & \cdots & \sum x_i^m \\
 \sum x_i & \sum x_i^2 & \cdots & \sum x_i^{m+1} \\ 
 \vdots & \vdots & \ddots & \vdots \\
 \sum x_i^m & \sum x_i^{m+1} & \cdots & \sum x_i^{2m}
\end{pmatrix} ,      b = \begin{pmatrix} \sum y_i \\
 \sum y_i x_i \\ \vdots \\ \sum y_i x_i^m \end{pmatrix}

Вектор коэффициентов: \vec c = \begin{pmatrix} c_0 \\
 c_1 \\ \vdots \\ c_m \end{pmatrix} = A^{-1}b = \begin{pmatrix} \sum 1 & \sum x_i & \cdots & \sum x_i^m \\
 \sum x_i & \sum x_i^2 & \cdots & \sum x_i^{m+1} \\ 
 \vdots & \vdots & \ddots & \vdots \\
 \sum x_i^m & \sum x_i^{m+1} & \cdots & \sum x_i^{2m}
\end{pmatrix}^{-1}  \begin{pmatrix} \sum y_i \\
 \sum y_i x_i \\ \vdots \\ \sum y_i x_i^m \end{pmatrix}

f(x) = \begin{pmatrix} c_0 \\
 c_1 \\ \vdots \\ c_m \end{pmatrix} \{ 1, x, \ldots, x^m \}  = c_0 + c_1x + c_2x^2 + ... + c_mx^m

На практике применима до степени m \leq 5 из-за быстрого нарастания ошибки при вычислении матрицы A^{-1}.

Приближение линейной зависимостью Править

Набор данных(точки) \{x_i, y_i\}_{i=1}^n, приближается линейной зависимостью(прямой):

f(x) = c_0 + c_1 x

Базисные функции линейной зависимости: \varphi_0(x) = 1, \ \varphi_1(x) = x

(далее производится суммирование \sum_{i=1}^n)


A = \begin{pmatrix} 
\langle\varphi_0, \varphi_0\rangle & \langle\varphi_1, \varphi_0\rangle \\
\langle\varphi_0, \varphi_1\rangle & \langle\varphi_1, \varphi_1\rangle \\
\end{pmatrix} 
= 
\begin{pmatrix} 
 \sum 1 & \sum x_i \\
 \sum x_i & \sum x_i^2 \\
\end{pmatrix}
,

 
b = \begin{pmatrix} 
 \langle y, \varphi_0\rangle \\
 \langle y, \varphi_1\rangle 
\end{pmatrix} = \begin{pmatrix}  
 \sum y_i \\
 \sum y_i x_i 
\end{pmatrix}

\sum_{i=1}^n 1 = n

\vec c = A^{-1} \vec b = \begin{pmatrix} 
 n & \sum x_i \\
 \sum x_i & \sum x_i^2 \\
\end{pmatrix}^{-1} \cdot \begin{pmatrix}  
 \sum y_i \\
 \sum y_i x_i 
\end{pmatrix}

A^{-1} = \frac{1}{n \sum x_i^2 - (\sum x_i)^2} \begin{pmatrix}
 \sum x_i^2 & -\sum x_i \\
 -\sum x_i  & n \\
\end{pmatrix}

\vec c = \begin{pmatrix} c_0 \\ c_1 \end{pmatrix} = 
A^{-1} \vec b =  \frac{1}{n \sum x_i^2 - (\sum x_i)^2} \begin{pmatrix}
 \sum x_i^2 & -\sum x_i \\
 -\sum x_i  & n \\
\end{pmatrix} \cdot \begin{pmatrix}  
 \sum y_i \\
 \sum y_i x_i 
\end{pmatrix}

\vec c = \begin{pmatrix} c_0 \\ c_1 \end{pmatrix} = 
\frac{1}{n \sum x_i^2 - (\sum x_i)^2} 
\begin{pmatrix}
 \sum x_i^2 \sum y_i - \sum x_i \sum y_i x_i \\ 
 n \sum y_i x_i  -  \sum x_i \sum y_i
\end{pmatrix}

Таким образом вычислили коэффициенты и нашли приближающую функцию:

c_0 = \frac{\sum x_i^2 \sum y_i - \sum x_i \sum y_i x_i}{n \sum x_i^2 - (\sum x_i)^2}

 c_1 = \frac{ n \sum y_i x_i  -  \sum x_i \sum y_i }{ n \sum x_i^2 - ( \sum x_i )^2 }

f(x) = c_0 + c_1 x

Ортогональность линейно-независимой системы функций Править

Чтобы упростить решение системы уравнений A \vec c = \vec b находя \vec c матричным методом \vec c = A^{-1} \vec b, нужно матрицу A привести к диагональному виду, для этого необходимо чтобы базисные функции \{\varphi_k(x)\}_{k=0}^m по которым разложена приближающая функция f(x) = \sum_{k=0}^m c_k \varphi_k(x) были ортогональны на определённом интервале, например на [x_1, x_n].

То есть их скалярное произведение: \langle\varphi_j, \varphi_k\rangle = \sum_{i=1}^n \varphi_j(x_i) \varphi_k(x_i)

обладало бы свойством ортогональности: \langle\varphi_j, \varphi_k\rangle =  \begin{cases}
\neq 0, \ j = k \\
\ \ 0, \ j \ne k
\end{cases} на интервале [x_1, x_n]

Тогда матрица Грама станет диагональной:


A = \begin{pmatrix}
 \langle\varphi_0, \varphi_0\rangle & 0 & \cdots & 0  \\
 0 & \langle\varphi_1, \varphi_1\rangle & \cdots & 0  \\ 
 \vdots & \vdots & \ddots & \vdots \\
 0 & 0 & \cdots & \langle\varphi_m, \varphi_m\rangle \\ 
\end{pmatrix}

а набор коэффициентов \vec c = \{c_k\}_{k=0}^m может быть легко вычислен:

\vec c = A^{-1} \vec b = \begin{pmatrix}
 \frac{1}{\langle\varphi_0, \varphi_0\rangle} & 0 & \cdots & 0  \\
 0 & \frac{1}{\langle\varphi_1, \varphi_1\rangle} & \cdots & 0  \\ 
 \vdots & \vdots & \ddots & \vdots \\
 0 & 0 & \cdots & \frac{1}{\langle\varphi_m, \varphi_m\rangle} \\ 
\end{pmatrix} \cdot \begin{pmatrix} 
 \langle y, \varphi_0\rangle \\
 \langle y, \varphi_1\rangle \\
 \vdots \\ 
 \langle y, \varphi_m\rangle \\
\end{pmatrix} = \begin{pmatrix} 
 \frac{\langle y, \varphi_0\rangle}{\langle\varphi_0, \varphi_0\rangle} \\
 \frac{\langle y, \varphi_1\rangle}{\langle\varphi_1, \varphi_1\rangle} \\
 \vdots \\ 
 \frac{\langle y, \varphi_m\rangle}{\langle\varphi_m, \varphi_m\rangle} \\
\end{pmatrix}

В данном случае, коэффициенты \{c_k\}_{k=0}^m приближающей функции f(x) называется коэффициентами Фурье

c_k = \frac{\langle y, \varphi_k\rangle}{\langle\varphi_k, \varphi_k\rangle}

А сама приближающая функция f(x) - обобщённым многочленом Фурье:

f(x) = \sum_{k=0}^m \frac{\langle y, \varphi_k\rangle}{\langle\varphi_k, \varphi_k\rangle} \varphi_k(x)


Однако, ортогональные функции, назовём их \{e_k\}_{k=0}^m (чтобы отличать от первоначально заданного линейно-независимых функций \{\varphi_k\}_{k=0}^m) заранее неизвестны и определяются относительно того или иного скалярного произведения \langle\varphi_k, \varphi_j\rangle.

Находятся \{e_k\}_{k=0}^m методом ортогонализации линейной-независимого набора векторов \{\varphi_k\}_{k=0}^m (процесс Грама-Шмидта):

e_0 = \varphi_0

e_1 = \varphi_1 - \frac{\langle \varphi_1, e_0 \rangle}{\|e_0\|^2} e_0

\ldots

e_m = \varphi_m - \sum_{k=0}^{m-1} \frac{\langle \varphi_m, e_k \rangle}{\|e_k\|^2} e_k

где \| e_k \|^2 = \langle e_k, e_k \rangle

Таким образом мы получили приближающую функцию наилучшего среднеквадратичного приближения, в виде разложения в ряд Фурье(обобщённого многочлена):

f(x) = \sum_{k=0}^m c_k e_k = \sum_{k=0}^m \frac{\langle y, e_k\rangle}{\| e_k \|^2} e_k

где c_k = \frac{\langle y, e_k\rangle}{\langle e_k, e_k\rangle} = \frac{\langle y, e_k\rangle}{\| e_k \|^2}

Однако, такой метод нахождение базисных функции \{e_k\}_{k=0}^m является численно неустойчивым, так как происходит накопление ошибки включающей в себя вычисления предыдущих функций.

Для того чтобы избежать накопления ошибки при вычислениях базисных функций \{e_k\}_{k=0}^m, нужно воспользоваться рекуррентным соотношением ортогональных многочленов:

{e_{k+1}\ =\ (\varphi_1-\beta_k)\ e_k\ -\ \alpha_k\ e_{k-1}},

где

\beta_n = \frac{\langle \varphi_1 e_k, e_k\rangle}{\langle e_k, e_k \rangle}, \qquad 	\alpha_k = \frac{\langle \varphi_1 e_k, e_{k-1}\rangle}{\langle e_{k-1}, e_{k-1}\rangle}.

Для случая когда приближение производится степенными функциями, т.е. \{\varphi_k\}_{k=0}^m = \{x^k\}_{k=0}^m рекуррентное соотношение имеет вид:


e_{k+1} = x e_k - \frac{\langle x e_k, e_k \rangle}{\langle e_k, e_k \rangle} e_k - \frac{\langle x e_k, e_{k-1} \rangle}{\langle e_{k-1}, e_{k-1} \rangle} e_{k-1}

где

k = 1,2,..; \ \ \ \  e_0 = 1, \ \ \ \    e_1 = x - \frac{\langle x, e_0 \rangle}{\langle e_0, e_0 \rangle} e_0

Иногда приведённые выше многочлены носят названия многочленов Чебышева, но не стоит их путать с классическими многочленами Чебышева.

Рекуррентное соотношение может быть упрощено[1], если область определения функции y, [x_1, x_n] отобразить на интервал симметричный относительно нуля [-1, 1]:

t:[x_1, x_n] \rightarrow [-1, 1]

Примером отображения t , для равноотстоящих точек, может быть функция:

t(x) = \frac{2x - (x_n + x_1)}{x_n - x_1}

Тогда рекуррентное соотношение будет иметь более простой вид [1]:

e_{k+1}(t)\ =\ t\ e_k(t)\ -\ \alpha_k\ e_{k-1}(t)

где

e_0 = 1; \ \ \ \ e_1 = t; \ \ \ \ \ \alpha_k = \frac{\langle t e_k, e_{k-1}\rangle}{\langle e_{k-1}, e_{k-1}\rangle};


Ссылки Править

  1. 1,0 1,1 "Вывод рекуррентного соотношения ортогональных многочленов из процесса ортогонализации Грама-Шмидта, а также схема применения полученного рекуррентного соотношения" Сухопаров C.Ю. http://vixra.org/pdf/1411.0072v1.pdf

Используемые материалы Править

1. "Вывод рекуррентного соотношения ортогональных многочленов из процесса ортогонализации Грама-Шмидта, а также схема применения полученного рекуррентного соотношения" Сухопаров С.Ю.[1]

2. http://solidbase.karelia.ru/edu/meth_calc/files/09.shtm «Аппроксимация функций методом наименьших квадратов»

Викия-сеть

Случайная вики