ПЕРЕОПРЕДЕЛЕННЫЕ СИСТЕМЫ ЛИНЕЙНЫХ УРАВНЕНИЙ. МЕТОД НАИМЕНЬШИХ КВАДРАТОВ
Работа позволяет исследовать и продемонстрировать особенности нахождения обобщенного решения переопределенных систем линейных уравнений. Рассматриваются переопределенные системы линейных уравнений, возникающие, например, при обработке результатов физических экспериментов. Исследуется влияние выбора базиса (степенные и тригонометрические многочлены, многочлены Лежандра), степени многочлена на свойства вычислительного алгоритма и обусловленность возникающей задачи.
Каноническая запись переопределенной системы линейных алгебраических уравнений имеет следующий вид:
........................................... n > s. (6.1)
Введем
пространства и
состоящие из
элементов вида
и имеющие размерности s и n > s соответственно. Обозначим через А прямоугольную матрицу
системы (6.1):
(6.2)
Тогда систему (6.1) можно записать в виде
(6.3)
Введем в «основное» скалярное
произведение, положив
(6.4)
Скалярное произведение в можно ввести множеством
других способов. Именно, произвольной симметричной и положительно определенной
матрице
т. е.
для любого
вектора
соответствует
скалярное умножение
(6.5)
Известно, что
любое скалярное произведение в пространстве можно записать формулой (6.5), подобрав
соответствующий самосопряженный оператор
Система (6.1),
как правило, не имеет классического решения, т. е. не существует такого
набора чисел который
обращает каждое из n уравнений (6.1) в тождество.
Определение.
Фиксируем В:
Введем
функцию от
положив
(6.6)
Примем за обобщенное решение
системы (6.1) вектор придающий наименьшее значение
квадратичной форме (6.6).
Замечание. Выбор зависит от исследователя. Матрица
В имеет смысл «весовой» матрицы и выбирается из тех или иных соображений
о том, какую цену придать невязке системы (6.1) при заданном
Теорема 1. Пусть столбцы матрицы А линейно независимы, т. е. ранг матрицы А равен s. Тогда существует одно и только одно обобщенное решение b системы (6.1). Обобщенное решение системы (6.1) является классическим решением системы уравнений
(6.7)
которая содержит s скалярных
уравнений относительно s неизвестных
В дальнейшем будем иногда использовать обозначение
Переопределенную
систему где
можно записать в виде:
где
— i-й столбец матрицы A,
а вектор
Требуется
найти коэффициенты линейной
комбинации
так, чтобы эта линейная
комбинация наименее отличалась от f:
Обозначим через подпространство
размерности s пространства
состоящее из всевозможных линейных комбинаций
векторов
…,
Пусть — обобщенное решение
переопределенной системы. Тогда линейная комбинация
— ортогональная в смысле
скалярного умножения
проекция
вектора
на
подпространство
так
как любой вектор из
имеет
вид:
Наименее
уклоняется от f элемент подпространства
имеющий вид
где
— решение системы (6.7) методом наименьших
квадратов (МНК).
В силу элемент
ортогонален любому
элементу
Если в пространстве вместо базиса
…,
выбрать какой-либо другой базис
…,
то система (6.7) заменится системой
(6.8)
с матрицей где
и правой частью i-я
компонента которой
Вместо решения
системы (6.7)
получим новое решение
системы
(6.8), но проекция f на
останется прежней.
Если нас интересует проекция заданного вектора f
на заданное подпространство то естественно стремиться к выбору базиса
…,
этого подпространства, по
возможности мало отличающегося от ортонормированного. Искомая проекция от выбора базиса в
не зависит, а система МНК (6.7) в
случае такого базиса будет иметь хорошо обусловленную матрицу.
Обусловленность линейной системы (6.7)
определяется числом которое определяет относительную погрешность
решения системы в зависимости от погрешности правой части f:
где векторы y и z определяются из решения следующих двух систем уравнений:
где — транспонированная матрица, е —
вектор с компонентами
выбираемые
так, чтобы обеспечить максимальный рост
на этапе обратной подстановки.
Реализованный в программе прямой
метод решения системы МНК (6.7) является вариантом метода Гаусса
последовательного исключения неизвестных с частичным выбором ведущего элемента
(по столбцу). Исключение по методу Гаусса состоит из двух этапов: прямого
хода и обратной подстановки. В прямом ходе на k-м шаге k-е
уравнение вычитается из оставшихся с целью исключения k-ого
неизвестного. Обратная подстановка состоит в решении последнего уравнения
относительно предпоследнего
— относительно
и
т. д. до
Это прямой метод, но может применяться и как итерационный. Подробнее см. работу [4] и библиографию к ней. Используется, если матрица линейной системы самосопряженная и положительно определенная. Матрица системы МНК как раз такая (см. п. 6.2). Метод неустойчив, если обусловленность системы достаточно велика. Обозначим матрицу системы МНК через C.
Существенным преимуществом метода является отсутствие необходимости знания границ спектра. В точной арифметике метод сходится не более чем за n итераций, где n — размерность матрицы.
Многочлены
…,
i = 3, 4, … (6.10)
называют многочленами
Лежандра. Они ортогональны на отрезке т. е.
k ¹ l;
В программе используются многочлены, ортогональные
на отрезке [a, b],
получающиеся из (6.10) соответствующей
линейной заменой переменных;
i = 1, …, N —
абсциссы точек исходных данных.
Изложение теоретических основ использования МНК ведется в соответствии с [1], см. также [8, 11, 27].
1. В окне «Данные и их аппроксимация» маркерами
отмечены точки с координатами k = 1, 2, …, 5, записанные
в файле данных. По умолчанию это файл DATA\LSQ1.DAT.
В этом файле
записаны результаты замеров какой-либо зависимости Расположение точек позволяет
приближенно считать, что зависимость
является линейной, т. е. имеет вид
Числа
желательно подобрать так, чтобы
при
получались
значения
Так
как измерений больше, чем неизвестных
то соответствующая система линейных уравнений
будет переопределенной и в общем случае не имеет классического решения,
поскольку не существует прямой, проходящей сразу через все экспериментальные точки.
Для
аппроксимации данной зависимости в меню «Параметры/Выбор системы
функций» выберите строку «Полиномы 1, x,
…» и в
строке «количество базисных функций» установите число 2.
Если считать все измерения равноправными, то весовую матрицу В следует выбрать единичной. Для этого в меню «Параметры/Выбор Евклидовой нормы» выберите строку «Матрица В Единичная».
Если Вы сочтете разумным придать каждому измерению свой вес, то в меню «Параметры/Выбор
Евклидовой нормы» выберите строку Матрица B произвольная
и строку «редактирование матрицы». Задайте диагональным элементам В
неравные значения. Запустите задачу на счет и сравните результаты расчетов.
Коэффициенты
можно посмотреть,
активизировав меню «Окна/Полная информация».
Если заранее
известен вид зависимости то тот же набор данных можно
интерпретировать как результат нескольких измерений для уточнения коэффициентов
2. Выберите
базис из полиномов Лежандра («Параметры/Выбор системы функций/Полиномы
Лежандра») для представления того же
набора экспериментальных данных т. е. будем искать зависимость y(x) в виде
Установите, зависит ли аппроксимирующий полином от выбора базиса
(см. п. 6.2).
3. Введите
набор данных из файла LSQ2.DAT («Параметры/Данные/Ввод данных
из файла»). Он отвечает более сложной, чем в предыдущем задании,
зависимости и
для хорошей аппроксимации нужно выбрать полином более высокой степени.
Установите число базисных функций 10 и базис-полином 1, x,
… в меню «Параметры/Выбор системы
функций».
Для решения системы МНК в меню «Метод» выберите «Метод сопряженных градиентов» и установите число итераций 8.
Следующий расчет выполните для базиса из полиномов Лежандра. Сравните полученные результаты. Обратите внимание на оценку числа обусловленности в обоих случаях. Опишите, как изменяются результаты, если уменьшить точность промежуточных вычислений. Для этого в меню «Параметры» в строке «Длина мантиссы» введите вместо установленной по умолчанию (53) цифру 24, а затем 12.
4. Данные в LSQ3.DAT соответствуют замерам зашумленного сигнала
Найти
Указание. Используйте базис из тригонометрических полиномов («Параметры/Выбор системы функций/Полиномы тригонометрические»).
5. Введите таблицу функции из файла LSQ4.DAT. Приблизьте эту функцию тригонометрическими многочленами степени m и m + r, r > 0, m + r < N, где N — число точек в таблице. Сравните первые m коэффициентов. Объясните полученный результат.
6. Используя мышь, введите произвольный набор данных с экрана терминала («Параметры/Данные/Ввод данных с экрана»). Аппроксимируйте их полиномами различной степени. Проанализируйте влияние выбора базиса, выбора матрицы В, задающей скалярное произведение, способы решения системы МНК («Метод Гаусса, Метод сопряженных градиентов»).