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

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И КОМПЬЮТЕРНЫЕ ТЕХНОЛОГИИ

УДК 519.248:681.5.001.3

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

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

Ростовский государственный университет Rostov State University of Transport

путей сообщения Communications

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

Ключевые слова: стохастический дискретный объект; нелинейные вероятностные критерии; оптимальное оценивание параметров.

The solution of the problem of the stochastic discrete plant's parameters optimal estimation problem based on the nonlinear probabilistic criterions is obtained in the paper. The optimal estimation algorithm of the unknown parameters of the discrete plant is offered using the criterion of the probability minimum of the error of estimation. The efficiency of the offered approach is demonstrated.

Keywords: stochastic discrete plant; nonlinear probabilistic criterions; optimal estimation of the parameters.

Задачи оптимального оценивания параметров стохастических систем составляют один из наиболее важных разделов все более расширяющейся в настоящее время проблематики, связанной с вопросами построения математических моделей и прогнозирования процессов самой различной природы. В подавляющем большинстве работ, посвященных решению существующих проблем данной области, в качестве объекта исследования в основном выступают системы или процессы, описывающиеся моделями с непрерывным временем [1—4]. Вместе с тем аналогичные проблемы, актуальные для дискретных стохастических моделей, остаются на сегодняшний день практически не освещенными.

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

тодов [5—9] (необходимости линейности используемых математических моделей, требования нормального вида распределения помех и т.д.).

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

хк = f (хк-^ Ак-1) + пк, (1)

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

Вектор наблюдения для к-го момента времени z к размерности М описывается следующим уравнением (в общем случае также нелинейным): zк = 8(хк ) + wк , (2)

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

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

I = | О (р(Хк | гк; Ак_1)) ) = №(Ак_1), (3)

X

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

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

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

Известно [11], что многомерная апостериорная плотность вероятности вектора состояния х для к-го момента времени р(хк | гк; Ак-1) определяется выражением

р(хк 14;А к _1) =

J ••• J p(xk-1 | zf-1; Ak-2) • p(x h( A f-1)

x A f -i)dx f -1 • p(zf 1 xf);

-x (4)

h(Ak_i) = J... J J... Jp(xk-i | zf-1; Ak-2)x

—^ —^ —^ —^

V

x p(xf | Xf-1; Af-1)dXf-1 • p(Zf | Xf) |dxt

где р(х к _1 | г!к-1; А к-2) — определенная на предыдущем — (к-1)-м, шаге апостериорная плотность вероятности вектора состояния объекта; А к _2 — полученная на (к—1)-м шаге оценка искомого вектора параметров объекта; р(хк | хк_!; Ак_!) — определяемая на текущем шаге алгоритма ^-мерная условная плотность вероятности вектора хк ; р(г к | хк) — определяемая на текущем шаге алгоритма функция правдоподобия для многомерного наблюдения.

Многомерная условная плотность р(х к | х к _ь А к _х) может быть получена из исходного уравнения объекта (1) при известном виде плотности распределения вероятности значений шума #(п) (в предположении их взаимной статистической независимости), как [12]:

p(xf I Xf-1; Af-1) = q(l(xf, Xf-1; Af-1)) • Y

4;

l(x f,x f-1; A f-1) =

xf(1) - f(1)(xf-1,A f-1) xf (N) - f(N)(x f-1,A f-1)

где Y/ — якобиан преобразования от вектора переменных п к к вектору хк (для частного случая аддитивного шума объекта Y/ = 1); %)(хк-1), I = 1... N — компоненты вектор-функции f(хк-1), входящей в (1); хщ), I = 1.. N — компоненты текущего вектора состояния объекта хк .

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

p(zf I Xf) = g(d(zf, xf)) • Y

d ;

d(zf, xf ) =

zf (1) s(1)(x f ) M

zf (m) - s(M)(xf)

где Yd — якобиан преобразования от вектора переменных wк к вектору г к (для частного случая аддитивного шума наблюдения Yd = 1); S(j )(х к), ] = 1.. М — компоненты вектор-функции s(x к), входящей в (2); 1к (j), j = 1.. М — компоненты вектора наблюдение г к .

Так как АПВ р(х к _1 11\ -1; А к-2) в правой части равенства (4) на к-м шаге является известной функцией (определенной на предыдущем шаге), рекуррентный алгоритм определения АПВ переменной состояния для к-го момента времени при наличии текущей последовательности дискретных отсчетов сигнала наблюдения гк принимает вид:

p(xk | zk; Ak-!> = A!Xk;Ak--l) h (Ak-i)

Л(х k,A k-i) =

(5)

= j" ••• j" I P(Xk-i | zk-1; Ak-2) • q(l(xk, xk-i; Akл)) x dxk-i • g (d(zk, xk));

*( a k-i) = j ••• Jmx k ,A k-i)dx k •

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

J = "ф

X

^(x k »A k-i) 4 h*(A k-i)

dx ь

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

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

'-шах '-тах

min J = min j••• jp(ek | zk)dek

'-шт '-шт

где ек = хк _ хк — вектор текущей ошибки оце-

Л

нивания; х к —вектор текущих оценок переменных состояния объекта; р(ек | хк) — многомерная апостериорная плотность вероятности ошибки оценивания.

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

Л

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

Р(ек I 4) = Р(ек + хк | 4;Ак_1) •

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

сшах сшах

min J = min j••• jp(ek | zk)dek =

'■шах '■шах

min j••• j p(ek + xk | zk; Ak-i)deA

(6)

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

min J = min Q(A k-i)

A t-i A k-,

(7)

где

шах шах

''шт ''шт

0(Ak-i) = j••• j p(ek + xk | zk; Ak-i)dek

Clin

А

Л(ek + x k, A k-i)

'■max '■max

сшт сшт

h (Ak-i)

de

k •

Здесь важно отметить, что в общем случае решения поставленной задачи вектор текущих

Л

оценок переменных состояния объекта хк , входящий в (7), представляет собой некоторый оператор X от многомерной апостериорной плотности распределения переменной состояния, т. е.

хк = х((хк | ; Ак-1)), и, следовательно, является нелинейной вектор-функцией от искомого

Л

параметра А к _1: х к = и(А к .

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

''шт ''шах

Q(A k-i) = j ••• j

'■шт '■шт

^Л(е k + U(A k-i), A k-i)Л h*(A k-i)

de

k • (8)

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

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

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

Хк = /(хк_ь ак_1) + пк = ак_1 ■ Хк_1 + пк , Х1 = ^

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

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

1к = ) + ^к = Ск ■ Хк + ^к ,

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

е

ьтах Л

^(ак_1) = | р(ек + Хк | г!; ак=

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

Для определения текущего значения оценки переменной состояния объекта использовался оптимальный фильтр Калмана [11], который в рассматриваемом случае имеет вид:

xf = U (af -1) = af -1 • xf -1 + Kf

/ f А ^ \

Zf - С • af -1 • xf-1

V /

А

x1

0,75;

Kf = CfD

Dw

Rf =

-1

lf-1 • Rf -1 + Dn

£f R.,

R1 = 1 .

С f 1 А А А

J P(xf-1 1 Zi- ; af-2) • l(ef + xf, xf-1; af^x^ • d(Zf, ef + xf)

h (af-1)

l (ek + xf, xf-1; af-1) = l (ef + U (af -1), Xf -1; af -1)

( +U (af-1)-af -1 •xf-1 )2 = 1 e 2Rn

^ Rn

А

d (Zf, Qf + xf) = d (Zf, ef + U (af-1)) =

( -c^ +U(af-1)))2

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

af-1 (P(ef + xf | Zf) = p(ef + U(a) | Zf) = V(af-1,ef)).

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

•^mm v -^max ^

dek ; (9)

2Dw

R

Априорную плотность вероятности для первой итерации алгоритма выберем нормальной с дисперсией В0 = 0,8 и математическим ожиданием 1,5. При этом интересно отметить, что в данном случае отклонения среднего значения априорной плотности от начального значения переменной со-

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

А

А

+

А

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

( етш. _ —0,2 , етах _ 0,6 )-

-

—W-

1.2

1.0 1.5 2.0

Рис. 2. Зависимость функции критериального выражения от оцениваемого параметра для к-го шага алгоритма (к=70)

Как показали результаты моделирования, вид приведенной на рис. 2 зависимости является характерным для критериальных выражений, получаемых на различных итерациях алгоритма. Здесь важно отметить, что свои наименьшие значения на различных шагах алгоритма (при к>65) критериальные выражения принимают в районе истинного значения искомого параметра ак—1 = 1,0 .

С учетом особенностей решаемой задачи оптимального оценивания параметров, выявленных в процессе проведении численных экспериментов, для минимизации выбранного критерия на очередном шаге алгоритма — т. е. для однозначного определения численного значения текущей оценки искомого параметра, предварительно задавался некоторый интервал возможных значений параметра ак—1 (в примере был выбран интервал

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

Результаты компьютерного моделирования процедуры нелинейного оценивания параметров показали, что при выборе количества дискретных значений сигнала наблюдения к > 69 отклонение оценки параметра объекта ак—1 от его истинного значения ак—1 = 1 не превышает 10 % от его величины. Наглядной демонстрацией этого является приведенная на рис. 3 зависимость значений оценок искомого параметра от номера шага алгоритма оценивания.

10 20 30 40 50 60 70 £

Рис. 3. Зависимость значений оценки искомого

параметра от номера шага работы алгоритма оценивания

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

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

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

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

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

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

Литература

1. Каргин A.B., Фатуев В.А. Об одном методе структурно-параметрической идентификации динамических систем // Автоматика и телемеханика. 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. №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. Vol. 2084. P. 466-512.

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

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

7. Красовский А.А. Справочник по теории автоматического управления. М., 1987.

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

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

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

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

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

13 сентября 2010 г.

Соколов Сергей Викторович — д-р техн. наук, профессор, Ростовский государственный университет путей сообщения. Тел. (863) 272-63-02.

Кучеренко Павел Александрович — канд. техн. наук, доцент, Ростовский государственный университет путей сообщения. Тел. (863) 272-63-02. E-mail: pavelpost83@mail. ru

Sokolov Sergey Victorovich — Doctor of Technical Sciences, professor, Rostov State University of Transport Communications. Tel. (863) 272-63-02.

Kucherenko Pavel Aleksandrovich — Cfndidate of Technical Sciences, assistant professor, Rostov State University of Transport Communications. Tel. (863) 272-63-02. E-mail: pavelpost83@mail. ru