Регуляризация неустойчивых задач наименьших квадратов на основе расширенных систем Текст научной статьи по специальности «Математика»

Вычислительные технологии

Том 12, № 6, 2007

РЕГУЛЯРИЗАЦИЯ НЕУСТОЙЧИВЫХ ЗАДАЧ НАИМЕНЬШИХ КВАДРАТОВ НА ОСНОВЕ РАСШИРЕННЫХ СИСТЕМ

А. И. ЖДАНОВ, Т. Г. ПАРЧАйКИНА Самарский государственный аэрокосмический университет им. акад. С.П. Королева, Россия e-mail: [email protected], [email protected]

The method of the solving of unstable least square problems by crude data is considered. This method is based on the transformation of least square problem to equivalent augmented system of linear algebraic equations with symmetric matrix.

Постановка задачи

Рассмотрим задачу наименьших квадратов в следующем виде [1]:

min ||/- Л«||2 + 2cTи, (1)

где A £ Rmxn, / £ c £ T — знак транспонирования, ||/||2 — евклидова векторная норма.

В [1] рассмотрен случай, когда m > n и rank A = n, а исходные данные d = {A, /, c} известны точно. В данной работе такие ограничения не вводятся. Случай c =0 был рассмотрен в [2].

Информация о задаче (1) содержит приближенные данные d = {A, /, с}:

||A - A||2 < h, ||/ -/||2 < i', ||c - с||2 < i", (2)

где h > 0, 8' > 0, 8" > 0 характеризуют погрешности задания приближенных данных d. В дальнейшем в качестве матричной нормы понимается спектральная матричная норма, т. е.

||A - a4|2 = sup ||Au - a4m|2. У«У2 = 1

Таким образом имеется семейство задач:

min ||/с - a4m|2 + 2сти, (3)

удовлетворяющих (2).

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2007.

Задача (1) эквивалентна решению системы линейных алгебраических уравнений (СЛАУ) [3]:

АтАи = Ат f - с, (4)

которая, в свою очередь, эквивалентна расширенной системе уравнений:

V + Аи = f, Ат V = с

(5)

или

Am A)(0 = (0 « gz = b- (6)

где Im — единичная матрица порядка m.

На основании теоремы Кронекера—Капелли система (4) совместна тогда и только

тогда, когда rank(ATA) = rank (ATA . ATf — c). В силу (5) условие совместности равносильно выполнению условия rank AT = rank (AT . c). Таким образом, можно сформулировать следующее утверждение.

Утверждение 1. Задача (1) разрешима тогда и только тогда, когда выполнено условие rankAT =rank(AT : c).

Решение задачи на основе расширенных систем

В общем случае под "решением" задачи (1) будем понимать единственное решение с минимальной нормой

u* = (AT A)+(AT f — c),

где (ATA)+ — псевдообратная матрица, или обобщенная матрица Мура—Пенроуза.

Как отмечалось в [4], решение системы (4) по приближенным данным d при h > 0 является некорректно поставленной по Адамару задачей, так как приближенное решение

и* = (AT A)+(AT f — c) (7)

неустойчиво к бесконечно малым возмущениям исходных данных. Задача (3) эквивалентна системе

AT Au = AT f — с (8)

и, в силу (6), расширенной системе уравнений

At A)(U) = (С) « Gz=ь. (9)

Преобразование задачи (3) к эквивалентной задаче (9) дает возможность существенно понизить число обусловленности исходной вычислительной задачи.

Утверждение 2. Единственное решение с минимальной нормой расширенной системы (9) находится как С* = G+b = (С, Uf)T, где U* определяется выражением (7)

и v* = f — AU* :

Для нахождения устойчивых решений задачи (3) воспользуемся методом регуляризации Тихонова [5].

Так как О = От, то ОтО = О2. Следовательно, регуляризованное решение системы (9) определяется как (единственное) решение уравнения Эйлера:

(О2 + а/т+„),г = О Ь . (10)

Возмущения в приближенной системе {О, Ь} удовлетворяют неравенствам

||G - G|h

0 A - A AT - AT 0

IIA - a||2 < h,

||b - b|2 = If - f ||2 + ||c - c|2 < ¿'2 + ¿"2 ^ ||b - b|2 < i,

где 6 = л/6'2 + 6''2. Из утверждений 1 и 2, теоремы 4 в [4] непосредственно следует теорема.

Теорема 1. Если в уравнении (10) положить а = Н, тогда ошибка уклонения

- z*||2 = O(h + ¿),

где fa — решение (10).

При численном решении уравнения Эйлера (10) для расширенных систем важно исследовать числа обусловленности этой системы. В [2] показано, что спектр матрицы

^m A

AT 0

имеет вид

Л

1 2

^/т + а2,

k = 1, 2,..., т;

1, имеет кратность m - т; 0, имеет кратность n - т,

где а > о >

> ат > ат+1 = ... = ап = 0 — сингулярные числа матрицы A, т.е. afc = \JЛк(ATA), Лк(ATA) — собственные числа матрицы ATA, k = 1, 2,... , n; т = rank (A).

Исследуем обусловленность регуляризованной системы (уравнения Эйлера) (10). Для этого оценим число обусловленности матрицы G.

Пусть матрица G — матрица полного ранга, т. е. f = rank (A) = n. Тогда, если f = n, минимальное сингулярное число матрицы A an = 5"mjn(A ) > 0 и

cond2(G2 + a/m+n) < cond.2(G2) = I 1 + + 2 ( О1 ) | < ^v/2cond2(A) + 1)2

Пусть матрица G — матрица неполного ранга, т. е. f = rank (A ) < n. Тогда получаем

• i-\ 2

1 Я

Лтах^ + a/m+n) + О? + а < 1 + О + О + a.

2

2

При а = к для максимального собственного значения матрицы С2 + а/т+п имеет место

а

оценка

, г \ ^ 1 , ~ , ~ 2

Атах(С + а/т+га) < 1 + а1 + + к.

Очевидно, что при Г < п Ат;п(С) = 0 и, следовательно, Ат;п(С2 + а/т+п) = а. Таким образом, при а = к, Г < п

со^2 (С2 + а/т+га) < 1 +

1 + д1 + д

к

Регуляризация систем на основе мнимого сдвига спектра

Используя свойство симметричности матрицы С, можно понизить число обусловленности регуляризованной задачи (10). Для этого воспользуемся методом мнимого сдвига спектра (В.Н. Фаддеевой).

В соответствии с этим методом вместо уравнения (10) рассмотрим уравнение

(С + г\/а1т+п)г = Ь, (11)

где г — мнимая единица, г2 = -1.

Пусть г = х + гу. Тогда (11) можно записать в виде расширенной вещественной СЛАУ: _ _

С -\/а1т+п \ ( х \ _ / Ь

- — • (12) С \у \0 '

Из (12) непосредственно следует, что х = Ке г является решением уравнения

(С2 + а1т+п, )х = СЬ .

Таким образом, задачу решения СЛАУ (10) можно заменить на задачу решения СЛАУ (11) с мнимым сдвигом спектра.

Исследуем число обусловленности матрицы Ш = С + Собственные значе-

ния Ак (Ш) матрицы Ш равны

Ак (Ш ) = Ак + г^а, к = 1, 2, ...,т + п,

где Ак = Ак (С) — собственные значения матрицы С.

Тогда сингулярные числа дк (Ш) матрицы Ш определяются как

дк (Ш) = у/Ак (С2 + а/т+га) = у]А2 + а.

Следовательно,

сопа2(С2 + а/т+га) = сопа2 (С + ^л/а1т+п)

и при Г = п

coпd2(С + < соп^(С) < ^/2coпd2(A) + 1.

Аналогично, если A неполного ранга, а a = h, то

cond2(G + ¿v^im+n) < ]J 1 + + ^ + 1.

Таким образом, из последних двух неравенств видно, что для плохо обусловленных матриц A (cond2(A ) ^ 1) вычисление регуляризованных решений za на основе СЛАУ (11) с мнимым сдвигом спектра существенно эффективнее, чем на основе СЛАУ (10).

Примеры

Рассмотрим задачу (1), если

2 -10 \ / 18 \ / 18 4=1 -1 1 2 J ' f = ( — J ' c = ( -9

Так как rank AT = rank (AT.c) = 2, то решение задачи (1) существует и имеет вид

-1

u = (ATA)+(ATf — c)

1

В элемент матрицы а13 введем возмущение такое, что |а13 — й13| = 0.0001. Очевидно, погрешности приближенных данных определяются величинами

||A — A||2 < 0.0001.

Приближенная система уравнений (8), которая определяет нерегуляризованную задачу:

5 —3 —0.9998 \ / u1 \ / —9 —3 3 2.9999 U2 = 9 | . (13)

—0.9998 2.9999 5.00000001 / \ щ / \ 9.0018

Система (13) невырожденная: det AT A = 0.00000001 = 0. Любым классическим (машинным) методом найдем ее решение:

— 1.26 ■ 106 — 2.52 ■ 106 1.26 ■ 106

Регуляризованные решения ua определялись из СЛАУ (12). На основании теоремы параметр регуляризации a был выбран a = h. Решая систему (12) любым классическим методом, получаем регуляризованное решение:

—0.99999 ' 0.99963 1.0006

h 10-3 10-5 10-6 10-7 10-8 10-9 10-10

||Ua - «*||2 6■10-2 6■10-4 6■10-5 6■10-6 6■10-7 6■10-8 6■10-9

||й* — U*||2 3 ■ 105 1.7 ■ 106 1.7 ■ 107 1.7 ■ 108 1.5 ■ 109 1.5 ■ 108 1.5 ■ 107

Сравнивая регуляризованное решение с точным решением и*, получаем

Ни« — и*||2 < 6 ■ 10-3,

что соответствует полученным в работе теоретическим результатам.

Приведем результаты исследования данного метода при различных к (см. таблицу).

Метод регуляризации на основе эквивалентных расширенных систем, предлагаемый в настоящей работе, позволяет решать (с гарантированной точностью) неустойчивые задачи наименьших квадратов по приближенным данным (3). Отметим, что за счет симметричности матрицы О удается понизить число обусловленности регуляризован-ной вычислительной задачи, при этом вычислительная трудоемкость решения СЛАУ (12) не намного выше трудоемкости решения СЛАУ (10). На основании теоремы, приведенной в данной работе, параметр регуляризации выбирается а = к. Таким образом, не требуется специальных правил выбора параметра регуляризации, которые существенно усложняют решение исходных некорректных задач.

Список литературы

[1] BjoRK A. Pivoting and stability in Augmented System Method // Proceedings of the 14th Dundee Conference. 1991. P. 1-16.

[2] ЖДАНОВ А.И. Регуляризация неустойчивых конечномерных линейных задач на основе расширенных систем // Журн. вычисл. матем. и математической физ. 2005. Т. 45, № 11. С. 1918-1926.

[3] Bjo RK A. Numerical stability of methods for solving augmented system / / Contemporary Math. 1997. Vol. 204. P. 51-60.

[4] МОРОЗОВ В.А. Алгоритмические основы методов решения некорректных задач // Вычисл. методы и программирование. 2003. Т. 45. С. 130-141.

[5] Тихонов А.Н. О приближенных системах линейных алгебраических уравнений // Журн. вычисл. матем. и математической физ. 1980. Т. 20, № 6. С. 1373-1383.

Поступила в редакцию 29 марта 2007 г., в переработанном виде — 28 сентября 2007 г.