Документ взят из кэша поисковой машины. Адрес оригинального документа : http://lnfm1.sai.msu.ru/grav/russian/lecture/mon/6a.doc Дата изменения: Sat Dec 20 13:03:58 2003 Дата индексирования: Mon Oct 1 23:35:46 2012 Кодировка: Windows-1251 Поисковые слова: р с р р с с с с р р с р с с р р

Глава 6

МЕТОД НАИМЕНЬШИХ КВАДРАТОВ

6.1. Средняя квадрати?еская аппроксимация функций

Нау?но-исследовательская практика ?асто нуждается в замене табли?ных
функций, полу?енных экспериментальныо или из астрономи?еских наблюдений,
аналити?еским выражением. Именно так решаются зада?и интерполирования,
?исленного дифференцирования и интегрирования. Однако, нали?ие погрешностей
в этих данных не позволяет подобрать аналити?еские выражения, дающие график
строго по экспериментальным то?кам (например, при "то?е?ном
интерполировании" или при построении интерполяционных сплайнов). Уменьшить
влияние этих погрешностей можно, если пользоваться для аппроксимации
(приближении) наблюденных зависимостей одной вели?ины от другой более
гладкой кривой, ?ем это дает опыт.
Допустим, ?то реально существует зависимость переменной Y от t, но
эта зависимость имеет m степеней свободы, т.е. содержит m неизвестных
параметров, которые нужно подобрать. Предположим, ?то наблюдения Y дают ряд
приближенных зна?ений [pic], полу?енных при n разли?ных зна?ениях t.
Неизвестные параметры обозна?им ?ерез [pic].
Итак, зада?а состоит в том, ?тобы определить оценки параметров [pic]
по заданному наблюдениями ряду [pic]. Запишем постановку зада?и в
математи?еских формулах
[pic].
При m ( n эта зада?а невыполнима, так как имеет бес?исленное множество
решений. При m = n система может иметь однозна?ное решение, а при m ( n -
?исло уравнений больше ?исла неизвестных и уравнения, вообще говоря,
несовместны. Именно к этому слу?аю применяется метод наименьших квадратов
(МНК), которому посвящается данная глава.
Допустим, ?то каким-либо способом мы выбрали набор неизвестных
параметров: [pic]. Подставим их в известную функцию зависимости y от t,
полу?им оценки [pic]:
[pic].
Из-за ошибок наблюдений [pic] - вы?исленное зна?ение - не будет совпадать с
наблюденным [pic]. Набор исходных данных у астрономов принято обозна?ать
буквой О (observatio, лат.), а вы?исленных по формуле - буквой С
(calculatio, лат.). Наилу?шей аппроксимацией будет та, которая
удовлетворяет условию
[pic],
где двойными скобками обозна?ена норма - в определенном смысле среднее
расстояние между О и С.
Как мы уже говорили, ?исло параметров функции есть ее ?исло степеней
свободы. С помощью этих параметров кривая - геометри?еский образ функции -
может менять свою форму и положение на плоскости. Для того, ?тобы провести
прямую, аппроксимирующую линейный тренд (ход) наблюдений, необходимо ввести
два параметра, определяющие параллельный перенос и вращение прямой. Чтобы
определить параболу второго порядка нужно иметь три параметра и т.д.
Синусоида имеет также три степени свободы: положение осевой линии,
амплитуду и фазу. Для определения орбит в зада?ах небесной механики
необходимо определить шесть элементов - это шесть параметров. Для их
определения необходимо минимум три пары координат ( и (, но для уменьшения
влияния погрешностей наблюдений, ?исло пар координат должно быть
существенно больше.
В 1792 году Лаплас в ка?естве нормы отклонений О-С принял сумму
абсолютных зна?ений остато?ных разностей
[pic].
Эта функция не дифференцируема и поиск минимальных зна?ений [pic]
затруднителен.
В 1905 году в книге "Новые методы определения орбит комет" Лежандр
предложил в ка?естве меры то?ности приближения теорети?еской кривой к
наблюдениям принять квадрати?ный критерий
[pic].
В этом слу?ае легко определить абсолютный минимум [pic] как функции
параметров [pic]. Условие минимума приводит к m уравнениям
[pic].
Эти уравнения полу?или название нормальных, а сам метод был назван методом
наименьших квадратов (МНК).
Принцип минимизации суммы квадратов отклонений (О-С) называют иногда
принципом Лежандра. Правда, немецкий у?еный Гаусс этот принцип использовал
раньше Лежандра на 10 лет. Поэтому было бы более правильным называть его
принципом Гаусса-Лежандра. Позже он полу?ил строгое математи?еское
обоснование в теории вероятностей и математи?еской статистики.

6.2. Метод наименьших квадратов (МНК) с независимыми наблюдениями

6.2.1. Применение МНК к линейным функциям

Рассмотрим ?астный слу?ай зависимости наблюдаемой вели?ины [pic] от
искомых параметров [pic]:
[pic].
В дискретных зна?ениях аргумента [pic] полу?им n зна?ений [pic], где [pic]
- ошибка измерения вели?ины l в момент [pic]. Обозна?им [pic], [pic],
[pic]. Теперь для определения искомых x, y, z будем иметь систему n
уравнений:
[pic].
Критерий [pic] будет выглядеть следующим образом:
[pic].
Построим систему нормальных уравнений
[pic],
[pic],
[pic].
Выполнив дифференцирование, полу?им
[pic],
[pic],
[pic].
Имеем систему, состоящую из трех линейных уравнений с тремя неизвестными,
которую легко тем или иным способом решить.
Для обозна?ения сумм произведений или квадратов Гаусс предложил
применять прямые скобки следующим образом
[pic], [pic], [pic] и т.д.
В этих обозна?ениях нормальные уравнения примут вид
[pic],
[pic],
[pic].
Основное свойство нормальных уравнений - симметри?ность матрицы системы:
[pic], [pic], [pic].
Решить эту систему можно используя, например, формулы Крамера, или один из
матри?ных способов.

6.2.2. Ковариационная матрица ошибок неизвестных

Поскольку суммы [pic], [pic], [pic] содержат ошибки наблюдений [pic],
то вы?исленные из нормальных уравнений неизвестные [pic] также
будут содержать ошибки. Следовательно, нашим решением будет лишь оценки
неизвестных [pic]. Запишем нормальные уравнения в матри?ной форме
[pic].
О?евидно, ?то погрешности неизвестных под?иняются той же системе уравнений
[pic],
где элементы матрицы-столбца в правой ?асти уравнений нужно понимать так
[pic], [pic], [pic].
Обозна?им матрицу нормальных уравнений для краткости буквой Ф. Полу?им
[pic],
[pic], [pic].
Пусть наблюдения являются равното?ными и независимыми. Тогда ковариационная
матрица ошибок измерений ряда [pic] является диагональной
[pic] .
Через [pic]мы обозна?или едини?ную матрицу (nxn). Теперь образуем
ковариационную матрицу ошибок оценок неизвестных:
[pic],
где [pic] - дисперсии ошибок неизвестных, [pic] - ковариации этих
неизвестных. Сразу же отметим, ?то если исходные данные представляют собой
независимые наблюдения, то ошибки неизвестных не будут независимыми и их
ковариации, вообще говоря, не равны нулю.
Используя решение нормальных уравнений для ошибок [pic], полу?им
[pic].
Вследствие симметри?ности матрицы Ф ее обратная матрица также симметри?на,
и их транспонирование не изменяет.
Поскольку матрица [pic] не содержит слу?айных погрешностей,
ковариационную матрицу ошибок неизвестных можно переписать так
[pic] .
Распишем первый элемент средней матрицы подробнее
[pic].
Но вследствие равното?ности и независимости измерений [pic] и [pic]
[pic]
Следовательно, в двойной сумме следует оставлять только те ?лены, у которых
i=j. Будем иметь
[pic].
То?но также можно показать
[pic],
[pic], ...
Итак, ковариационная матрица ошибок неизвестных имеет вид
[pic].
Мы полу?или, таким образом, важную формулу, позволяющую вы?ислить не только
дисперсии ошибок неизвестных, но и их ковариации

[pic] .
Сравним полу?енную формулу с формулой для дисперсии ошибки едини?ного
измерения [pic], имеющего вес [pic] (см.гл.5)
[pic].
Видим полную аналогию. Роль дисперсии ошибки измерений играет
ковариационная матрица ошибок неизвестных, а роль веса - матрица нормальных
уравнений. Иными словами: матрица нормальных уравнений есть матри?ный вес
вектора неизвестных.

6.2.3. Вы?исление ошибок неизвестных

Обозна?им элементы обратной матрицы [pic] ?ерез [pic], т.е.
[pic].
Тогда
[pic].
Здесь диагональные элементы - дисперсии ошибок неизвестных:
[pic], [pic], [pic];
а недиагональные - их ковариации:
[pic], [pic], [pic].
Постоянные q11, q22, q33 называют также весовыми коэффициентами (в отли?ие
от весов неизвестных px, py, pz). Чтобы полу?ить веса неизвестных, нужно
найти обратную вели?ину
[pic], [pic], [pic].
В подавляющем большинстве слу?аев нам нужно определить средние
квадрати?еские ошибки МНК-оценок неизвестных. В этом слу?аае вы?исляют СКО
"единицы веса" по формуле, которую мы введем несколько позже:
[pic],
где vk - остато?ные разности: [pic], а m - ?исло неизвестных (у нас
m=3).
Вели?ина (0, как правило, заменяет стандартное отклонение (0 и
является его несмещенной оценкой. Поэтому СКО неизвестных определяют по
формулам
[pic], [pic], [pic].
Для системы трех уравнений веса неизвестных легко вы?ислить, образуя
матрицу Ф. Диагональный элемент q11 обратной матрицы [pic] равен отношению
определителя матрицы, полу?енной вы?еркиванием первой строки и первого
столбца, к определителю системы нормальных уравнений D
[pic],
Аналоги?но
[pic],
[pic].
Веса px, py, pz - соответствено - отношения определителя системы к
определителю матрицы, полу?енной вы?еркиванием соответствующих строки и
столбца, где расположен диагональный элемент обратной матрицы нормальных
уравнений.
Численный пример.
Решим систему уравнений [pic]. Численные зна?ения [pic] сведены в
таблицу

| [pic]| [pic]| | | | |
| | |[pic] |[pic] |[pic] |[pic] |
| 1 | 0 | 2 | 7 | 6.8 | +0.2 |
| 0 | 3 | -2 | 1 | 0.8 | +0.2 |
| -1 | 2 | 0 | 3 | 3.9 | -0.9 |
| 2 | -1 | 1 | 2 | + 2.2 | -0.2 |
| 3 | 2 | -2 | 1 | 0.6 | +0.4 |
| -2 | -1 | 3 | 6 | 5.4 | +0.6 |
| 0 | 3 | -2 | 1 | 0.8 | +0.2 |
| 4 | 1 | 0 | 5 | 5.2 | -0.2 |


Составим нормальные уравнения. Первый элемент матрицы нормальных
уравнений есть сумма квадратов элементов первого столбца ([pic]), второй
элемент первой строки есть сумма произведений элементов второго и первого
столбцов [pic], и т.д. Полу?им
[pic].
Решим полу?енную систему методом обращения матриц:
[pic].
В пятом столбце таблицы приведены вы?исленные зна?ения [pic], а в шестом -
разности [pic]. Минимальная сумма квадратов разностей равна
[pic] , т.е. [pic].
Определим веса неизвестных
[pic], [pic], [pic].
Отсюда
[pic], [pic], [pic].

Решение следует записать в виде
x=0.72(0.10, y=2.29(0.15,
z=3.06(0.16.




6.3. МНК для линейных уравнений. Матри?ный подход

6.3.1. Матри?ная МНК-оценка параметров

Рассмотрим теперь общий слу?ай применения МНК к линейным уравнениям.
Пусть наши уравнения содержат m неизвестых параметров, которые будем
обозна?ать [pic]. Система линейных уравнений, коэффициенты которой будем
с?итать то?ными вели?инами, а наблюдения [pic] содержат погрешности [pic],
имеет вид:
[pic],
[pic],
. . . . . . . . . . . . . . . . . . . . . . . . . .
[pic].
Эти уравнения можно записать компактно в матри?ной форме:
[pic],
где
[pic], [pic], [pic], [pic].
Будем с?итать, ?то вектор l определяется независимыми измерениями, но эти
измерения, вообще говоря, не являются равното?ными. Допустим, ?то [pic]
имеет вес [pic], [pic] - вес [pic], ..., [pic] - вес [pic]. Образуем
квадрати?ескую форму следующим образом:
[pic],
где
[pic], [pic] - оценка вектора параметров х.
Здесь снова l - вектор наблюдений (О), а [pic] - вы?исленный вектор l (С).
Принцип Гаусса-Лежандра в данном слу?ае рассматривается расширенно: будем
минимизировать взвешенную сумму квадратов остато?ных разностей [pic]:
[pic],
[pic].
Для того, ?тобы определить МНК-оценку вектора х нужно найти минимум
квадрати?ной формы [pic]
[pic],
где дифференицироваие по вектору [pic] озна?ает дифференцирование
скалярной функции по каждой из составляющих этого вектора, а результат
располагается в виде матрицы-столбца:
[pic].
Выведем простое правило дифференцирования скаляра по вектору.
Допустим, ?то z - скалярная функция вектора х, а y - известный вектор той
же размерности, ?то и х. Рассмотрим их скалярное произведение
[pic].
Дифференцируя z по [pic], полу?им
[pic], [pic], . . . , [pic].
Располагая результаты в виде столбца, имеем
[pic].
Таким образом, дифференцирование произведения [pic] по х даст y, а
дифференцирование [pic] даст также y. В данном слу?ае это -
транспонированный коэффициент, стоящий перед х.
Итак, при дифференцировании по вектору х результатом будет векторный
коэффициент, стоящий после транспонированной переменной х или
транспонированный коэффициент, стоящий перед векторной переменной х.
Пользуясь этим правилом, продифференцируем квадрати?ную форму [pic]
по
[pic]:
[pic].
Следовательно, нормальные уравнения имеют вид
[pic] .
Это уравнение отли?ается от полу?енного нами в 6.2.1. только тем, ?то оно
у?итывает веса наблюдений. Если наблюдения равното?ны, то [pic] - едини?ная
матрица размера (n x n) и нормальные уравнения совпадают с "класси?ескими"
[pic].
МНК-оценкой вектора параметров х будет
[pic] .
Покажем, ?то это оценка несмещенная. Пусть X - то?ные зна?ения
вектора параметров, а L - наблюдения, не содержащие ошибок. Тогда,
о?евидно, должно выполняться строгое равенство
[pic].
Умножим слева обе ?асти равенства на [pic], полу?им
[pic].
Вы?ислим то?ное зна?ение вектора Х:
[pic].
У?итывая формулу для МНК-оценки х, полу?им
[pic].
Если l содержит только слу?айные ошибки, т.е.
[pic],
то выполняется также равенство
[pic].
Другими словами, [pic].
Вывод: МНК-оценка х является несмещенной при любой матрице весов
наблюдений Р.

6.3.2. Ковариационная матрица МНК-оценки

Вы?ислим теперь ковариационную матрицу ошибок неизвестных параметров.
Для упрощения записи снова обозна?им матрицу нормальных уравнений буквой Ф.
Теперь она имеет вид
[pic].
Вектор ошибок [pic] можно записать так:
[pic],
[pic].
По определению, ковариационная матрица ошибок равна
[pic].
Здесь мы у?ли, ?то Р - диагональная матрица и P=PT.
Ковариационная матрица неравното?ных независимых наблюдений имеет
диагональный вид
[pic].
Перемножим эту матрицу на Р, например, слева. Полу?им
[pic].

Но [pic] есть "дисперсия единицы веса". Следовательно,
[pic].

Вернемся к формуле, определяющей ковариационную матрицу
[pic].
Мы вновь полу?или тот же результат, ?то и при равното?ных наблюдениях (см.
6.2.2): ковариационная матрица неизвестных с то?ностью до множителя равна
обратной матрице нормальных уравнений. В данном слу?ае [pic] - есть
дисперсия ошибки наблюдений, вес которому мы приписали равным единице
(дисперсия единицы веса).
Ковариационная матрица ошибок неизвестных содержит ис?ерпывающую
информацию о то?ности определения каждого параметра из [pic] и ковариации
их ошибок. Остается неопределенной вели?ина [pic]. Ее, как и в слу?ае
одного ряда неравното?ных наблюдений, можно априорно задать или вы?ислить
по результатм эксперимента.

6.3.3 Апостериорная оценка дисперсии единицы веса

Теперь полу?им формулу, позволяющую вы?ислить несмещенную оценку
дисперсии единицы веса по данным наблюдений. Поступим также, как и в слу?ае
определения дисперсии единицы веса при неравното?ных наблюдениях.
Рассмотрим взвешенную сумму квадратов остато?ных разностей
[pic].
Поскольку [pic] является несмещенной оценкой, то справедливо равенство
[pic].
Вы?ислим среднее зна?ение для [pic]
[pic].
Выполним преобразования в каждом слагаемом:

[pic]
[pic]

Но [pic], поэтому
[pic].
Ранее мы полу?или, ?то [pic], поэтому
[pic].
Обозна?им временно
[pic], [pic], [pic].
Теперь [pic].
Поскольку а - вектор-столбец (n x 1), то приведенное выше равенство -
скалярная вели?ина. Заметим, ?то
[pic],
где [pic] - след (сумма диагональных элементов) квадратной матрицы размера
(n х n). Последнее утверждение легко проверить на простом примере:
[pic].
С другой стороны
[pic].
След полу?енной матрицы, равный сумме ее диагональных элементов, есть
[pic]. Таким образом, [pic].
В нашем слу?ае [pic], отсюда следует справедливое утверждение, ?то
[pic].
Но [pic], ибо [pic],.

Умножая полу?енную матрицу справа на Р, полу?им
[pic].
Теперь имеем
[pic].
Третий ?лен равенства для [pic] представляет собой транспонированное
выражение второго ?лена. А так как [pic] - скалярная вели?ина, то оба эти
?лена совпадают.
Наконец, определим последнее, ?етвертое слагаемое
[pic]
Итак, окон?ательно
[pic].
Следовательно,
[pic].
Отбрасывая треугольные скобки, полу?им несмещенную оценку дисперсии единицы
веса - квадрат средней квадрати?еской ошибки единицы веса
[pic],
где [pic], - остато?ные разности.

6.4. Практи?еские методы вы?исления МНК-оценок и их средних квадрати?еских
ошибок

6.4.1. Метод определителей

Этот метод применяют на практике, если ?исло неизвестных невелико и
не превосходит трех. Метод ориентирован прежде всего на "ру?ные
вы?исления" с использованием микрокалькуляторов. В основе метода лежат
известные формулы Крамера.
Пусть система трех уравнений имеет вид
[pic],
[pic],
[pic].

Вы?ислим определитель системы
[pic].
Заменяя последовательно первый, второй и третий столбцы матрицы нормальных
уравнений на столбец правой ?асти, вы?ислим определители
[pic], [pic], [pic].
Теперь неизвестные определяются следующим образом:
[pic], [pic], [pic].
Вы?ислением неизвестных не закан?ивается метод наименьших квадратов.
Необходимо вы?ислить средние квадрати?еские ошибки этих неизвестных.
Весовые коэффициенты суть диагональные элементы [pic] обратной
матрицы системы нормальных уравнений, т.е.
[pic].
Для определения [pic] нужно заменить первый столбец матрицы нормальных
уравнений на первый столбец правой ?асти, вы?ислить определитель полу?енной
матрицы и поделить на D:
[pic].
Для определения [pic] нужно заменить второй столбец матрицы нормальных
уравнений на второй столбец правой ?асти и выполнить аналоги?ные
вы?исления:
[pic].
Поступая то?но также, полу?им
[pic].
Теперь нужно определить квадрат средней квадрати?еской ошибки единицы веса
[pic],
где x, y, z - зна?ения наших неизвестных, полу?енные из решения системы
нормальных уравнений.
Последнюю формулу можно преобразовать так, ?то необходимость
вы?исления остато?ных разностей типа [pic] отпадет. Запишем [pic]следующим
образом
[pic]
Используя скобки Гаусса, перепишем полу?енный результат
[pic]
Легко видеть, ?то последние три выражения в квадратных скобках равны нулю,
так как x, y, z удовлетворяют системе нормальных уравнений. Отсюда следует
[pic]
Поскольку [pic] мы уже вы?ислили, составляя нормальные уравнения, осталось
довы?ислить сумму квадратов правых ?астей исходных уравнений (наблюдений)
[pic]. Результатом обработки данных будут зна?ения x, y, z и их ошибки:
[pic] , [pic] , [pic] .


6.4.2. Метод последовательных исклю?ений (схема Гаусса)

Этот метод можно применять к любому коли?еству неизвестных. Схема,
предложенная Гауссом первона?ально, была расс?итана для ру?ного с?ета с
использованием специальных вы?ислительных бланков. Схема предусматривает:
1. составление нормальных уравнений с параллельным контролем правильности
вы?ислений,
2. последовательное исклю?ение неизвестных, на?иная с первого. Вы?исления
сопровождаются контролем,
3. параллельное вы?исление весовых коэффициентов,
4. определение неизвестных, на?иная с последнего,
5. вы?исление суммы квадратов остато?ных разностей, вы?исление ошибок
неизвестных.
Несколько слов о терминологии.
Исходные уравнения, надлежащие математи?еской обработке косвенных
наблюдений методом наименьших квадратов, в астрономи?еской традиции
называют условными уравнениями. Этот термин заимствован из зада? геодезии,
где МНК применяется довольно широко. Однако, использование его в
астрономи?еской практике не имеет основания. Дело в том, ?то в геодезии
дело имеют с линейными уравнениями, в которых на неизвестные накладываются
дополнительные условия. Например, если x, y, z - есть три угла
треугольника, то какими бы они ни были, должно выполняться равенство
x+y+z=1800. В противном слу?ае отрезки прямых не образуют замкнутую
геометри?ескую фигуру: образуется невязка - несовпадение первой и последней
то?ки замкнутого треугольника (или полигона). Термин невязка также
применяется в астрономи?еской практике, но несколько в ином смысле: это
разность (О-С):
[pic].
Во избежание недоразумений, вели?ины [pic] следует называть остато?ными
разностями, т.е. так же, как их называют в англоязы?ной литературе
(residual).
Для иллюстрации схемы Гаусса ограни?имся всего тремя неизвестными,
хотя ее можно распространить на любое их коли?ество.
Составляется таблица исходных уравнений:
| x | y | z | | | | |
| a1 | b1 | c1| l1 | s1| | v1 |
| | | | | |[pic] | |
| a2 | b2 | c2| l2 | s2| | v2 |
| | | | | |[pic] | |
| a3 | b3 | c3| l3 | s3| | v3 |
| | | | | |[pic] | |
| a4 | b4 | c4| l4 | s4| | v4 |
| | | | | |[pic] | |
| ...| ...| | ...| | | ...|
| | |... | |... |... | |
| an | bn | cn| ln | sn| | vn |
| | | | | |[pic] | |


В таблицу исходных уравнений помещают :
1) коэффициенты перед неизвестными [pic];
2) правые ?асти lk несовместных "уравнений" [pic];
3) контрольные суммы [pic];
4) зна?ения [pic]. Этот столбец вы?исляется после того, как неизвестные
будут определены;
5) остато?ные разности [pic] в последнем столбце таблицы.
Вы?исление контрольных сумм (пятый столбец) необходимо в слу?ае
ру?ного с?ета (на микрокалькуляторе). Эти суммы предохранят от возможных
промахов при составлении нормальных уравнений и их решении. Однако, ?исла
[pic] и lk должны быть выписаны с одинаковым ?ислом знаков после запятой. В
слу?ае, когда они определены ина?е, прибегают к приему "выравнивания"
порядков коэффициентов соответствующим выбором неизвестных. Поясним
сказанное на примере. Пусть
[pic].
Как видно из этого уравнения, все ?исла ?етырехзна?ные. Следовательно,
можно пользоваться ?етырехзна?ными вы?ислениями. Следует ожидать, ?то и
решения будут иметь ?етыре зна?ащих цифры . Однако, вы?исление суммы
[pic]
не имеет смысла, т.к. первое ?исло имеет только одну цифру после запятой, а
третье - шесть. Дополнение "нужных мест" нулями (?то возможно) приведет к
резкому увели?ению ?исла знаков для s
[pic].
Мы полу?им девятизна?ную контрольную сумму. Чтобы этого избежать,
перепишем уравнение следующим образом: