Стохастическое оценивание параметров нелинейных дискретных объектов на основе обобщенных вероятностных критериев Текст научной статьи по специальности «Электротехника, электронная техника, информационные технологии»

ISSN 0868-5886 НАУЧНОЕ ПРИБОРОСТРОЕНИЕ, 2011, том 21, № 2, c. 112-119 МАТЕМАТИЧЕСКИЕ МОДЕЛИ, ОБРАБОТКА ДАННЫХ -

УДК 519.248: 681.5.001.3

© П. А. Кучеренко, С. В. Соколов

СТОХАСТИЧЕСКОЕ ОЦЕНИВАНИЕ ПАРАМЕТРОВ НЕЛИНЕЙНЫХ ДИСКРЕТНЫХ ОБЪЕКТОВ НА ОСНОВЕ ОБОБЩЕННЫХ ВЕРОЯТНОСТНЫХ КРИТЕРИЕВ

Показана актуальность исследования альтернативных (по отношению к традиционным) методов стохастического оценивания параметров нелинейных дискретных объектов. Впервые предложен алгоритм оптимального оценивания параметров математических моделей дискретных нелинейных стохастических объектов на основе критерия минимума вероятности ошибки оценивания. Рассмотрен модельный пример, иллюстрирующий эффективность реализации предлагаемого алгоритма.

Кл. сл.: нелинейный дискретный объект, стохастическое оценивание параметров, обобщенные вероятностные критерии, критерий минимума вероятности ошибки оценивания

ВВЕДЕНИЕ

Стохастическое оценивание параметров математических моделей объектов и систем является на сегодняшний день одним из важнейших направлений современных исследований для целого ряда практически значимых областей научной деятельности человека — управления различными системами и процессами, имитационного моделирования, метрологии, приборостроения и пр. Существующие в настоящее время исследования и разработки в данной области рассматривают в основном непрерывные объекты или процессы, протекающие в непрерывном времени [1-4], при этом аналогичные вопросы, касающиеся такого широкого класса моделей, как дискретные модели стохастических систем, остаются на данный момент практически не освещенными. Кроме того, разработанные методы стохастического оценивания параметров требуют, как правило, для своей удовлетворительной реализации принятия ряда упрощающих ограничений — линейности математических моделей объектов и/или наблюдателей, нормального вида распределения помех, их аддитивности и пр. [5-9]. Ниже предлагается подход к решению задачи стохастического оценивания (идентификации) параметров моделей нелинейных дискретных объектов, позволяющий избавиться от существующих ограничений разработанных методов, а также повысить потенциальную точность процедуры оценивания за счет использования обобщенных (нелинейных) вероятностных критериев.

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

Пусть вектор состояния нелинейного дискретного объекта хк задан в самом общем случае разностным уравнением

хк = Г Ак-^ "к ), (1) где Г (...) — известная нелинейная вектор-функция с компонентами, допускающими обращение; хк-1 — ^-мерный вектор переменных состояния

на (к -1)-м шаге времени; пк — ^-мерный вектор шума с известной ^-мерной плотностью распределения вероятности д(пк); Ак-1 — вектор (или матрица) параметров объекта соответствующей размерности, в общем случае нестационарный.

Наблюдение для к-го момента времени z к размерности М описывается следующим уравнением (в общем случае также нелинейным):

z к = , V к ), (2)

где 8(...) — известная нелинейная вектор-функция наблюдения с компонентами, допускающими обращение; V к — М-мерный вектор шума с известной М-мерной функцией плотности распределения вероятности g^к). Для сокращения дальнейшей записи набор векторов сигналов наблюдения zi (/ = 1,..., к) на текущем интервале времени обозначим через z1k.

В рассматриваемом общем нелинейном стохастическом случае задача стохастического оптимального оценивания текущего неизвестного векторного параметра Лк1 может быть сформулирована как задача нахождения его значения, удовлетворяющего некоторому вероятностному критерию оптимальности J. В качестве таких критериев, обеспечивающих максимальную (потенциально возможную) точность процедуры оценивания, целесообразно использовать нелинейные вероятностные критерии, зависимость которых от апостериорной плотности вероятности (АПВ) переменных состояния р(хк | z1k; Лк_1) в общем случае является нелинейной:

J = { Q (р(х>*; Лк _1)) dxk = W (Лк _1),

(3)

где Q — известная нелинейная аналитическая

функция; X — заданная область пространства состояний.

Здесь важно заметить, что различные вариации вида критериальной функции Q позволяют охватить достаточно широкий класс вероятностных критериев оптимальности [10].

Таким образом, задача в данной постановке сводится к отысканию текущей апостериорной плотности р(хк | z1k; Лк_1) и последующему определению текущего значения искомого векторного параметра объекта Лк1, удовлетворяющего выбранному критерию оптимальности J .

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

Известно [11], что многомерную апостериорную плотность вероятности вектора состояния х для к-го момента времени р(хк | z1k; Лк_1) можно представить в виде

Р(Хк^1к; Лк_1) = Л^х;дЛк )1);

и (Лк_1)

Л(Хк, Лк _1) =

= { ... |™^р(Хк_1 ^1к_1; Л к _2) • q(^(xk, Хк_1; Л^ 1 х, (4)

х dxk _1 •

к, Хк

И\Лк_1) = {...|Л(Хк,Лк,

где р(хк_1 | z1k 1; Лк 2) — определенная на предыдущем (к - 1)-м шаге апостериорная плотность

вероятности вектора состояния объекта; Лк_2 — полученная на (к - 1)-м шаге оценка искомого вектора параметров объекта; 1(хк, хк_1; Лк_1) — N мерный вектор с компонентами, полученными в результате обратного преобразования соответствующих компонентов f (...); Yl — якобиан преобразования от вектора переменных пк к вектору хк; d(zк,хк) — М-мерный вектор с компонентами, полученными в результате обратного преобразования соответствующих компонентов s(...);

— якобиан преобразования от вектора переменных wк к вектору zк .

Тогда с учетом (4) обобщенный вероятностный критерий оптимальности идентификации параметров нелинейной дискретной системы (3) можно окончательно представить следующим образом:

J = |ф

Л(Хк, Лк _1) . Л*(Лк_1)

Л

dx,

В качестве одного из возможных вариантов критерия оптимальности идентификации J используем далее критерий минимума отклонения апостериорной плотности вероятности текущей ошибки оценивания переменных состояния объекта ек от заданной функции на выбранном интервале ее предельно допустимого изменения —

от е„

д° еп

т. е.

тш J = тш [ ... [ (гк(ек) _р(ек 1 ))2de

л , Л,. ,

где ек = хк _ хк — вектор текущей ошибки оценивания, хк — вектор текущих оценок переменных состояния объекта; р(ек | z1k) — многомерная апостериорная плотность вероятности ошибки оценивания; гк (ек) — эталонная функция, выбираемая исходя из физического смысла и особенностей задачи идентификации конкретного объекта.

Учитывая линейную зависимость значений текущей ошибки ек от переменной состояния хк :

А

ек = хк _ Хк, выразим АПВ ошибки оценивания р(ек | z1k) через определенную ранее АПВ переменной состояния объекта р(хк | z1k; Лк_1):

р(ек1 ) = р(ек + хк1 ; Лк _1)

к

е

е

В этом случае процедуру параметрической идентификации на основе минимизации предложенного критерия можно представить следующим образом:

тт J = тт [ ... [ (гк(ек) - Р(ек \ А))2dek =

Ак-1 Ак-1 * ■> ет1п ет1п

етах етъХ { А \ 2

Ь | ... | I гк(ек) -Р(ек + х \ zlk;Ак-1) I dek. (5)

= тт

Ак-1

Производя соответствующую замену переменных в полученном ранее выражении для АПВ параметров состояния (4) и обозначив критериальное выражение через 0(Ак-1), процедуру минимизации критерия (5) можно представить в следующем компактном виде:

тт J = тт 0(А к-1)

(6)

где

П(Ак-1) = | ...] ( Гк(вк) -р(ек + xk\zk; Ак-1) I de

Л(ек + хк, Ак-1) Ь'(Ак -1)

dek.

Здесь важно отметить, что в общем случае решения поставленной задачи вектор текущих оценок переменных состояния объекта хк, входящий в (6), представляет собой некоторый оператор L от многомерной апостериорной плотности распределения переменной состояния, т. е. хк = = L (р(хк \ z1k; Ак-1)) и, следовательно, является нелинейной вектор-функцией от искомого пара-

А

метра Ак-1 : Хк = и (Ак-1 ) .

Тогда критериальное выражение в (6) окончательно можно представить в следующем обобщенном виде:

0(Ак -1) = {...{

Л(вк + и(Ак -1), А к-1) . ^*(А к-1)

Л

dek . (7)

Как было отмечено выше, процедура оптимального оценивания текущего значения параметра, удовлетворяющего предлагаемому критерию, предполагает минимизацию функции (7) на основе известных методов оптимизации функций многих переменных, вариант применения которых рассмотрим ниже в процессе численного моделирования конкретной процедуры оценивания.

ПРИМЕР

Проиллюстрируем эффективность предложенного метода следующим примером с последующим сравнительным анализом полученных результатов с результатами альтернативного подхода к решению поставленной задачи идентификации. Пусть стохастический дискретный объект задан нелинейным разностным уравнением

хк = /(хк -1, ак-1)+пк = ак-1 • хк-31 + пк, х =1:

(8)

где ак-1 — неизвестный искомый параметр (для рассматриваемого далее модельного примера исходное значение этого параметра ак-1 = 1 для всех к (к = 1,..., ^набл, ^набл — общее количество доступных измерений); пк — шум объекта с плотно. л 1 ( дь

стью вероятности д(п) = — -

п ^ (п - да)2 + дЬ2

(плотностью вероятности Коши) с параметрами да = 0, дЬ = 0.03 .

Наблюдение переменных состояния заданного объекта осуществляется дискретным измерителем, описываемым уравнением

2к = -К хк) + ^ = ск • 45 + ^,

где ск — известный параметр наблюдателя (для рассматриваемого далее модельного примера значение параметра ск = 1 для всех к,); ^к — шум измерителя с плотностью вероятности g(м>) =

1 ( gЬ ^

п

(п - ga)2 + gЬ2

(плотностью вероятности

Коши) с параметрами ga = 0 , gЬ = 0.04 .

Соответственно функция 0,(ак-1) для к-го шага алгоритма примет вид

^(ак-1) = ](Гк (ек) - Р(ек + хк \ ^; ак-1)

de =

е

е

к-1

к-1

е

е

е

е

е

е

2

е

гк (ек) _

{р(Хк_1 | г^1;ак_2) • q^ 1(ек + Хк,Хк_1;ак_1)^^ • g^й(гк,ек + Хк)^

И'(ак _1)

\

q ^1 (ек + хк, хк _1; ак _1) ^ = q (1 (ек + и (ак _1), Хк _1; ак _1)) = q (ек+и (ак _1) _ ак _1 • Хк_31):

g Iй (гк, ек + Хк ) I = g (й (гк, ек + и(ак_1 )) ) = g ( гк _ С • (ек + и(ак_1 ))1 5 ) .

de,

(9)

В качестве эталонных функций гк (ек), к = 1,..., ^абл выберем численно заданные на интервале [ етш =_1, етах = 1] зависимости (рис. 1), удовлетворяющие очевидным требованиям к виду АПВ текущей ошибки оценивания для случая корректного (близкого к истинному) определения значения искомого параметра, а именно требованию положительной определенности, требованию расположения максимальных значений этой функции вблизи нулевого значения ошибки и одновременно достаточно небольшой величины ее интегрального значения на выбранном интервале. Отметим также, что предлагаемый вариант критерия

позволяет дополнительно учитывать динамические особенности используемого метода получения оценок переменных состояния (выбор которого для данного примера обсуждается ниже), например динамику абсолютной и/или относительной величины получаемых ошибок оценивания на различных временных интервалах.

Как показали приведенные ниже результаты численных экспериментов, подобная функциональная зависимость гк (ек) достаточно эффективно обеспечивает необходимую точность и скорость сходимости алгоритма идентификации для исследуемого объекта.

е

1

¡ежей

ШВШВююшИввШ^^

1.0

-0.5

Рис. 1. График зависимости эталонных функции гк (ек) от номера к шага алгоритма идентификации.

ек — ошибка оценивания переменных состояния объекта

Рис. 2. Зависимость АПВ текущей ошибки оценивания V (ак-1, ек) от идентифицируемого параметра ак -1 для к-го шага алгоритма (к = 22)

Априорная плотность вероятности для первой итерации алгоритма задавалась плотностью вероятности Коши с параметром сдвига 1.3 и параметром масштаба 0.2.

Определение оценок было произведено с использованием дискретного нелинейного алгоритма калмановской фильтрации [11] при следующем выборе параметров алгоритма фильтрации: дисперсии шумов объекта и наблюдателя принимались равными Dn = 0.01 и Dw = 0.05 соответст-

А

венно; начальное значение оценки — Х1 = 0.85; начальное значение элемента корреляционной матрицы ошибок — R1 = 1.

На рис. 2 представлен полученный в результате моделирования график входящей в формулу (9) АПВ текущей ошибки оценивания для к-го шага алгоритма (к = 22), которая, являясь функцией текущей ошибки ек , зависит в то же время от значе-

А

ний искомого параметра ак-1 (р(ек + Хк \ г1к) =

= р(ек + и (ак-1)Кк) = V (ак-1, ек)).

Определение интегралов, входящих в выражение (9), здесь и в дальнейшем производилось численно с использованием квадратурных формул с шагом А = 0.05 . Бесконечные пределы интегрирования по переменной состояния Х были заменены на конечные значения, удовлетворяющие точностным требованиям к алгоритму оценивания

( Хтп = ^ Хтах = 4).

В качестве одной из интересных особенностей графика рис. 2 можно отметить сравнительную простоту распознаваемости диапазона наиболее вероятных значений текущих оценок искомого параметра объекта, удовлетворяющих выбранному критерию, — за счет резкого изменения формы функции V (ак-1, ек) вдоль оси оцениваемого параметра ак-1 . При этом сечение построенной функции вдоль оси ек (на выбранном интервале [ етт, етах ]), наименее отличающееся от эталонной функции гк (ек) (а следовательно, и минимум критериального выражения (9)), располагается вблизи истинного значения искомого параметра (ак-1 = 1).

Рис. 3. Зависимость функции критериального выражения 0,(ак-1) от искомого параметра ак-1 для различных шагов к алгоритма

Рис. 4. График зависимости функции критериального выражения 0.(ак_1) от номера к шага работы алгоритма идентификации

Наглядным подтверждением этого является приведенный график самой функции 0(ак-1) на этом же шаге алгоритма (см. рис. 3).

Как показали результаты моделирования, вид приведенных на рис. 3 зависимостей является характерным для критериальных выражений, получаемых на различных итерациях алгоритма. Здесь важно отметить, что, имея многоэкстремальный характер, критериальные выражения на различных шагах алгоритма (при к > 17) принимают свои наименьшие значения в районе истинного значения искомого параметра ак -1 = 1.0 (рис. 4).

Для минимизации определенной на очередном шаге алгоритма целевой функции 0(ак-1) использовался метод численной минимизации, основанный на прямой подстановке набора значений искомого параметра (заданных в интервале его возможных значений с шагом 0.05) и обеспечивающий в рассматриваемом случае достаточную оперативность и удовлетворительную точность получаемых результатов.

Результаты моделирования процедуры нелинейной параметрической идентификации показали, что при выборе количества дискретных значений сигнала наблюдения к > 17 отклонение оценки

Рис. 5. Зависимости значений ак х

оценки искомого параметра от количества наблюдений для двух рассмотренных подходов

параметра объекта ак -1 от его истинного значения ак-1 = 1 не превышает 5-10 % от его величины (рис. 5).

Для анализа сравнительной эффективности предложенного подхода выберем в качестве альтернативного подхода к идентификации параметров заданного объекта (8) один из наиболее известных и широко используемых методов параметрического оценивания — метод наименьших квадратов (МНК) и сопоставим полученные результаты рассмотренных методов.

Критерий МНК, подлежащий минимизации по переменным состояния объекта и искомому параметру, в рассматриваемом случае формулируется следующим образом [13]:

е=еГ( хк -/ (хк _l, а) )2+(^-хк) )2

Соответственно, расширенный вектор оценок переменных состояния и искомого параметра А, определяемый по методу наименьших квадратов, определится как

А =

= а^тт( Е).

А

(10)

Процедура минимизации в (10) производилась на основе одного из достаточно распространенных методов прямой оптимизации — метода Нелде-ра—Мида, показавшего наилучшие результаты из всех исследованных минимизационных процедур для данного примера. Указанный метод был реализован с использованием стандартного набора методов оптимизации среды математических вычислений МаШетайса.

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

Сравнительный анализ приведенных зависимостей показывает, что, несмотря на менее уверенную работу предлагаемого алгоритма в начальный интервал времени (при к < 17), значения оценок искомого параметра, получаемые в дальнейшем, оказываются значительно более близкими к его истинной величине, чем при использовании альтернативного подхода на основе МНК. Кроме того, приведенный график наглядно демонстрирует качественное расхождение в тенденциях изменения получаемых оценок параметров с увеличением времени работы алгоритма идентификации.

Представленные результаты исследований подтверждают принципиальную возможность эффективной реализации предложенного метода параметрической идентификации с использованием обобщенных вероятностных критериев, и в частности с использованием критерия минимума отклонения АПВ ошибки оценивания от заданной функции.

2

ЗАКЛЮЧЕНИЕ

Таким образом, в работе получено решение задачи стохастического оптимального оценивания параметров нелинейных дискретных объектов, обладающее рядом новых свойств:

— более высоким по сравнению с традиционными методами уровнем потенциальной точности процесса идентификации параметров за счет использования обобщенных вероятностных критериев, зависящих в общем случае нелинейно от АПВ распределения вектора состояния и позволяющих охватить широкий класс условий оптимальности;

— возможностью его использования при различных видах плотности распределения вероятности шума как объекта, так и наблюдателя;

— отсутствием ограничений на использование предлагаемого подхода также и для нестационарных неизвестных параметров стохастических объектов;

— принципиальной возможностью применения метода для нелинейных дискретных объектов и наблюдателей, в том числе при нелинейной зависимости функции объекта от искомого параметра.

Работа выполнена при поддержке РФФИ (грант № 10-07-00158).

СПИСОК ЛИТЕРАТУРЫ

1. Каргин А.В., Фатуев В.А. Об одном методе структурно-параметрической идентификации динамических систем // Автоматика и телемеханика. 2006. № 4. С. 116-125.

2. Мамай В.И., Сотников В.И., Щербань О.Г. Субоптимальная параметрическая идентификация нели-

нейных динамических систем // Известия вузов. Радиоэлектроника. 2005. № 3. С. 15-23.

3. Narayanan M., Narayanan S. Parametric identification of nonlinear systems using multiple trials // Nonlinear Dynamics. 2007. N 4. P. 341-360.

4. Blasco X., Herrero J.M., Martinez M. Nonlinear parametric model identification with Genetic Algorithms. Application to a Thermal Process // Lecture Notes in Computer Science. 2001. V. 2084. P. 466-512.

5. Льюнг Л. Идентификация систем. Теория для пользователя. М.: Наука, 1991. 432 с.

6. Гроп Д. Методы идентификации систем. М.: Мир, 1979. 302 с.

7. Справочник по теории автоматического управления / Под ред. А.А. Красовского. М.: Наука, 1987. 712 с.

8. Штейнберг Ш.Е. Идентификация в системах управления. М.: Энергоатомиздат, 1987. 87 с.

9. Сейдж Э., Мелса Дж. Идентификация систем управления. М.: Мир, 1974. 248 с.

10. Хуторцев В.В., Соколов С.В., Шевчук П.С. Современные принципы управления и фильтрации в стохастических системах. М.: Радио и связь, 2001. 808 с.

11. Тихонов В.И., Харисов В.Н. Статистический анализ и синтез радиотехнических устройств и систем. М.: Радио и связь, 1991. 608 c.

12. Левин Б.Р. Теоретические основы статистической радиотехники. М.: Радио и связь, 1989. 656 с.

13. Эйкхофф П. Основы идентификации систем управления. М.: Мир, 1975. 680 с.

Ростовский государственный университет путей сообщения, г. Ростов-на-Дону

Контакты: Кучеренко Павел Александрович, [email protected]

Материал поступил в редакцию 24.10.2010.

THE STOCHASTIC ESTIMATION PROCEDURE OF THE PARAMETERS OF THE NONLINEAR DISCRETE OBJECTS BASED ON THE GENERALIZED PROBABALISTIC CRITERIONS

P. A. Kucherenko, S. V. Sokolov

Rostov State University of Transport Communication, Rostov-on-Don

The relevance of researching of the alternative (against to the traditional) methods of the stochastic estimation of parameters for the nonlinear discrete objects is shown. For the first time the algorithm for the optimal estimation of the parameters for the mathematical models of the discrete nonlinear stochastic objects on the ground of the minimum probability criterion of the estimation error for the proposed algorithm is suggested. The model example, illustrating the efficiency for realization of the suggested algorithm is considered.

Keywords: nonlinear discrete plant, stochastic parameters estimation, generalized probabilistic criteria, minimum criterion of the estimation error probability