Дж. Деммель - Вычислительная линейная алгебра, страница 26

Однако, когда значение степени достигает 18, норма невязки неожиданно и очень резко возрастает. На левом рисунке можно видеть, как хаотически ведет себя график многочлена 19-й степени (сплошная линия). Мы увидим позднее, что зто связано с плохой обусловленностью задачи. Обычно, чтобы избежать плохо обусловленных задач, для аппроксимации используют многочлены сравнительно невысоких степеней [6Ц.

В Ма11аЬ'е полиномиальная аппроксимация данных осуществляется функцией ро1ууйс. Полиномиальная аппроксимация допускает альтернативу. В более общем случае, имеется система независимых функций /»(у),..., /„(У), действующих нз в««в К, и набор точек (ум 61),...,(У,„,Ь,„), где у; Е 11" и 6 Е К. Нужно найти наилучшее приближение этих точек вида 6 = 2 " х./ (у). Другими словами, нужно выбрать х = [хы..., х„]т так, чтобы минимизировать невязки т; = 2 ", хз/з(у«) — Ь; для 1 < «< т.

Полагая а, = /1(у«), можем записать это в виде т = Ах — 6, где А — матрица размера т х п, х — вектор размерности и, а 6 и т — векторы размерности т. Удачный выбор базисных функций /,(у) может дать лучшее приближение и лучше обусловленную задачу, чем при использовании многочленов [33, 84, 168].

О ПЗ 3.1. Введение Исходныедвниые (помечеиы символом о) и полиномивльиые приближения з Норма иевлзки в задаче наименьших квадратов !О 10 1 О-и -Э -то -о о в 1о 0 О 1О 10 20 Степень аппроксимирующего миогочленв Рис. 3.1. Полиномивльная аппроксимация кривой Ь = в1п(ху/5) + у/5 и нормы не- вязок (см. цветной вариант рисунка на обложке книги. — Перев.).

среднего балла (ат), который абитуриент имел в старших классах школы, и результатов двух тестов готовности: устного (аз) и письменного (аз). Основываясь на прошлых данных для ранее принятых студентов, можно построить линейную модель вида Ь = ~ ., а.х . Наблюдениями являются четверки чисел з ап, аип агз и Ьи по одной для каждого из т студентов в базе данных.

Итак, желательно минимизировать величину 61 Ь Гт Гг а11 а12 а13 а21 а22 а23 [" 1- :— А х — 6, аш1 ашз а,„г что можно рассматривать как задачу наименьших квадратов. Пример 3.3. Задача наименьших квадратов впервые была поставлена и сформулирована Гауссом при решении практической задачи по заказу германского 1 Стандартные обозначения в статистике отличаются от принятых в линейной алгебре: статистики пишут Х)З = у вместо Ах = Ь. Существует статистическое обоснование метода наименьших квадратов, который статистики называют линейной регрессией. Предположим, что числа а, известны точно, так что шум присутствует лишь в числах Ьб при этом шумы в различных Ь, независимы и нормально распределены с нулевым средним и одним и тем же средним квадратичным отклонением о. Пусть х — решение задачи наименьших квадратов, а хт — вектор истинных значений параметров.

Тогда х называется оценкой наибольшего правдоподобия для хт, а ошибка х — хт нормально распределена и имеет нулевое среднее и ковариационпую матрицу аг(АтА) 1. Мы встретимся с матрицей (А А) 1 снова при решении задачи наименьших квадратов методом нормальных уравнений. Более подробные сведения о статистических аспектах задачи' можно найти в (ЗЗ, 259). О 114 Глава 3. Линейные задачи наименьших квадратов правительства. Существуют важные экономические н юридические причины для точного определения границ между участками земли, принадлежащими разным хозяевам.

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

Задача Гаусса состояла в том, чтобы предложить метод обновления координат вех из государственной базы данных, основываясь на этих более точных измерениях. Для этой цели он и изобрел метод наименьших квадратов, который мы вскоре изложим [33]. Задача, решенная Гауссом, не потеряла своей актуальности и к ней периодически приходится возвращаться. В 1974 г.

Национальная Геодезическая Служба США начала работу по обновлению геодезической базы данных США, состоящей из приблизительно 700 000 точек. Возросшие стимулы для такой работы включали в себя предоставление достаточно точных данных инженерам и местным планировщикам для строительных проектов и геофизикам для изучения движения тектонических плит земной коры (они могут перемещаться на расстояние до 5 см за год). Соответствующая задача наименьших квадратов была на тот момент самой большой из когда-либо решавшихся: примерно 2.5 миллиона уравнений относительно 400 000 переменных. Кроме того, задача была очень разреженной, что сделало возможным ее решение в 1978 г. на имевшихся в то время компьютерах [164].

Теперь мы кратко обсудим постановку задачи. Она, в действительности, нелинейная и решается посредством аппроксимации ее последовательностью линейных задач; каждая из них есть линейная задача наименьших квадратов. База данных представляет собой список точек (геофизических вех), снабженных координатами: широтой, долготой и, возможно, высотой. Для простоты изложения будем считать землю плоской и сопоставим каждой точке ь' линейные координаты в, = (хм у,)т.

Для всех точек нужно вычислить поправки бв; = (бхь,бу;)т так, чтобы измененные координаты в[ = (х'„у,')т = в;+ба; лучше соответствовали новым, более точным измерениям. Эти измерения включают в себя как расстояния между избранными парами точек, так и углы между отрезками, соединяющими точку 1 с точками у и и (см. рис.

3.2). Чтобы понять, как превратить эти новые измерения в уравнения, рассмотрим треугольник на рис. 3.2. Кто углы помечены соответствующими (поправленными) координатами; показаны также углы и' и длины сторон Ь. С помощью этих данных и простых тригонометрических тождеств легко выписать необходимые уравнения. Например, улучшенное измерение для д; можно использовать в соотношении [(вз в[) (вй в[)! где соз в; выражен посредством скалярных произведений соответствующих сторон треугольника. Предполагая, что величины бв; малы в сравнении с вп это соотношение можно линеаризовать таким образом: умножим обе части на знаменатель дроби, выполним все необходимые умножения, что приводит к много- члену четвертой степени относительно «б-переменных» (вроде бх;), и отбросим все члены, содержащие в качестве множителя более чем одно 5-переменное.

Это 115 3.2. Матричные разложения г' „-- ( х', у) ) х'), = ( х'), у), ) х')=(х'(,у) ) Рис. 3.2. Вывод уравнений при обновлении геодезической базы данных. дает уравнение, в которое б-переменные входят линейно. Если собрать подобные линейные уравнения для всех новых измерений углов и расстояний, то получится переопределенная линейная система относительно всех б-переменных.

Мы хотели бы найти наименьшие поправки, т. е. наименьшие значения величин бт, и т.д., которые бы наилучшим образом соответствовали этим поправкам. Это и есть задача наименьших квадратов. О Позднее, когда будет развит соответствующий аппарат, мы покажем, что и сжатие изображений может интерпретироваться как задача наименьших квадратов (см. пример 3.4). 3.2. Матричные разложении для решения линейной задачи наименьших квадратов Линейная задача наименьших квадратов может быть решена несколькими явными способами, к обсуждению которых мы теперь переходим.

Эти способы; 1. нормальные уравнения, 2. ЯК-разложение, 3. ЯЧП, 4. преобразование в линейную систему (см, вопрос 3.3). Первый метод наиболее быстр', но и наименее точен; он пригоден для задач с малым числом обусловленности. Второй метод является наиболее стандартным; его стоимость может превышать стоимость первого метода в два раза.

Третий метод наиболее полезен для плохо обусловленных задач, т. е. задач, в которых матрица А не имеет полного ранга; он в несколько раз дороже предыдущих методов. Последний метод позволяет проводить итерационное уточнение приближенного решения в случае плохо обусловленной задачи. Все методы, кроме третьего, могут быть адаптированы для эффективной обработки разреженных матриц (33]. Мы последовательно обсудим все эти методы. Для методов 1 и 2 мы поначалу предположим, что А имеет полный столбцевой ранг п.

116 Глава 3. Линейные задачи наименьших квадратов 3.2.1. Нормальные уравнения Для вывода нормальных уравнений будем искать точку х, в которой градиент функции ЙАх — ЬЯ =- (Ах — Ь)т(Ах — Ь) обращается в нуль. Итак, мы хотим, чтобы (А(х+ е) — Ь)т(А(х+ е) — Ь) — (Ах — Ь)т(Ах — Ь) О= 1пп еео ~~е!1г 2ет(АтАх АтЬ) + етАтАе = 1пп е-ео !)е11г Второе слагаемое — ( ~~ — ~~х — л стремится к нулю при е э О, поэто- ))е()г — ~!е)!» му множитель А Ах — АтЬ в первом слагаемом также должен быть нулем, т.е.