Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
вектор уИ1'3,5) , который как своей нормой, так и нормой соответствующей невязки очень походит на указанные ранее векторы хиг*3'.
Дальнейшее поведение пошаговой регрессии критически зависит от деталей конкретного алгоритма. Это связано с очень плохой обусловленностью задач, отвечающих множествам из четырех или пяти столбцов.
ГЛАВА 27
МОДИФИКАЦИЯ РАЗЛОЖЕНИЯ ПРИ ДОБАВЛЕНИИ ИЛИ УДАЛЕНИИ СТРОКИ
(С ПРИЛОЖЕНИЯМИ К ПОСЛЕДОВАТЕЛЬНОЙ ОБРАБОТКЕ ЗАДАЧ С БОЛЬШИМИ ИЛИ ЛЕНТОЧНЫМИ МАТРИЦАМИ КОЭФФИЦИЕНТОВ)
В этой главе мы адаптируем ортогональные преобразования к последовательной обработке данных для задачи НК. Подобная адаптация является средством экономии машинной памяти для определенного типа задач, связанных с обработкой очень большого объема информации.
Описываемые нами методы имеют важные приложения и к задачам, где нужно получить последовательность решений для массива данных, к которому последовательно добавляется (или из него изымается) информация. Потребность в вычислениях такого рода возникает в классе задач, называемых последовательным оцениванием или фильтрацией. Методы окаймления, представленные в данной главе, имеют отличную численную устойчивость. Не так обстоит дело во многих опубликованных алгоритмах для этого класса задач, которые основаны на модификации обратной матрицы.
Удаление данных может быть внутренне неустойчивой операцией при определенной связи удаляемых строк со всем информационным массивом. Мы предложим некоторые методы, в которых сделана попытка избежать внесения в задачу дополнительной численной неустойчивости.
Запишем задачу наименьших квадратов обычным образом:
Ах^Ъ, А =А,пХ„, Ь = Ь-шХ ,. (27.1)
В § 1 мы рассмотрим задачу, в которой число строк А велико сравнительно с числом столбцов. В § 2 внимание сосредоточено на случае, когда А
160
имеет ленточную структуру. Методы, представленные в этих двух параграфах, применимы к задаче последовательного оценивания, так как при необходимости на каждом шаге последовательного накапливания легко можно вычислить вектор решения или отвечающую ему ковариационную матрицу. В § 3 демонстрируется применение этих методов на примере выравнивания посредством линейных сплайнов. В § 4 описана процедура среднеквадратичной аппроксимации посредством кубических сплайнов с равноудаленными узлами. Это еще один пример использования метода ленточного последовательного накапливания из § 2. Методы удаления данных рассматриваются в § 5.
§ 1. Последовательное накапливание
Опишем алгоритм преобразования матрицы [А : Ь] к верхнему треугольному виду, где не требуется, чтобы вся эта матрица находилась единовременно в памяти машины. Формальное описание этой процедуры дает алгоритм SEQHT (см. 27.10).
Договоримся вначале об обозначениях и приведем неформальное описание процесса. Запишем матрицу А н вектор Ь в блочном виде:
А,'
А =
Ьг
(27.2)
(27.3)
Здесь Aj — матрица размера mt X п, а Ь — вектор порядка mt. Разумеется, т - тх + ... + mq. Числаmt могут быть, в частности, равны 1, что позволяет получить большую экономию памяти.
В алгоритме строится последовательность треугольных*) матриц
[Rt : d{], i = 1.....q, с тем свойством, что задача наименьших квадратов
Rtx s dj имеет то же множество решений и ту же норму невязки, что и задача
Важное обстоятельство, позволяющее экономить память, заключается в том, что матрица [Rj : di ] может быть построена и хранима в массиве, где прежде находилась окаймленная матрица
[At:bt] =
ч
ft
(27.4)
*) При этом матрицы [Л,-: dt\, вообще говоря, неквадратшые. (Примеч. пер.) 11 Ч.Лоусон 161
Таким образом, максимальное число хранимых строк равно max {/И/ +
/-« /
+ гшп[л + 1, L ntf]} .
Удобно обозначить через [r0:d0] (пустую) матрицу, не имеющую строк. В начале /'-го шага алгоритма мы располагаем Щ-\ X (л+1)-матрицей [/?;_, :?//_!), вычисленной на /' — 1-м шаге, и новой матрицей данных \at:bi] размера т;Х(л+1). Положим ~mi=mi-l+mi. Составим окаймленную щ X (л + 1)-матрицу (27.4) и приведем зту матрицу к треугольному виду методом Хаусхолдера:
L 0:0 J} m,-m,
л _ ,
Здесь т{ = min {и + 1, mt).
Теперь 1-й шаг алгоритма завершен. Для] последующих ССЫПОК обозначим результат а шагов алгоритма через
л 1 ' .
Г r <Л} л , , .... ¦
Легко видеть, что существует ортогональная матрица Q такая, что
q[a
b] = \0 e\. Lo Oj
Следовательно, задача наименьших квадратов
(27.7)
эквивалентна исходной задаче (27.1) в том смысле, что обе имеют одинаковое множество решений и один и тот же минимум нормы невязки. Кроме того, матрица r имеет те же сингулярные числа, что и матрица а из (27.1).
Опишем теперь вычислительный алгоритм, который формализует зту процедуру. С этой целью положим
О, / = 0,
V, =
2 т„ />0, / = 1
(27.8)
ц = max {т/ + гшп[л + 1,Vi-i]} . (27.9)
к/<<?
Весь процесс обработки будет происходить в массиве W машинной памяти, состоящем по крайней мере из ц строк ил + 1 столбцов. Обозначение W(ii : t2, /i :/г) будет относиться к подмассиву массива W, образованному пересечением строк z'i.....i2 и столбцов /i, ... ,/2. Через W(i,j)
