[Все статьи]

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

Наш мир не идеален, ни в чем нельзя быть уверенным с абсолютной точностью. Кто помнит лабораторные работы по физике, тот должен знать, что измерение какой-либо физической величины обычно проводят несколько раз при различных условиях, а найденный результат записывают в виде 20,3±2,3. Это необходимо для того, чтобы нейтрализовать погрешности приборов, трясущиеся руки экспериментатора, вспышки на солнце и так далее.

Метод наименьших квадратов (далее МНК), о котором пойдет речь в этой статье, является одним из способов противостоять ошибкам измерений.

Общая формулировка метода выглядит так:

Пусть имеется система уравнений:

Формула 1

Здесь f(x) эта некая функция, конкретный вид которой не известен, известен лишь ее общий вид, например известно что это прямая, или многочлен, или синусоида и так далее. МНК позволяет зная общий вид функции найти ее конкретный вид (коэффициенты) который наилучшим образом вписывается в экспериментальные данные.

Особенностью МНК является то, что число уравнений превышает число неизвестных коэффициентов в функции f(x). Таким образом, в общем случае точного решения системы не существует.

Обратите внимание, что система решается не относительно xn, а относительно неизвестных коэффициентов функции f(x).

Основная идея МНК состоит в том, чтобы при нахождении конкретного вида функции минимизировать сумму квадратов ошибок во всех исходных уравнениях.

Иными словами нужно свести к минимуму функцию:

Формула 2

Может возникнуть вопрос почему сумма квадратов? Дело в том, что во-первых, квадрат любого числа всегда неотрицателен, и следовательно сумма квадратов всегда не отрицательна, т.е. ограничена снизу, а следовательно у нее есть минимум; во-вторых, при нахождении минимальной суммы квадратов мы уменьшаем максимальную ошибку.

Согласитесь, что иметь в двух точках ошибку в 5 единиц, лучше, чем в первой точке иметь нулевое отклонение, зато во второй точке иметь отклонение 10. Сумма отклонений в обоих случаях будет одинаковой, а вот сумма квадратов отклонений в первом случае будет меньше.

Чтобы визуально понять что такое МНК, можно поиграть с функцией «Линия тренда» в Microsoft Excel.

Для примера решим небольшую задачку. Пусть имеются результаты некого физического эксперимента. Представим, что мы измеряли значение параметра Y при различных значениях параметра X.

X Y
1 1.1
3 1.98
5 3.2

Пусть, мы так же знаем, что теоретически эти параметры связаны линейной зависимостью:

Формула 3

Параметров k и b мы не знаем, их предстоит найти.

Если построить график этой зависимости, то будет видно, что точки лежат примерно на одной прямой.

График 1

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

Придется искать прямую, проходящую примерно через данные точки. В чем нам и поможет МНК.

Составим дополнительную таблицу.

X Y Ошибка Квадрат ошибки
1 1.1 k * 1 + b - 1.1 (k * 1 + b - 1.1)2
3 1.98 k * 3 + b - 1.98 (k * 3 + b - 1.98)2
5 3.2 k * 5 + b - 3.2 (k * 5 + b - 3.2)2
Сумма (k * 1 + b - 1.1)2 + (k * 3 + b - 1.98)2 + (k * 5 + b - 3.2)2

Таким образом, задача состоит в том, чтобы найти такие k и b, чтобы сумма представленная в последней строке (назовем ее функцией ошибки) таблицы была минимальной.

Как известно из курса математического анализа, экстремум (минимум или максимум) функция многих переменных достигает в точке, где ее частные производные равны нулю.

Таким образом:

Формула 4

Приравнивая производные к нулю и получим систему уравнений:

Формула 6

Решение которой находится тривиально. Таким образом получим, что k=0.525, b=0.518. Вообще говоря, в общем случае, нужно еще доказать, что это именно минимум, и что это минимум глобальны (в случае если система имеет несколько решений). Однако, если мы в заданный набор данных пытаемся вписать полином, то найденное решении и будет глобальным минимумом.

Тогда уравнение прямой запишется:

Формула 7

На графике это будет выглядеть следующим образом:

График 2

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

Бывает, что общий вид уравнения кривой неизвестен, а нужно просто подобрать какую-то функцию, чтобы промоделировать тот или иной процесс. Тут начинается полностью творческий процесс по подбору. Иногда общий вид угадывается сразу, иногда приходится поломать голову и повспоминать, что проходили на математическом анализе.

Например, в следующий набор экспериментальных данных хорошо впишется кривая вида

Формула 8
График 3

А для следующих данных:

График 4

может подойти:

Формула 9

или

Формула 10

Первый вариант дает меньшие отклонения, но в этом случае, чтобы найти минимум функции ошибки придется решить систему нелинейных уравнений и выделить глобальный минимум. Наиболее целесообразно делать это численными, а не аналитическими методами.

Второй же случай приведет к системе линейных уравнений, которая решается тривиально (не смотря на большое число неизвестных).

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

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

Рассмотрим квадратную систему линейных уравнений.

Формула 11

Ее можно переписать в матричном виде:

Формула 12

Введем дополнительные обозначения:

  • A – матрица коэффициентов (или главная матрица системы);
  • x – вектор-столбец неизвестных;
  • B – вектор-столбец свободных членов.

Тогда систему можно переписать следующим образом:

Формула 13

Умножим обе части уравнения слева на матрицу обратную матрице А.

Формула 14

Умножение матрицы на обратную ей дает в результате единичную матрицу (или матрицу идентичности, как ее еще называют), таким образом, имеем, что:

Формула 15

Таким образом, все что нужно сделать чтобы решить квадратную систему уравнений, это найти матрицу обратную матрице коэффициентов и умножить полученный результат на вектор-столбец свободных членов. Получившийся результат даст нам вектор-столбец значений неизвестных.

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

Таким образом решение системы уравнений требует всего двух трех строчек программного кода (не считая затрат на вычисление коэффициентов в уравнениях).

Описанный выше подход применим когда мы имеем дело с квадратной системой уравнений.

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

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

Однако, здесь на помощь приходит понятие обратной матрицы Мура-Пенроуза, или, как ее еще называют, псевдообратной матрицы.

Псевдообратная матрица в общем случае вычисляется по формуле:

Формула 16

Здесь А* - это Эрмитова матрица, или самосопряженная матрица. Так как в нашем случае мы имеем дело с полем действительных чисел, то ее можно заменить на траснпонированую.

Таким образом, имеем

Формула 17

Соответственно приближенное решение системы линейных уравнений будет иметь вид:

Формула 18

В теории линейной алгебры доказывается, что решение полученное подобным образом, совпадет с решением, полученным согласно общего алгоритма МНК.

Стоит отметить, что операции транспонирования, обращения и перемножения матриц как правило входят в библиотеки функций/классов для операций с матрицами. Таким образом решение подобной задачи программным способом не представляет трудностей.

В последующих статьях будет рассмотрено практическое применение МНК.

Цветков Павел.
Апрель, 2008 г.

[На верх страницы]
Loading... Загрузка...