CC BY
169
УДК 517.518.8
РЕГУЛЯРИЗАЦИЯ РЕШЕНИЯ ЗАДАЧИ ИДЕНТИФИКАЦИИ ПРИ ИСПОЛЬЗОВАНИИ АЛГОРИТМА ЧУВСТВИТЕЛЬНОСТИ
А.В. Димаки, А.А. Светлаков
Томский государственный университет систем управления и радиоэлектроники E-mail: dav18@yandex.ru
Предложен способ регуляризации плохо обусловленных систем линейных алгебраических уравнений, возникающих при использовании алгоритма чувствительности. Способ основан на применении изменяющегося от итерации к итерации значения параметра регуляризации, пропорционального ошибке между известными значениями функции и аппроксимирующим ее дифференциальным уравнением. Для реализации способа не требуется априорная информация о характеристиках известных значений. Показано, что скорость и область сходимости алгоритма чувствительности при использовании описанного способа не хуже, чем для метода регуляризации А.Н. Тихонова.
Ключевые слова:
Идентификация, выбор параметра регуляризации, алгоритм чувствительности.
1. Введение
Известна задача идентификации параметров обыкновенного дифференциального уравнения, используемого в качестве модели некоторого динамического объекта [1, 2]. При решении такой задачи зачастую возникает ситуация, когда получаемая система линейных алгебраических уравнений (СЛАУ) оказывается плохо обусловленной. Это приводит к неустойчивости и чрезмерной чувствительности решения системы к малым ошибкам задания коэффициентов при неизвестных и правой части СЛАУ [1]. Такая ситуация характерна, например, при использовании алгоритма чувствительности для решения задачи идентификации. Известным выходом из данной ситуации является применение того или иного метода регуляризации, в частности, метода, предложенного А.Н. Тихоновым [3], или других методов подобного назначения [4-8]. К недостаткам таких методов относятся либо необходимость подбора параметра регуляризации, либо требование наличия априорной информации о статистических характеристиках ошибок задания коэффициентов, в первую очередь о дисперсиях этих ошибок. Данные недостатки усугубляются еще и тем, что при использовании итерационной процедуры для решения задачи статистические характеристики ошибок могут изменяться от итерации к итерации, причем весьма существенно. В данной работе предложен, эвристически обоснован и апробирован способ регуляризации решения СЛАУ, свободный от указанных недостатков, т. е. не требующий априорной информации о характеристиках ошибок задания коэффициентов СЛАУ. Апробация способа проведена на примере решения одной из широко распространенных задач идентификации, представляющей значительный самостоятельный интерес.
2. Постановка задачи и описание
алгоритма чувствительности
Рассмотрим объект, описываемый обыкновенным дифференциальным уравнением, имеющим следующий вид:
^ХУ = I (х, а, у, у', у",...), (1)
сХд
где х - независимая переменная, а - «-мерный вектор неизвестных параметров, а у', у", ..., у(д-1) -первые д-1 производные переменной у по х. Здесь д и « - некоторые ограниченные натуральные числа такие, что д<п. Отметим, что дифференциальное уравнение (ДУ) д-го порядка (1) можно путем известных алгебраических преобразований свести к системе из д уравнений первого порядка. Как известно [9], и как нетрудно убедиться непосредственной проверкой, для этого необходимо и достаточно ввести д функций £;=,г;(х), определив их равенствами вида:
^/ёх = г+1,' = ^^;
ёг /
"Ах = 1(ха г, -1)-
Задача идентификации ДУ (1) сводится к нахождению таких оценок неизвестных компонент вектора а на основе имеющихся измерений х, у,, 7=1...т переменных хи у, чтобы решение полученного ДУ аппроксимировало результаты измерений наилучшим в некотором смысле образом [1]. Здесь т - натуральное ограниченное число такое, что т>п, а у определяется соотношением у=у,+е], где у, - истинное неизвестное нам значение переменной у, в] - неизвестное значение случайной ошибки е измерений, которая является случайной величиной, удовлетворяющей условиям:
М{е}=0; М{е2}=0 <«>, (2)
где М{.} - операция математического усреднения, а о2 - дисперсия величины е.
Одним из наиболее эффективных алгоритмов, позволяющих оценивать порядок и неизвестные параметры обыкновенного дифференциального уравнения, является алгоритм чувствительности [1]. В основу данного алгоритма положен хорошо известный в численном анализе метод линеаризации нелинейных функций, а также использование
так называемых функций чувствительности Wi по параметрам а, г=1...п.
Пусть у=у(а,х) - решение идентифицируемого дифференциального уравнения (1). Как известно [1], функции Щ определяются равенствами вида
ду(а, х)
^ (а, х) =-
да,
I = 1, п,
(3)
где а=(а1,а2,...,ап)Т - вектор-столбец неизвестных параметров размерности п; Т- символ операции транспонирования векторов и матриц. Из данных равенств видно, что 1-ая функция чувствительности Щ является частной производной решения у=у(а,х) дифференциального уравнения по г-му параметру а,, зависящей от вектора а и переменной х. Алгоритм чувствительности является итерационным. Число итераций, выполняемых при решении той ли иной задачи, заранее не фиксировано, а определяется по ходу решения задачи оценивания в зависимости от желаемой и достигнутой точности получаемых оценок.
Для описания последовательности действий, выполняемых на каждой итерации, будем считать, что: 1) даны т>п значений y¡ решения у=у(а,х), соответствующие значениям х переменной х; 2) уже выполнено к-1 итераций и 3) получен вектор оценок ак-1 параметров а,, г=1...п. Представим вектор ак новых оценок а,к параметров а, г=1...п, которые нам необходимо получить, равенством вида
а* = а*-1 + Да*. (4)
Из данного равенства видно, что для вычисления вектора ак необходимо и достаточно вычислить вектор поправок Dak. Будем считать, что компоненты вектора Dak являются достаточно малыми величинами и имеют место приближенные равенства
у (а*, х]) « у (а*-1, х]) + (а*-1, х] )Да* , ] = 1...т. (5) /=1
Конкретные численные значения поправок Дак нам неизвестны. Поэтому каждое из этих равенств является линейным алгебраическим уравнением относительно вектора Дак, а вся совокупность этих равенств - СЛАУ относительно данного вектора. Воспользовавшись векторами и матрицами, всю совокупность этих уравнений можно представить в следующем, традиционном для линейной алгебры виде:
W* Да* «Ду *,
где Wk и Дук - (тхп)-матрица и т-мерный вектор, элементы и компоненты Дук которых определяются следующими равенствами:
дУ(а*,х.)
да.
. = 1...т, I = 1...п;
Ду. = У(а,х.)-у(а",х.), } = 1...т.
При любом тхп данная система уравнений является, очевидно, переопределенной и, вообще говоря, несовместной. Для ее решения воспользуемся хорошо известным методом наименьших квадратов [10], полагая при этом, что все имеющиеся у
нас значения ¡=1...т являются равноточными, и, соответственно, содержащиеся в них ошибки измерения s¡ удовлетворяют условиям (2). Действуя в соответствии с данным методом и проведя ряд хорошо известных арифметических преобразований, получаем так называемую нормальную систему уравнений, которая в векторно-матричной записи имеет следующий, предельно компактный и традиционный в линейной алгебре вид:
В * Да* = Дг*. (6)
Здесь Дак - п-мерный вектор-столбец, компонентами которого являются искомые нами поправки, а п-мерный вектор-столбец Дък и (пхп)-матрица Вк вычисляются в соответствии с равенствами
Дг* = (W* )Т Д у*
и
В * = (W* )Т W*.
Решая систему (6) тем или иным способом и подставляя найденный вектор поправок Дак в равенство (4), получаем вектор ак новых оценок. На этом к-ая итерация реализованного алгоритма чувствительности заканчивается. Итерации, как уже отмечено выше, повторяются до тех пор, пока не будет достигнута требуемая точность аппроксимации имеющихся у нас значений ~ вычисляемыми значениями ~(ак,х,) решения у(а,х), либо изменение оценок параметров не станет меньше некоторой заранее заданной величины.
Важнейшей проблемой, возникающей при использовании алгоритма чувствительности, является тот факт, что СЛАУ (6) может быть плохо обусловленной, что негативно влияет на сходимость алгоритма чувствительности и не позволяет получить приемлемые по точности оценки ак вектора параметров а в ур. (1). Для решения указанной проблемы авторами был предложен рассматриваемый ниже способ регуляризации решения системы (6). Его работоспособность и пригодность для решения прикладных задач проиллюстрирована примером оценивания неизвестных параметров уравнения Пирсона, широко используемого при построении плотностей распределения случайных величин по их выборочным значениям [11].
3. Предлагаемый способ регуляризации
плохо обусловленных СЛАУ
Для регуляризации плохо обусловленных СЛАУ в настоящее время предложен целый ряд методов, базирующихся на различных идеях и подходах [3-8]. Отличительной особенностью всех этих методов является то, что применение любого из них требует привлечения той или иной априорной информации о решении исследуемой СЛАУ. Одним из наиболее простых и известных способов регуляризации является метод, предложенный А.Н. Тихоновым [12]. В рамках данного метода вместо решения СЛАУ (6) ищется решение другой СЛАУ, имеющей следующий вид:
(Вк + /Е)Аак = Агк
(7)
где у- параметр регуляризации - некоторое достаточно малое положительное число, а Е - единичная (ихи)-матрица.
Как показывают результаты исследований и накопленный опыт решения прикладных задач, чем меньше значение параметра у, тем ближе решения систем (7) и (6), и тем менее устойчиво это решение. Обратно, чем больше значение у, тем устойчивее решение системы (7), и тем больше отличается оно от решения системы (6). Отсюда вытекает, что: 1) значение /необходимо выбирать таким образом, чтобы система (7) оказалась достаточно хорошо обусловленной, и, вместе с тем, ее решение не очень сильно отличалось от решения системы (6); 2) выбор такого значения у является основной проблемой, возникающей при использовании рассматриваемого способа регуляризации плохо обусловленных СЛАУ. Одним из наиболее широко используемых способов решения данной проблемы является способ, в соответствии с которым выбор желаемого значения у сводится к решению нелинейного уравнения вида
Агк - ВкАак = И,
(8)
где ||.|| - евклидова норма векторов и матриц, а к -заданное достаточно малое положительное число, в качестве которого может быть использована и часто используется дисперсия о2. Данный способ является весьма трудоемким с точки зрения вычислительных затрат, что ограничивает возможности его применения, в частности, при решении задач управления динамическими объектами и т. д.
Таким образом, становится совершенно очевидной необходимость создания и применения какого-либо иного способа регуляризации, не требующего априорной информации о решении имеющейся СЛАУ и экономичного с точки зрения вычислительных затрат. Один из таких способов предлагается ниже, его сущность заключается в следующем.
Введем функцию определив ее равенством следующего вида:
1
^ =-Х (У у - У (ак, х, ))2
т
(9)
У=1
ной структуры дифференциального уравнения у=у(а,х). Данные особенности натолкнули авторов на идею замены постоянного значения параметра регуляризации у в уравнении (7) на вычисляемое на каждой к-й итерации значение функции (9). Таким образом, в рамках предложенного способа уравнение (6) приобретает следующий вид:
(Вк + Бк Е)Аак = Агк.
Важнейшей особенностью данного способа является использование на каждой итерации алгоритма чувствительности переменного значения параметра регуляризации /к, удовлетворяющего равенству вида
7к = б" , к = 1,2,...
Использование данного способа избавляет от необходимости решать на каждой итерации достаточно сложное и трудоемкое уравнение (8) относительно параметра уи, вместе с тем, как будет видно ниже, позволяет получить вполне приемлемые для практических приложений результаты.
4. Результаты тестирования предлагаемого
способа регуляризации
В качестве тестовой задачи для исследования влияния различных способов регуляризации решения СЛАУ (6) была использована задача аппроксимации плотности стандартного нормального распределения при помощи решения уравнения Пирсона [11]. Данная задача, как известно [11], сводится к нахождению оценок неизвестных параметров уравнения Пирсона, которое имеет следующий вид: ёу (х) = х - ао
У( х)
Ь0 + Ь1х + Ь2 х
-ёх,
где к - номер итерации алгоритма чувствительности.
Очевидно, что данная функция характеризует меру близости между измеренными значениями ~ функции у=у(а,х) и ее значениями ~(ак,х;), вычисленными на основе модели (5) в точках х на к-ой итерации. По мере того, как оценки компонент вектора параметров а будут приближаться к истинным значениям, значение функции У будет уменьшаться и стремиться к к.
При этом следует отметить, что значение функции ¡$к зависит не только от найденных на данной итерации оценок ак параметров а, но и от получен-
где у(х) - плотность распределения вероятности случайной величины X. Определив функции чувствительности Щ в соответствии с равенством (3), и, находя решения получившихся дифференциальных уравнений в каждой точке х, осуществляем итерации алгоритма чувствительности до выполнения условия
1 т
Б = -1 (Уу - У (ак, ху ))2 <£.
В настоящей тестовой задаче значение £ полагалось равным 10-5.
Известно [11], что для стандартного распределения Гаусса с параметрами (0; 1) истинные значения параметров уравнения Пирсона определяются равенствами:
а0=0; ¿0=1; ¿!=0; Ь2=0.
При решении тестовой задачи идентификации были выбраны следующие начальные значения оценок а0 параметров уравнения Пирсона:
а0=0,1; Ь0=0,9; ^=0; Ь2=0. (10)
Отличительной особенностью решаемой задачи является то, что от итерации к итерации численные
значения элементов матрицы Б^ могут изменяться на несколько (6 и более) порядков. Соответственно, при использовании метода регуляризации А.Н. Тихонова будем варьировать значение параметра регуляризации у так, чтобы оно оставалось «достаточно малым». Следует иметь в виду, что «достаточно малое» значение y на некоторой итерации может стать «весьма большим» на следующей и наоборот, что связано с изменением численных значений элементов матрицы Б*.
Оценки параметров уравнения Пирсона, найденные при различных значениях у с использованием метода регуляризации (7), приведены в табл. 1. Там же приведены значения чисел обусловленности C матрицы Б^, вычисленные на последней итерации алгоритма чувствительности. Число обусловленности C определялось в соответствии с равенством
C = | |Bk||| |(Bk )-1||.
Таблица 1. Результаты оценивания параметров уравнения Пирсона при различных значениях у. Начальные приближения: 3=0,1; Ь0=0,9; Ь=0; Ь2=0
у 30 b0 b b C Nter
0 -1,2210-5 0,971 2,7610-5 9, 6710-3 27300,31 5
10-4 -1,26.10-5 0,971 2,8410-5 9,6710-3 27206,91 5
10-3 -1,8210-5 0,971 4,0310-5 9,6810-3 26393,53 5
10-2 -2,6510-4 0,970 5,4310-4 1 , 01.10-2 20225,91 5
0,1 -2,1110-4 0 ,963 4,0810-3 1,310-2 5777,03 9
1 -3,7510-3 0,957 7,1410-3 1,5410-2 767,52 28
10 -3,9610-3 0,957 7,5210-3 1,57.10-2 83,93 235
100 5,3510-4 0,957 9,2310-4 1,3210-2 9,26 2550
Как видно из табл. 1, при идентификации плотности гауссова распределения, которое является одним из распределений Пирсона, значение параметра регуляризации влияет на сходимость алгоритма весьма слабо. Иными словами, алгоритм сходится и при у=0; при у>>0 число обусловленности матрицы Бк уменьшается, а число итераций, необходимых для достижения заданной точности, увеличивается. В целом, можно сделать вывод, что в данном случае применение каких-либо методов регуляризации не требуется.
Попытаемся теперь решить ту же самую задачу с другими начальными приближениями, например, заданными равенствами
а0=0,1; Ь0=0,9; Ь1=—0,1; Ь2=-0,1. (11) При этом вновь будем использовать метод регуляризации А.Н. Тихонова. Результаты решения данной задачи приведены в табл. 2.
В целом, приведенные результаты говорят о том, что в данном случае сходимость алгоритма чувствительности в большой степени определяется выбором начальных приближений параметров уравнения Пирсона. В том случае, когда начальные приближения далеки от истинных значений, алгоритм чувствительности начинает «зацикливаться» на некоторой промежуточной стадии решения, и дальнейшего увеличения точности не происходит,
и, как правило, монотонное уменьшение ошибки аппроксимации не наблюдается. Найденные на этом этапе оценки параметров могут значительно (на порядки) отличаться от истинных значений. Тот факт, что при у=0,01 и у=10 решение было получено, можно, по-видимому, считать случайным. Аналогичное поведение алгоритма чувствительности при использовании метода регуляризации (7) наблюдалось при решении ряда других тестовых задач, в которых варьировалось количество измеренных значений т, начальные приближения и тип аппроксимируемой плотности распределения.
Таблица 2. Результаты оценивания параметров уравнения Пирсона при различных значениях у. Начальные приближения: а0=0,1; Ь0=0,9; Ь1=-0,1; Ь~0,1
у 30 b0 b b'i C Nter
0 0 о Не удалось достичь заданной точности
10-2 -2,76.10- 4 0,970 5,17.10-4 0,01 20223,89 224
0,1 1 Не удалось достичь заданной точности
10 3,79.10-3 0,963 -4,6210-3 0,01 84,06 239
100 Не удалось достичь заданной точности
Применим теперь при тех же начальных приближениях способ регуляризации, изложенный нами в разделе 3. Результаты решения задачи приведены в табл. 3.
Таблица 3. Результаты оценивания параметров уравнения Пирсона при использовании предлагаемого способа регуляризации. Начальные приближения: 30=0,1; Ь0=0,9; Ь=-0,1; Ь=-0,1
30 b0 b b C Nter
-1,2710-5 0,971 2,2210-5 9,6710-3 27413,92 24
На начальных итерациях алгоритма чувствительности значение невязки (9), используемое нами в качестве параметра регуляризации, колеблется и незначительно убывает (рисунок). Аналогично ведет себя и число обусловленности С матрицы Бк. Резкий рост числа обусловленности на 11-й итерации (рис. 1, б) вызван тем, что произошла смена типа распределения, описываемого решением уравнения Пирсона. После этого как число обусловленности С, так и значение параметра регуляризации & монотонно убывает.
Как видим, достичь заданной точности аппроксимации (^=10-5) удалось за весьма небольшое количество итераций N¡¡¿=24. Найденные оценки значений параметров уравнения Пирсона близки к истинным значениям, а также к тем, что были найдены при начальных приближениях (10) (табл. 1) с использованием метода регуляризации (7).
Для оценки скорости сходимости алгоритма при использовании данного способа регуляризации проведем еще один численный эксперимент, в котором получим оценки параметров уравнения Пирсона при начальных приближениях (10). Результаты приведены в табл. 4.
Номер итерации к а
Номер итерации к б
Рисунок. Зависимость от номера итерации: а) значения параметра регуляризации 9; б) числа обусловленности матрицы Зк при решении задачи с начальными приближениями (11)
Таблица 4. Результаты оценивания параметров уравнения Пирсона при использовании предлагаемого способа регуляризации. Начальные приближения: ао=0,1; Ьо=0,9; Ь=0; Ь-=0
ao bo b Ь-2 C Nter
-1,22-105 0,971 2,75 10-5 9,67'10"3 27300,31 5
Полученные в результате данного теста (табл. 4) оценки параметров уравнения Пирсона близки к тем, что приведены в табл. 1, и также близки к истинным значениям параметров уравнения Пирсона для стандартного нормального распределения. Количество итераций, затраченных на получение данных оценок, не превышает количества итераций, которые потребовались для нахождения решения с использованием метода регуляризации А.Н. Тихонова, в тех случаях, когда его применение позволяло получить решение задачи.
Выводы
Предложен способ регуляризации плохо обусловленных СЛАУ, возникающих при использова-
нии алгоритма чувствительности. Способ основан на применении изменяющегося от итерации к итерации значения параметра регуляризации, пропорционального ошибке между известными значениями функции и аппроксимирующим ее дифференциальным уравнением.
Результаты тестирования предложенного авторами способа регуляризации показали, что: 1) данный способ позволяет расширить область сходимости алгоритма чувствительности; 2) для его применения не требуется никакой априорной информации о характеристиках ошибок измерения значений идентифицируемой функции; и 3) скорость сходимости алгоритма при использовании описанного способа регуляризации не хуже, чем при использовании метода регуляризации А.Н. Тихонова.
Следует заметить, что эффективность данного способа подтверждается лишь результатами его апробации на серии тестовых примеров. В дальнейшем необходимо разработать теоретическое обоснование эффективности предложенного способа регуляризации.
СПИСОК ЛИТЕРАТУРЫ
1. Рубан А.И. Идентификация нелинейных динамических объектов на основе алгоритма чувствительности. - Томск: Изд-во ТГУ, 1975. - 272 с.
2. Кондрашин А.В., Хорьков В.И. Исследование и идентификация управляемых технических систем. - М.: ИспоСервис, 2000. - 220 с.
3. Тихонов А.Н., Уфимцев Н.В. Статистическая обработка результатов экспериментов. - М.: Изд-во МГУ, 1988. - 174 с.
4. Воскобойников Ю.Е., Мухина И.Н. Локальный регуляризи-рующий алгоритм восстановления контрастных сигналов и изображений // Автометрия. - 2000. - № 3. - С. 45-53.
5. Hoang N.S., Ramm A.G. Solving ill-conditioned linear algebraic systems by the dynamical systems method // Inverse Problems in Science and Engineering. - 2008. - V. 16. - № 5. - P. 617-630.
6. Воскобойников Ю.Е. Численная реализация и сравнение четырех способов выбора параметра регуляризации в устойчивых
алгоритмах деконволюции // Научный вестник НГТУ. - 2004.
- № 2 (17). - C. 27-44.
7. O'Leary D.P. Near-optimal parameters for Tikhonov and other re-gularization methods // SIAM J. Sci. comput. - 2001. - V. 23. -№4. - P. 1161-1171.
8. Майстренко А.В., Светлаков А.А., Старовойтов Н.В. Регуляризация простейшего алгоритма цифрового дифференцирования сигналов // Научный вестник НГТУ. - 2006. - № 4. - С. 53-65.
9. Смирнов В.И. Курс высшей математики. Т. 2. - М.: Наука, 1974. - 656 с.
10. Ильин В.А., Позняк Э.Г. Линейная алгебра. - М.: Физматлит, 2002. - 318 с.
11. Кендалл М., Стьюарт А. Теория распределений. - Пер. с англ.
- М.: Наука, 1966. - 588 с.
12. Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. - М.: Наука, 1974. - 186 с.
Поступила 16.04.2009 г.
