Метод наименьших квадратов: вывод формул и пример с 5 точками

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

Почему не интерполяция#

У вас уже есть инструмент для восстановления функции по точкам — интерполяция: многочлен, проходящий через все узлы. Но интерполяция честна только тогда, когда данные точные — значения функции в узлах вычислены, а не измерены. Измерения содержат шум, и интерполяционный многочлен с рвением воспроизводит каждую случайную кочку: график начинает змеиться между точками, а прогноз за пределами данных улетает в космос. Пять зашумлённых точек интерполировать многочленом четвёртой степени — значит построить модель шума вместо модели сигнала.

Тут и появляется метод наименьших квадратов (МНК). Параметров у модели мало (у прямой — два), данных много (пять, десять, миллион точек). Системы «уравнений больше, чем неизвестных» точного решения не имеют, зато имеют наилучшее компромиссное: минимизируем сумму квадратов отклонений.

Смотрите, что произошло с алгеброй. Прямая имеет два параметра, данных пять — система пяти уравнений $y_i=a{x}_i+b$ с двумя неизвестными переопределена: общего решения у неё нет, третья точка почти наверняка не ляжет на прямую через первые две. Раньше такие системы браковали как несовместные. МНК меняет вопрос: вместо «существует ли точное решение» спрашивается «какое приближение наименее виновато». Это тот же дух, что в численных методах решения СЛАУ: там мы учились жить с ошибками округления, здесь — с ошибками измерений.

Сводим задачу к минимуму суммы квадратов#

Для прямой $y=ax+b$ и точек $(x_i,y_i)$ остаток (отклонение) $i$-й точки равен $y_i-a{x}_i-b$. МНК выбирает коэффициенты, при которых сумма их квадратов минимальна:

Почему квадраты, а не модули? Квадрат гладкий — минимизируется производными, даёт явные формулы, удобен для анализа ошибок и связан с дисперсией из статистики. Минимум суммы модулей тоже существует (это медианная регрессия), но там формулы замысловатее. Нелинейность по $a$ и $b$ здесь только кажущаяся: $S$ — многочлен второй степени от $a$ и $b$, а у такого параболоида минимум ищется в лоб, приравниванием частных производных к нулю.

Нормальные уравнения: вывод для прямой y = ax + b#

В матричной записи обе строки сворачиваются в одни нормальные уравнения: $A^{T}Ac=A^{T}y$, где столбцы $A$ — это $x_i$ и единицы, а $c=(a;b)$. Уравнение так и называется «нормальным», потому что содержит условие ортогональности: остатки $y-Ac$ перпендикулярны столбцам $A$, то есть вектор $Ac$ — это проекция данных на плоскость моделей. Решать нормальные уравнения можно тем же методом Гаусса или LU-разложением.

Картинка, которая всё объясняет: векторы данных живут в пятимерном пространстве (по одному измерению на точку), а всевозможные прямые $ax+b$ заметают в нём двумерную плоскость. МНК опускает из точки данных перпендикуляр на эту плоскость — основание перпендикуляра и есть предсказание модели. Квадрат длины перпендикуляра — это в точности сумма квадратов остатков, а условие перпендикулярности — в точности нормальные уравнения. Геометрия превращает «подгонку» в точную конструкцию.

Считаем руками: пять точек от и до#

Берём конкретные данные: точки $(1;2)$, $(2;3)$, $(3;5)$, $(4;4)$, $(5;6)$. Все суммы считаем честно, в столбик. На экзамене схема та же: сначала пять сумм, потом две формулы. Ошибаются чаще всего в числителе наклона — держите в голове правило «$n$ умножить на перекрёстное произведение минус произведение сумм»:

Подставляем в выведенные формулы. Определитель $D=n\sum x^2-(\sum x{)}^2=5\cdot 55-15^2=275-225=50$. Наклон: $a=\frac{n\sum xy-\sum x\sum y}{D}=\frac{5\cdot 69-15\cdot 20}{50}=\frac{345-300}{50}=\frac{45}{50}=0,9$. Сдвиг: $b=\frac{\sum y-a\sum x}{n}=\frac{20-0,9\cdot 15}{5}=\frac{20-13,5}{5}=1,3$.

Проверка добросовестности — остатки. Значения прямой: $2,2$; $3,1$; $4,0$; $4,9$; $5,8$. Отклонения данных от модели: $-0,2$; $-0,1$; $+1,0$; $-0,9$; $+0,2$ — сумма ровно ноль (МНК обязан балансировать остатки), сумма квадратов $S_{\min }=0,04+0,01+1+0,81+0,04=1,9$. Сдвиньте мысленно прямую вверх или увеличьте наклон — сумма квадратов вырастет: вы стоите в минимуме параболоида. Прогноз для $x=6$: $\hat{y}=6,7$.

1. 1
   Собрать данные в таблицу и посчитать пять чисел: $\sum x$, $\sum y$, $\sum xy$, $\sum x^2$ и $n$
2. 2
   Посчитать определитель $D=n\sum x^2-(\sum x{)}^2$ — он должен быть заметно больше нуля, иначе все $x$ совпадают и наклон не восстановить
3. 3
   Найти $a$ и $b$ по формулам вывода
4. 4
   Проверить баланс: сумма остатков должна выйти нулём (в пределах округлений)
5. 5
   Посмотреть на график: модель обязана согласовываться со смыслом данных, а не только с суммами

Полиномиальный МНК и связь со статистикой#

Прямая — частный случай. Если точки выстраиваются в дугу, берут модель $y\approx c_0+c_{1}x+c_2{x}^2$ и повторяют ту же процедуру: минимизируют сумму квадратов, приравнивают производные по $c_0,c_1,c_2$ к нулю — выходят нормальные уравнения размера 3×3. Механика та же, растёт только размер системы. Соблазн увеличивать степень огромен: с многочленом достаточно высокой степени можно пройти сколь угодно близко к каждой точке. Но это возврат к интерполяции с её змеящимся графиком: модель начинает описывать шум. Степень выбирают по смыслу данных и следят за поведением между узлами; для высоких степеней нормальные уравнения вдобавок становятся плохо обусловленными — та самая история с источниками погрешностей.

У МНК двойное гражданство: в численных методах это аппроксимация, в статистике — регрессия. Если шум нормальный и одинаковый по разбросу, оценка МНК совпадает с оценкой максимального правдоподобия — глубинная причина, почему метод так живуч. Качество подгонки измеряют коэффициентом детерминации: $R^2=1-\frac{S_{\min }}{S_{\text{полн}}}$, где $S_{\text{полн}}$ — разброс данных относительно среднего $\bar{y}$. В нашем примере $\bar{y}=4$, $S_{\text{полн}}=4+1+1+0+4=10$, значит $R^2=1-1,9/10=0,81$: модель объясняет 81% изменчивости данных. Для прямой есть бонус: $R^2$ равен квадрату коэффициента корреляции Пирсона — две характеристики одной природы. Сам термин метод наименьших квадратов разобран в словаре.

Две оговорки про честность метода. Первая — выбросы: квадрат растёт быстрее отклонения, поэтому один дикий замер тянет прямую на себя сильнее ста нормальных. Данные с подозрением на сбои сначала чистят или переходят к робастным вариантам, где крупные остатки штрафуются мягче. Вторая — вычислительная: у матрицы $A^{T}A$ обусловленность равна квадрату обусловленности $A$, и на неудачных данных нормальные уравнения теряют точность. Поэтому промышленные пакеты решают ту же задачу через QR-разложение — ответ тот же, чувствительность к округлениям ниже.

Когда МНК уместен#

• Данных много, параметров мало — прямая, парабола, экспонента, подгоняемая к облаку точек
• Шум случайный и примерно одинаковый по всему диапазону измерений
• Нужен прогноз или сглаженная зависимость, а не точное прохождение через узлы
• Выбросы проверены: один дикий замер тянет прямую на себя, квадраты усиливают его влияние
• Данные точные по определению (таблица значений функции) — тогда честнее интерполяция

Пример неудачного применения — экстраполяция. Прямая, подогнанная к продажам за квартал, не имеет права предсказывать продажи через десять лет: модель описывает локальный тренд, а не закон природы. Хорошее применение — калибровка: пять эталонных температур, пять замеров сопротивления, прямая МНК — и датчик готов к работе в измеренном диапазоне. Разница между двумя сценариями не в математике, а в вопросе, который вы задали данным.

Закрепите руками: возьмите любые пять точек из своего задачника, посчитайте четыре суммы и два коэффициента по формулам выше, а затем сверьте график в симуляторе. Дальше по курсу та же идея минимизации квадратов появится в роли движка: в численном решении ОДУ ошибки шага тоже складывают в квадраты, когда сравнивают методы. Термин система линейных уравнений из нормальных уравнений — в словаре.