CC BY
81
УДК 621.396
А. Л. Лебедев
РЕШЕНИЕ НЕКОРРЕКТНЫХ ЗАДАЧ МЕТОДАМИ МНОГОКРИТЕРИАЛЬНОГО МАТЕМАТИЧЕСКОГО ПРОГРАММИРОВАНИЯ
Рассмотрено решение некорректных задач, основанное на методах многокритериального программирования. В методах регуляризации, применяемых для решения некорректных задач, дополнительно определяется параметр регуляризации, для оценивания которого в настоящее время не существует однозначно формализованных процедур. Методы математического программирования позволяют избежать этой проблемы и учесть дополнительные условия задачи. В качестве примера рассмотрена задача о пеленгации точечных источников радиоизлучения с помощью антенной системы, состоящей из слабонаправленных элементов (вибраторов). Приведены алгоритмы вычисления точечных оценок амплитуд сигналов, пеленгов и углов места источников радиоизлучения; интервальные оценки получены методом статистических испытаний и аналитически.
При решении многих прикладных задач возникают ситуации, когда небольшие изменения исходных данных могут привести к значительным отклонениям в решении. Такие задачи относятся к классу некорректных [1]. Они возникают, например, при решении интегрального уравнения типа свертки; восстановлении входного сигнала по выходному сигналу и импульсной характеристике системы; определении интервалов между импульсами и их амплитуд при неточно заданном ядре интегрального уравнения; при сведении уравнения теплопроводности к интегральному уравнению; при решении задач пеленгации и другие.
Общепринятый метод решения некорректных задач — метод регуляризации. В функционал, точка минимума которого определяет точечную оценку решения, вводится компонента 71, отвечающая за близость решения к заданным результатам наблюдений (точкам), например 71 = ||Аг — мЦ^, и сглаживающий функционал (оператор) 72, который учитывает априорную информацию о возможном решении (например, в методе регуляризации Тихонова 72 = ||г||2). Здесь
п
||а||Г = |а»Г, ^ — размерность вектора а. Решение должно удовле-
г=1
творять двум целевым функциям 71 и 72 [1], которые в методе регуляризации объединены с помощью параметра регуляризации А в один функционал
7 (г) = 71 (г) + А7 (г), (1)
безусловный минимум которого надо найти. Метод регуляризации можно рассматривать как частный случай метода взвешенных сумм
— одного из способов перехода от задачи векторной оптимизации к задаче математического программирования. При этом вопрос о выборе оптимального параметра Л (как и весовых коэффициентов в методе взвешенных сумм) остается открытым, так как к настоящему времени не существует однозначно формализованных процедур его нахождения. В ряде случаев в методах регуляризации небольшие изменения параметра регуляризации могут привести к различным решениям. Таким образом, параметр регуляризации Л определяет выбор компромиссного решения и этот выбор не единственен.
Поскольку метод взвешенных сумм — только один из способов решения задач векторной оптимизации, то можно выделить и другие подходы, где не требуется определять Л. К ним, в частности, относятся [2, 3]:
1) сжатие области допустимых решений (метод е-ограничений);
2) целевое программирование (архимедова модель и модель с приоритетами).
При этом задачи векторной оптимизации позволяют учитывать любые дополнительные физические условия задачи, которые известны о возможном решении, например неотрицательность решения, добавлять другие целевые функции (минимизация энтропии и др.), чего нельзя сделать в методах регуляризации. Кроме того, после сведения задачи векторной оптимизации к однокритериальной задаче математического программирования, последнюю можно решать непосредственно или через построение и решение двойственной к ней задачи.
В качестве примера применимости методов векторной оптимизации для решения некорректных задач рассмотрим решение задачи пеленгации. Главная цель — показать, что параметр Л, который подлежит определению в методах регуляризации, здесь находить не нужно.
Задача состоит в определении пеленгов и местоположения источников радиоизлучения (ИРИ) в процессе мониторинга радиообстановки. Интерес представляет решение задачи, когда на одной частоте одновременно функционирует несколько ИРИ, т.е. когда надо найти пеленги каждого ИРИ. Один из возможных методов решения подобной задачи — метод регуляризации. Но, как уже было сказано, слабое место методов регуляризации — определение параметра регуляризации, для нахождения которого к настоящему моменту не существует строго формализованных процедур.
Рассмотрим пеленгацию нескольких неподвижных ИРИ, излучающих на одной частоте амплитудно- или частотно-модулированный сигнал при помощи антенных систем (АС), состоящих из набора слабонаправленных вибраторов (в общем случае вид сигнала не создает
принципиальных трудностей для решения задачи пеленгации). Требуется определить амплитуды и пеленги сигналов Zj (t), j = 1, 2,..., K, приходящих по направлениям (9j, 0j) к линии отсчета углов расположения вибраторов АС (здесь Oj — азимут, 0j — угол места). Рассмотрим математическую модель задачи пеленгации [4]
имхi (ti) = Amхк (O, 0) ZKxi (ti) + Пмxi (ti), (2)
где M — число детекторов; K — число источников радиоизлучения; u (t) — вектор сигналов с выходов элементов АС; z (t) — вектор сигналов модулирующих функций;
A (O, 0 )= [ a (Oi, в1) a №,£2) ... а (OK ,0K ) ] ;
O = [ Oi O2 ... Ok ]т; 0 = [ 0i 02 ... вк ]т.
В матрице A(0,0) m-й элемент, m =1, 2, ..., M, вектора-столбца a (Oj, 0j) имеет вид
am (Oj, 0j) = [exp {j (uo/c) (m - (M + 1)/2) dcos (Oj) cos (0j)}],
j = 1, 2, ..., K,
— для линейной АС;
am (Oj, 0j) = [exp {j (uo/c) (2пЛ/Л) cos (Oj - Ym) cos (0j)}],
j = 1, 2, ..., K
— для круговой АС.
Здесь u0 = 2nf0, f0 — несущая частота; c — скорость распространения электромагнитных волн; d — расстояние между соседними элементами АС; R — радиус окружности, вдоль которой расположены элементы АС; Л — длина волны сигналов ИРИ; Oj — пеленг j -го излучателя; 0j — угол места j-го излучателя; n (t) — вектор шума. Предполагается, что помехи являются постоянным в пространстве и времени белым шумом, не коррелирующим с излучателями, c нулевым математическим ожиданием и ковариационной матрицей а2/ (а
— среднеквадратическое отклонение (СКО); / — единичная матрица).
Введем новую переменную Xj = cos (Oj) cos (0j), отождествляя ее с cos (Oj): в нелинейных системах, зная величину Xj и сравнивая набег фаз на разных вибраторах, по простым формулам вычисляются отдельно и cos (Oj), и cos (0j). В линейной системе для этого достаточно хотя бы один из вибраторов сместить с линии расположения вибраторов. В примерах будет рассмотрена линейную АС.
Исходная система уравнений (2) относительно углов O и 0 и амплитуд, описывающих математическую модель задачи пеленгации, — нелинейная. Решить такую систему можно только в случае близких к
решению начальных приближений, иначе, как правило, итерационный процесс прекращается в локальном минимуме.
Сведем систему (2) к линейной путем расширения базиса, для чего в заданном секторе введем сетку по = , #2,..., #N]• Тогда определению подлежат только амплитуды z, соответствующие элементам сетки [4, 5]:
u (ti) = А z (ti) + n (ti). (3)
Тем самым сразу решается задача о числе ИРИ в заданном секторе. Размерность вектора z определяется шагом дискретизации по углу пеленгации.
Задача (3) — некорректная, так как сколь угодно малые изменения исходных данных могут приводить к произвольно большим изменениям решений. Некорректность задач объясняется тем, что уравнения, описывающие поведение систем, — "близки", для систем линейных алгебраических уравнений это выражается большим отношением значений максимального и минимального собственных чисел матрицы АтА. Обратная матрица (АтA)- ведет себя как фильтр верхних частот с высоким коэффициентом усиления, и сколь угодно малые отклонения в компонентах вектора u могут приводить к сколь угодно большим изменениям вектора решения z.
Пользуясь понятием регуляризующего функционала и считая его частным случаем метода взвешенных сумм для задачи векторной оптимизации, имеем две целевые функции, которые требуется одновременно минимизировать на некотором множестве D, определяемом условием u = Az (векторный критерий):
Ji (z) = ||Az — u||2 ^ min;
(4)
J2 (z) = ||z||p ^ min, 0 < p < то. В данной формулировке это задача двухкритериальной оптимизации.
Чтобы сформировать вторую целевую функцию для задачи двухкритериальной оптимизации, требуется определить вид сглаживающего функционала J2 (z). При дискретизации по углу пеленгации получают разреженный вектор решения z, большинство компонент которого равны нулю. Для такого решения целесообразно использовать lp-регуляризацию [5, 6]: J2 (z) = ||z||p, в частности при p =1 имеем
N
J2 (z) = £ |z |.
i=i
Задача двухкритериальной оптимизации для задачи пеленгации в данном случае имеет вид:
J1 = ||Az — u||2 ^ min;
z
1 ь ■ > 0 (5)
J2 = zi ^ min, zi ^ 0
i= 1
при ограничениях
Az = u.
Метод пороговой оптимизации (метод е-ограничений) приводит к различным возможным комбинациям целевых функций и ограничений
(p = 1):
1) min ||Az — u|2 при ||z|1 ^ (6a)
z2
2) min ||Az — u||2 при 11z11 < (66)
z2
3) min ||z||1 при ||Az — u|2 ^ (6в)
z
4) min ||z||1 при ||Az — u||2 ^ \fe. (6г)
z
Каждая задача, связанная с тем или иным физическим явлением или процессом, требует своего подхода. В настоящей статье рассмотрены задачи (6a) и (6в). Задача (6a) может быть решена методами квадратичного программирования [7, 8], задача (6в) — методами нелинейного программирования [9-11].
Для решения двухкритериальной задачи (5) применено целевое программирование [12]. В общем случае полученная таким способом однокритериальная задача решается методами нелинейного математического программирования.
В рассматриваемых примерах вычисляются точечные и интервальные оценки решения. Дисперсии оценок амплитуд получены методом статистических испытаний и аналитически. Аналитически ковариационная матрица оценок или дисперсий определяется путем обращения матрицы вторых производных логарифма функции правдоподобия, полученной на основе исходных данных для системы нелинейных уравнений, при найденных значениях точечных оценок [3].
В приведенных далее примерах на АС поступают два сигнала с пеленгами 9[ = 56025' и #2 = 128°57' и амплитудами 10 и 12 мВ соответственно (K = 2). Полагаем, что число элементов АС M = 16, d = А/2, f0 = 20 МГц. Подставив эти данные в исходную нелинейную систему уравнений (2), находим вектор u. Для каждого метода решения многокритериальной задачи выбрана своя оптимальная сетка по На
Пеленгационная панорама, полученная различными методами:
a — 1Р-регуляризация; б — квадратичное программирование; в — нелинейное программирование; г — архимедова модель; д — модель с приоритетами; е — минимизация энтропии
компоненты вектора и действует аддитивный гауссов шум с нулевым математическим ожиданием и СКО а = 0,5 мВ. Вектор и везде один и тот же. На рисунке (а) для сравнения с другими методами приведена пеленгационная панорама, полученная методом 1Р -регуляризации (р = 1, А = 1). Вычисленные по результатам статистических испытаний СКО решения приведены в сводной таблице.
Пример 1. Рассмотрим решение задачи (6a), т.е. функционал J2 (z) переведен в ограничения. При неотрицательных компонентах вектора z имеем задачу квадратичного программирования, методы решения которой хорошо разработаны [7, 8]:
J1 (z) = IIAz — u||2 ^ min;
z
J2 (z) = ||z|1 < S = J2опт (z) = 23,94;
z G D, z > 0,
где S — некоторая априорная оценка нормы решения (например, максимальная мощность сигнала, на которую рассчитан приемник). Результат решения задачи квадратичного программирования приведен на рисунке (б). Шаг по 9' выбран равным 1, 82°.
Точечные оценки решения (z1 и z2 — амплитуды сигналов, имеющих пеленг 91 и 92 соответственно) равны: z1 = 7,62 мВ; 9' = 55, 80°; z2 = 7,58мВ; 92 = 127,80°.
Ограничение на неотрицательность переменных было введено исходя из физического смысла задачи — амплитуды сигналов по определению не могут быть отрицательными. Если это ограничение снять, то результат решения окажется хуже — появятся дополнительные "шумовые" пики. Поэтому условие z ^ 0 всегда будем использовать в дальнейшем. Методы регуляризации этого сделать не позволяют.
Пример 2. Рассмотрим задачу (6в). Из двух функционалов J1 (z) = = ||Az — u||2 и J2 (z) = ||z||p функционал J1 (z) переводится в ограничения, что приводит к задаче нелинейного программирования [9-11]:
J2 (z) = ||z||p ^ min;
z
J1 (z) = || Az — u||2 < e = Jlопт (z) = 10,86;
z G D, z > 0,
где параметр e позволяет учитывать шум, присутствующий в компонентах вектора u. Примем e = 3a vM. Для p =1 и шага A9' = 1,21° получены следующие точечные оценки решения (см. рисунок (б)): z1 = = 8,60 мВ; 9' = 56,36°; z2 = 11,60 мВ; 92 = 128,49°.
Воспользуемся для решения той же задачи методом целевого программирования. При этом возможны два подхода: архимедова модель и модель с приоритетами [12]. Выбор значений S и e аналогичен их определению в задачах квадратичного и нелинейного программирования.
Пример 3. Архимедоба модель. Все целевые функции переводятся в ограничения и осуществляется максимизация взвешенной суммы меры их отклонений от ограничений:
max(w1d1 + w2d2) при J1 (z) + d1 ^ e, J2 (z) + d2 ^ S,
di,d,2
где ^ — весовые коэффициенты, ад^ = 1; — отклонения от
г= 1
ограничений.
На рисунке (г) приведен результат решения задачи при w 1 = 0,2; = 0,8; е = 10,86; 8 = 23,94; Д0' = 1, 45°. Выбор весовых коэффициентов влияет на решение. Увеличение w1 и соответственно уменьшение w2 приводит к более точному определению амплитуд сигналов, но при этом возникают "ложные" пики. Выбранное соотношение весовых коэффициентов позволяет усилить влияние функционала 72.
Получены точечные оценки решения: = 6,66 мВ, 0^ = 55,96°; г2 = 11,71 мВ, ¿2 = 128,49°.
Найденные отклонения: = 0; = 2,01.
Пример 4. Модель с приоритетами. Проводится последовательный перевод целевых функций в ограничения и максимизация отклонения значений целевых функций от ограничений. При этом найденное на данном шаге значение отклонения используется как оптимальное отклонение на следующем шаге:
71 (г) + ¿1 < е;
71 (г) + ¿1сит|й1оп1=й1 = ^ 72 (г) + ¿2 < 8.
На рисунке (д) приведен результат решения задачи при е = 10,86; 8 = 23,94; шаг Д0' = 0,81°.
Точечные оценки решения: ¿1 = 6,72 мВ; = 55,96°; г2 = 11,54 мВ; 02 = 128,49°. Найденные отклонения: =0; = 2,14.
Пример 5. Отличительная особенность задач математического программирования — это возможность добавлять к исходной задаче любые необходимые ограничения и целевые функции или переводить последние в ограничения. В задачу (5) введем дополнительную це-
N
левую функцию — энтропию 73 = 1п г |, получив следующую
задачу:
i=1
Ji = ||Az—иУ2 ^ min;
z
N
J2 = S Zi ^ min; (7)
i=i z v >
N
J3 = S |Zi ln Zi| ^ min
i=i z
при ограничениях
Az = u, zi ^ 0.
Рассмотрим случай, когда функционалы 7 и 72 переведены в ограничения (метод пороговой оптимизации), т.е.
N
штУ^ 1пг| при ||Аг — и||2 ^ е, Ц^к ^
г=1
Поскольку 1п0 не существует, то в качестве нижней границы для м выбрано значение 0,01, т.е. г ^ 0,01. Результат решения приведен на рисунке, е, шаг Д0' = 0,81°. Получены точечные оценки решения: г1 = 6,61 мВ; в[ = 55,96°; г2 = 11,91 мВ; ^ = 128,49°.
Пример 6. Рассмотрим случай, когда приходит один сигнал с амплитудой 5 мВ, пеленгом 55° и углом места 30°. В линейной АС, состоящей из 16 вибраторов, сдвинем последний вибратор в сторону на 60° — угол между вертикальной плоскостью линейной АС и вертикальной плоскостью, проведенной через фазовый центр и смещенный вибратор. По-прежнему считаем, что на вектор и действует аддитивная помеха с нулевым математическим ожиданием и СКО а = 0,5 мВ. При решении методами нелинейного программирования приходим к системе
cos в cos в = 0,5; ( Л
cos (в + y) cos в = -0,3746. (8)
Решая систему (8), получаем искомые пеленг в = 55,27° и угол места в = 28,65°.
Для определения ковариационной матрицы оценок, а иногда и для уточнения точечных оценок, полученных методами математического программирования, можно использовать исходную нелинейную систему уравнений (2), взяв в качестве начального приближения найденные точечные оценки решения. Такое начальное приближение близко к оптимальному решению и приводит в точку глобального минимума. Действительно, поскольку на вектор u действует аддитивный гауссов шум с СКО, равным а, и нулевым математическим ожиданием, то u подчиняется многомерному нормальному закону распределения вида
f (u) = 1
u-A(e')z 2
/г) 2\М/2
(2па2) '
где А (0') г — математическое ожидание вектора и. Объединим эту информацию в виде функции правдоподобия. Вычислив матрицу, обратную матрице вторых производных логарифма функции правдоподобия при найденных точечных оценках решения системы (3), согласно работе [3] можно найти ковариационную матрицу решений. Например,
для метода квадратичного программирования (пример 1) аналитический расчет дает следующие оценки СКО:
^ан (¿1) = 0,05 мВ; аан (¿2) = 0,04 мВ; ^ан (5) = 0,05°; аан (<%) = 0,06°.
В таблице приведена сводка основных характеристик полученных решений, из которой видно, что решение методом /р-регуляризации, хорошо согласуется с решением, найденным методами многокритериального математического программирования. Полученные оценки СКО характеризуют точность только математических алгоритмов. Предлагаемый метод решения задачи пеленгации показал свою эффективность и лишен главного недостатка, присущего методам регуляризации, — необходимости определения параметра регуляризации.
По результатам данных, приведенных в таблице, можно сделать выводы:
1) минимальные интервальные оценки для пеленгов сигналов обеспечивает метод сжатия области допустимых значений с функционалом
— и||2 в ограничениях;
2) самым быстродействующим является метод сжатия области допустимых значений с функционалом ||г||р в ограничениях (является более быстродействующим, чем метод /р-регуляризации);
3) минимальные интервальные оценки для амплитуд сигналов дает архимедова модель;
Сводная таблица СКО
Метод град ), град мВ мВ Время работы (отно сительные единицы)
¿р -регуляризация 0,10 0,13 0,18 0,22 7
Квадратичное программирование 0,12 0,15 0,21 0,36 3
Нелинейное программирование 0,04 0,04 0,24 0,32 17
Архимедова модель 0,11 0,14 0,11 0,12 14
Модель с приоритетами 0,09 0,15 0,23 0,37 50
Минимизация энтропии 0,07 0,12 0,20 2,33 11
4) задача минимизации энтропии позволяет учесть все требования и дает средние показатели.
СПИСОК ЛИТЕРАТУРЫ
1. Тихонов А. Н., А р с е н и н В. Я. Методы решения некорректных задач.
- М.: Наука, 1979. - 285 с.
2. Плотников В. Н., Зверев В. Ю. Принятие решений в системах управления. Ч. 1. Теория и проектирование алгоритмов принятия оперативных решений: Учеб. пособие. - М.: Изд-во МГТУ, 1993. - 172 с.
3. Грешилов А. А. Математические методы принятия решений: Учеб. пособие для вузов. - М.: Изд-во МГТУ им. Н.Э. Баумана, 2006. - 584 с.
4. Грешилов А. А., Назаренко Б. П., Плохута П. А. О пеленгации источников излучений // Вестник МГТУ им. Н.Э. Баумана. Сер. "Естественные науки". 2007. - № 3.
5.Malioutov D. M. A sparse signal reconstruction perspective for source localization with sensor arrays / Master of Science thesis, Massachusetts Institute of Technology. - 2003. - P. 35-59.
6. C e t i n M., Karl W. C. Feature-enhanced synthetic aperture radar image formation based on nonquadratic regularization // IEEE Trans. Image Processing.
- 2001. - Vol. 10, No. 4. - P. 623-631.
7. Coleman T. F.,Li Y. A reflective Newton method for minimizing a quadratic function subject to bounds on some of the variables // SIAM Journal on Optimization.
- 1996. - Vol. 6, No. 4. - P. 1040-1058.
8. G i 11 P. E., M u r r a y W., Wright M. H. Practical optimization. - London: Academic Press. UK, 1981.
9. C o l e m a n T. F., L i Y. An interior, trust region approach for nonlinear minimization subject to bounds // SIAM Journal on Optimization. - 1996. - Vol. 6. -P. 418-445.
10. C o l e m a n T. F., L i W. On the convergence of reflective Newton methods for large-scale nonlinear minimization subject to bounds // Mathematical Programming.
- 1994. - Vol. 67, no. 2. - P. 189-224.
11.Han S. P. A globally convergent method for nonlinear programming // Journal of Optimization Theory and Applications. // - 1997. - Vol. 22. - P. 297.
12. ШтойерР. Многокритериальная оптимизация. Теория, вычисления и приложения / Пер. с англ. - М.: Радио и связь, 1992. - 504 с.
Статья поступила в редакцию 1.04.2008
Алексей Леонидович Лебедев родился в 1983 г., окончил в 2007 г. МГТУ им. Н.Э. Баумана. Аспирант кафедры "Высшая математика" МГТУ им. Н.Э. Баумана. Специализиуется в области методов решения некорректных задач.
A.Ye. Lebedev (b. 1983) graduated from the Bauman Moscow State Technical University in 2007. Post graduate of "Higher Mathematics" department of the Bauman Moscow State Technical University. Specializes in the field of methods of solving ill-posed problems.
