Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
и рассмотреть задачу наименьших квадратов
(25.53) (25.54)
(25.55)
где р связано с х ортогональным линейным преобразованием
х = Vp. (25.56)
Задача (25.55) эквивалентна задаче Ах=*Ъ ъ том смысле, как эта эквивалентность была определена в гл. 2 для произвольных ортогональных преобразований задачи наименьших квадратов.
Поскольку S — диагональная матрица (5 = Diag (ss„}), то влияние каждой компоненты р на норму невязки видно непосредственно. Вводя в решение компоненту ру со значением
Р1 = Ф1, (25.57)
мы снижаем квадрат нормы невязки на величину ?2 .
Предположим, что сингулярные числа упорядочены так, что $* > sk+il * = 1,...,л—1. Тогда естественно рассмотреть "пробные" решения зада-, чи (25.55) следующего вида:
л
Pi
,(*) =
Рк
О
О J
* = 0,1.....н
(25.58)
мша pj шшсделяюня ширмуламИ (23.3 }). прооньпгвекторр^"7 является нормальным псевдорешеннем задачи (25.55), если сингулярные числа Sj,) > к, считаются нулевыми.
Из пробных векторов р(Ас* получаем пробные решения дг*** задачи Ах=* Ь:
*(*) = урС)= ? к = 0,...,п. (25.59)
Здесь -/-й столбец V. Заметим, что (
|| г = ||p(ft) || 2 = | fta = ? /ii-V (25.60)
/=1 /=1 \ Sj J
и, следовательно, Идг(** II - неубывающая функция от *. Квадрат нормы невязки, отвечающей вектору х , равен
m
р2 = II6 - Лх<*> В2 = ? gj. (25.61)
Исследование столбцов матрицы V, ассоциированных с малыми сингулярными числами, — очень эффективное средство выявления почти линейно зависимых наборов столбцов А (см. (12.23), (12.24)).
Наша практика показывает, что вычисление и вывод на печать матрицы V и величин sk,sll ,рк,х{к), Вдс(Л) И,Р* и Рк для всех значений к (от 0 или 1 до и) чрезвычайно полезны при анализе трудных практических задач метода наименьших квадратов. Для некоторых статистических интерпретаций (см. (12.2)) интерес представляют также величины
(25.62)
Предположим, что матрица А плохо обусловлена. Тогда некоторые из сингулярных чисел с большими номерами будут существенно меньше предшествующих. В этом случае некоторые значения р, с большими номерами могут быть слишком велики. В типичной ситуации пытаются найти такой индекс к, чтобы все коэффициенты р,-, / <к, были достаточно малы, все сингулярные числа х,-, / <к, достаточно велики, а норма невязки р* достаточно мала. Если такой индекс к существует, то в качестве приемлемого решения можно взять пробный вектор ****.
Эта техника успешно применялась во многих приложениях (см., например, [87], где рассматривается приложение к численному решению фредгольмовых интегральных уравнений первого рода).
Как только вычислены сингулярные числа и вектор g, уже просто вычислить для ряда значений \ числа II II и cj\, нужные в методе стабилизации Левенберга-Марквардта (см. (25.47) и (25.48)). Эти величины представляют интерес из-за свойства оптимальности, выражаемого теоремой 25.49.
В следующей главе рассматривается пример, в котором все упомянутые здесь величины вычисляются и интерпретируются для конкретного набора входных данных.
152
ГЛАВА 26
ПРИМЕРЫ НЕКОТОРЫХ МЕТОДОВ АНАЛИЗА ЗАДАЧИ НАИМЕНЬШИХ КВАДРАТОВ
Рассмотрим задачу наименьших квадратов Лхэг ft, матрица коэффициентов которой приведена в табл. 26.1. Мы предположим, что имеется неопределенность порядка 0,5 X 10"8 в элементах А и порядка 0,5 X 10~4 в элементах ft.
Рассмотрим вначале сингулярный анализ этой задачи, как он описан в § 6 гл. 25. Результаты, полученные при выполнении программы сингулярного анализа на машине UNIVAC 1108, воспроизведены на рис. 26.1.
В общем случае программа производит замену переменных х = Dy и решает задачу (AD)y э< ft. В данном примере мы положили D = 7, поэтому символ у на рис. 26.1 нужно отождествить с х.
Напомним, что вектор g (столбец выдачи, стоящий под шапкой
"G") вычисляется по формуле g = UTb, где U - ортогональная матрица. Так как четвертая и пятая компоненты g меньше предполагаемой неопределенности в ft, то мы хотели бы рассматривать эти две компоненты как нулевые. Это свидетельствует в пользу того, чтобы считать третий пробный вектор х(3) (столбец выдачи под шапкой" SOLN 3 " ) наиболее удовлетворительным решением.
Другой метод состоит в сравнении чисел ак (см. (25.62)), стоящих в столбце под шапкой " N.S.R.C.S.S. " (что является сокращением от Normalized Square Root of Cumulative Sum of Squares, т.е. нормированный квадратный корень из накопленной суммы квадратов), с предполагаемой неопределенностью 0,5 X 10"4 в ft. Замечаем, что аг ( = 1,1107 X 10~2) значительно больше этой предполагаемой неопределенности, в то время как а3 ( = 4,0548 X 10"5) чуть меньше ее. Это может служить основанием для
того, чтобы выбрать х*3* в качестве предпочтительного пробного решения.
Если предполагаемая неопределенность в ft интерпретируется статистически как среднеквадратичное отклонение погрешностей в ft, то произведение этой величины (0,5 X 10~4) на обратное к сингулярному числу (эти обратные числа находятся в столбце с шапкой "1/синг.число") дает среднеквадратичное отклонение соответствующей компоненты вектора р (столбец с шапкой "Р"). Как видим, первые Три компоненты р превосходят по абсолютной величине соответствующие среднеквадратичные отклонения, в то время как последние две компоненты меньше своих отклонений. Это также может быть причиной, чтобы предпочесть пробное
