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

Доклады БГУИР

2011 № 7 (61)

УДК 621.376

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

А.В. ОВСЯННИКОВ

Белорусский государственный университет пр. Независимости, 4, Минск, 220030, Беларусь

Поступила в редакцию 10 октября 2011

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

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

Введение

Метод максимального правдоподобия (ММП) в связи с его применением к оценке параметров известен давно (Р. Фишер, Г. Крамер и др.). Получающиеся на основе ММП алгоритмы обработки обладают простотой, приемлемым качеством, не требуют существенных априорных сведений об информационных параметрах сигналов. В то же время, использование только белого гауссовского шума при построении алгоритмов обработки на основе соответствующего функционала (К. Ито, Р.Л. Стратонович) сужает применение метода. Подход с использованием функции потерь (Я.З. Цыпкин и др.) и соответственно выборочный эмпирический функционал (ВЭФ) позволяют учитывать разнообразие возможных распределений шумов, строить робастные и адаптивные алгоритмы, однако выбор функции потерь (или "наихудшего" в классе распределения) во многом остается субъективным.

Применение функционалов марковских процессов для реализации фильтров максимального правдоподобия (далее просто фильтр) [1-3] хотя и является формально ограниченным, в то же время позволяет существенно выходить за рамки марковской теории [4]. Так в более общих случаях имеет смысл рассмотрение фильтра недельтакоррелированного шума (или класса шумов с экспоненциальной огибающей корреляционной функции), заданного СДУ [2]. В сочетании с наличием в СДУ флюктуирующих параметров и малых отношений сигнал-помеха (ОСП) на входе, синтез такого фильтра представляется интересным с теоретической и практической точки зрения.

Цель статьи - получить обобщенный функционал и реализации фильтров максимального правдоподобия, отличающихся от классических оптимальных фильтров для СДУ, описывающих медленные флюктуации интенсивности помеховой обстановки.

Теоретические положения

Алгоритм ММП реализуется соотношением Vх ln P (r | X) = 0, (1)

где Р(г | X) = Р(г - s) = Р(Ъ) - функция правдоподобия, г^) = s(t,X) + Ъ(0, Vх = (й / йХ) -оператор, X- вектор информационных параметров сигнала s(t,X), Ъ(0 - марковский процесс, в общем случае с произвольной корреляционной функцией. В частности, модуль огибающей корреляционной функции может иметь экспоненциальный характер [2]. Переходя к многомерной функции правдоподобия в (1) введем обозначения: Ъ = [Ъ^—,Ъп], г = [1,...,гп], s = [к^...,кп],

Ъ = г - , , = 1, п . При постоянном коэффициенте диффузии Ь=сош^ связь СДУ марковского процесса Ъ(0 с его одномерной плотностью вероятности (далее плотностью) Р(Ъ) имеет вид:

Ь(Ъ) = Ъ (0 + / (Ъ) = тп(/), / (Ъ) = bZ (Ъ) / 2, (2)

где обозначено Ъ (t) = )/ й, 2 (Ъ) = -й 1п Р(Ъ)/ й Ъ - нелинейное преобразование, осуществляемое над процессом Ъ(0 с плотностью [4] Р(Ъ) = (С /4ъ )ехр(-2 Г, / (х)йх / Ь), у = ^2ЬГы , N

- односторонняя спектральная плотность белого шума п(0 с нулевым средним, С - константа, Ь - обеляющий оператор. Плотность перехода марковского процесса (2) имеет вид (Д = t¡+1 -1¡

малый интервал времени, ) = Ь = / 2)

Р(Ъ+1 I Ъ,) = ехр |- -1

Ъ — Ъ

+ / (Ъ,)

Д

Д I ^2%ЬД

Методика синтеза функционалов и фильтров максимального правдоподобия

Рассмотрим функцию 1п^п 1Р(Ъ,+1 I Ъ,). Выделим сигнальную часть и введем обозна-

чение

1 п-1

*(Ъ)=-2ь §

Ъ,+1 -Ъ,

Д

+ / (Ъ,)

Д

(3)

Для того чтобы в дальнейшем оперировать симметризированным стохастическим интегралом, проведем дифференцирование (3) по информационному вектору X = [X i, Xi+1] в точке

X,,,+1 = (X* + X,0+1)/ 2 , где X- вектор экстраполированных на (, +1)-ый шаг оценок. Получим систему уравнений

1 п-1 -1 §§

Ь £

п-1 -1 §

ь £

Ъ ,+1,,+1 Ъ /',/'+1 Д

Ъ/ + 1,/+1 Ъ,,г

+ Г (Ъ,,,+1)

V xS , (XiMl)

Д

+ / (Ъ,,,+1)

Д

V xSi+1( X ,,,+1) Д

+ V XГ (Ъ,,,+1)

Д = о,

(4)

Д = 0;

где обозначено Ъj,,+1 = (Ъ0 +Ъ*)/2, Ъ0 = Г -к,-(XI0+l), Ъ* = Г -Ку(X*), у =, ,, +1. Возможность представления (4) определяется ошибкой, возникающей за счет отбрасывания последующих членов разложения в ряд (X,,,+1) в точках [X *, XI'+1]. В частности, необходимым условием справедливости записи системы (4) являются следующие:

'К (X,*)| > (X,0+1 - X*)т 8ир(| VX VXs, (X;)|)(Xi0+l - X*)/ 2, К+ДX,0+1) > (Xi0+l -X,*)т 8ир(|VXVXs,+l(X,°l)|)(X,°l -X*)/2.

Заметим, что необходимое условие не является достаточным, поскольку не учитывает влияния следующих членов аппроксимации.

Просуммируем полученные выражения системы (4):

2

1 n-1

1 z btt

— i+1,i+1 — i, А

■ + f C^i-i+t)

V xsr+t( Xi,i+t) -V xsr (Xi,i+t) А

-V xf (-■,,+1)

А = 0.

Поскольку входящие в последнее уравнение переменные оказываются симметризированными, переходя к пределу при А ^ 0 получим:

1 т 1 т

1 (X \ т) = V ^ (X *,т) = - -1 ¿©V =1V X* + ff¡'V х^ - V X = 0. (5)

Здесь и далее, там где это не вызывает сомнений, зависимости функций от аргументов опускаются.

Описанная выше методика симметризации позволяет использовать ее для синтеза функционала Р(£) непосредственно из (3). Такой же результат можно получить, формально интегрируя уравнение (5) по области изменения Х:

1 T

P(-) = hT exp{-F(,,T)} = hT exp{-—{L2(,)dt} .

2b i

(6)

где ^ - нормирующая величина, не зависящая от ^ .

Представляет интерес внутренняя структура функционала (6), следующая из (5)

P(,) = Иа exp

J_ 2b

> + ,

lo Jü

Jo'f=

Tt

P, (0) - у P, (0)

(7)

где hА - нормирующая величина. При выводе (7) учтено, что Т = 1/ М (

f-

мени процесса (2), ст? = b / 2M (f-

1 гт ■■

р.. (0) = ^ 2 ^ - нормированные производные автокорреляционной функции (АКФ).

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

) - постоянная вре-1 ст„

1 С '

) - стационарная дисперсия, р, (0) = JTTJ0 -—dt и

V(T) = exp(-T / T,) / T,. Таким образом, полагая И = lim ИА exp--(р, —- р,) , функционал (7)

T —м А—0

T

2

2

принимает вид P(,) = И exp I--

2b

'I + ,0

ff2 dt

, где составляющая , - характеризует дина-

мическую площадь процесса ,(t).

В случае однозначного преобразования — = f^l(ini) имеем P(—) = P(f)|df /d . Полагая, что n = 1 n(t )dt / А и M (n2) = N /2А, запишем одномерную плотность

JA

P(—) = А exp(- f2A /2b)/>/2rcbA . С учетом Pf (,) = lim Пn_ P(—) получаем функционал слу-

f n—=1

чайного процесса, формируемого нелинейным уравнением f (,) = yn(t):

Pf (,) = Ит exp^2b[f2 (,)dt),

(8)

где обозначено кт = ЛА ехр(т /2т - т / А), т = пА . При т ^ да и т >> А функционалы (7), (8) совпадают. Время достижения стационарного состояния определяется постоянной времени системы (2) т .

Ограничившись в качестве достаточных статистик первыми двумя моментами аппроксимации (6), качество оценки X (Т) фильтром (5) определится величиной

KT ={vXF(X*,T)}- = ¿|£Т(vxZ(^*))2 + L(^*)VXL(^)Jdt} \

dt} . (9)

В табл. 1 приведены общие уравнения фильтров для трех видов сигналов. Там же, введены обозначения 5*(0 = г (?) - s(t, Х*(?)), f '(5*) = df (5)/d 5|5=5* .В уравнениях фильтров следует учесть, что подынтегральные составляющие, содержащие производную по времени, определяют динамическую разность между граничными значениями на концах интервала интегрирования. Так, общее уравнение фильтра (№1), (№2) принимает вид

b Л [/(О + Ж)/(Г)] dt = 1 (лЖ) + £7(Г)/'£>) = 0, (10)

где А/(5*) = f (5*(Т)) - f (5*(0)). Переходя к функции нелинейного преобразования из (10) получим (А2 (5*) = 2 (5*(Т)) - 2 (5*(0))):

ЛZ + Ь?т

j0Z (Г) Z ^')dt = 0. (11)

2 4J°

Качество фильтра в соответствии с (9)

T

Kt =-ь{Л/'0Т) + j0T [(/'(О)2 + /(Г)/"(О (12)

Таблица 1. Уравнения фильтров максимального правдоподобия

№ Вид сигнала Уравнение фильтра

1 s(t, X) = X = const 1 Т -Г[ г(?) + / (5*) ]/ '(5*)^ = 0 Ь 0

2 s(t, X) = X (t) 1 Т -{[5 *(?)+/(5*) ]/ '(5*> = 0 Ь 0

3 s (t, X) = s0(t) X (t) 1 Т -Г[5 *(?)+/ (5*)][ ^) + / '(5*) ^) ] Л = 0 Ь 0

Обобщенная структура фильтра приведена на рис. 1. Отличие от известных схем демодуляторов [5] - наличие блока фиксации динамического диапазона А/(5 ) и блока /'(5 ) .

Рис. 1. Структурная схема фильтра (10)

Пример. Рассмотрим работу фильтра (11) при изменении структуры шума. В частности, шум может иметь характер процесса с медленными нестационарными изменениями. В этом случае модель шума можно описать одномерным распределением [6]:

р =аГ((1 + к)/2)^*к = р/а . 4ьГ (к / 2)

(13)

Такая модель шума позволяет учесть плавные изменения интенсивности и амплитуды процесса Ъ(0. При параметре а ^ 0 распределение стремиться к гауссовскому, а при а ^ да к лапласовскому (см. рис. 2). Нелинейное преобразование в случае (13) I(Ъ) = рй(аЪ).

0.5 0.4 0.3

т 0.2

0.1

о

-10

Л а=10

7

ГЧ <1-11.5

а=0„1

/

10

О-в 0.6

0.4

0.2 Т(|) о -0.2 4.4 -0.6

А уи

г < 1......л

1к 01=0.5 /................

^.5

2.6

Рис. 2. Распределение Р(Ъ) (13)

Рис. 3. Характеристика фильтра (14)

Пусть сигнал имеет форму №2 табл. 1. Уравнение фильтра (11) принимает вид А1 (Ъ*)/2 + Ь^ /4 = 0, (14)

где ¥(Ъ) = ар2?Л(аЪ)*А2(аЪ) - характеристика фильтра. Вид формы характеристики приведен на рис. 3 (р = 1). Такая форма нелинейной характеристики не фиксирует уровень ошибки Ъ (/) оценки X (/) (жесткий или мягкий ограничитель), а адаптивно уменьшает его [7]. Текущая

оценка качества фильтра, определяемая из (12), имеет вид Кт =-2 где ¥ '(Ъ) = -а2р2(сЛ(2аЪ) - 2)эес^(аЪ).

|а! ' (Ъ*) + V '(Ъ)Л/2}\

20 15 10

Косп 5 0 -5 -10

л

/2 •3

У / . 5

4

0.5

1.5

2.5

X 105

Рис. 4. Коэффициент амплитудного подавления Рис 5 Выигрыш в ОСП фильтра (14) в сравнении с шума с распределением (13) построенным на основе ВЭФ

Дискретную реализацию фильтра можно получить применяя к (5) стохастическую градиентную пр°цедуру ^ = К -(Vх 1п(1 С))) \ где Ък = rk-s(х¡), Ък+1 = Гк+1 -s(X'k) -экстраполированное значение. Конкретизация такой процедуры для фильтров №1 и №2 табл. 1 имеет вид Xk+l = Xk-Кк+1{2(ЪЭ+1)-1(Ъ0) + Ь1 (Ък)1 '(Ък)/2), К+1 = К- +1'(ЪЭ+1>/2.

Эффективность фильтров максимального правдоподобия

Нелинейное преобразование текущего значения процесса 2(5, t) = -£(5)У.£(.)/ Ъ изменяет ОСП на его выходе (дех ■ 1) де^1х = ц., где Ц. = а^М(22_) - матрица коэффициентов амплитудного подавления шума. Элементы матрицы 1. = М(22) (информационная матрица Фишера) на основании определенных ранее плотностей Р(..+1 1.) и Р(..) для /(5,.) = а5,-

имеют

следующий

вид

151+1,1+1 = (ЪА)-1, 15,+1,, =-(ЪА)-1(1 -аА))

51+1,1

151,1+1 151+1,1,

I ^ = (ЪА) 1(1 - аА))2. Вычисляя величину дисперсии а2; = Ъ / Аа2, получаем базовую матрицу коэффициентов подавления (где обозначено г = ехр(-аА), аА ■ 1):

Ц.

(1 - г Г2 -г (1 - г )-2 -г(1 - г)-2 -г 2(1 - г У2

(15)

Ц

.0

51+1,1+1 2

Для количественной оценки подавления достаточно рассматривать только величину = (1 - г)2. Такую же оценку можно получить для конкретной функции / (.■), вычисляя

ст^ и 15. Для ориентировочной оценки ц. с произвольной /(. ■) можем принять

Ц5 ~ Ц/ Ц0 .

(16)

А), а

В данном случае коэффициент, входящий в (16), равен величине г = ехр(-М(/0'

Ц/ =оМ (252) - коэффициент амплитудного подавления, численные значения которого для различных функций / (5) легко вычисляются.

Рассмотрим процесс с одномерной плотностью (13). Информация Фишера в этом случае 15 = ар2 /(а + Р) и с учетом определенной в [7] стационарной дисперсии (а,Р) зависимость коэффициента Ц/ от параметра k = Р / а имеет вид (см. рис. 4):

2к+2 Г

Ц /(к) = -

к +1

к (1 + к Г

4 Ъ

к к к к 2, 2, 2,

ккк — +1,— +1,— +1 222

-1

Потери эффективности фильтра с характеристикой 2 (5, ?), оценивают реальным коэффициентом подавления [5] ц5р=а2 {м. (2 )| /М. (2 .). Используя неравенство Коши-

Буняковского-Шварца и, ограничившись первым членом разложения /' в ряд Маклорена, получим оценку верхней границы коэффициента подавления ц5р <Ц.М.(/.2)/М.(£2). Здесь под £ = ) понимается процесс, в общем случае отличный от дельтокоррелированного. Например, для линейного фильтра с характеристикой /(5 ) = а. можно получить точные значения ц. . Так, в случае действия на входе фильтра процесса с нормированной корреляционной функцией рс = ехр(-р|т|) оценка ц= 1/(1 + Р/а), а для р^ = зт(А^т)/А^х приходим к оценке

цр = агС§(О) / О, О = А^ / а .

Моделирование

Результаты моделирования фильтра (14) в сравнении с обычным фильтром максимального правдоподобия, построенным в соответствии с ВЭФ [8] (1) представлены на рис.5. Выиг-

рыш в ОСП оценивался коэффициентом Косп = 10lg(Р0 / Р1), где Р1, Р0 - оценки мощности шума на выходе указанных выше фильтров. Отмеченным экспериментальным кривым соответствуют значения параметров плотности (13), приведенные в табл. 2 (р = 1, РХвх - -9дБ - мощность

сигнала). Полученный результат поясняется наличием в фильтре (10) величины /'(Ъ), которая

). Такое отличие от

формально определяет постоянную времени системы (2) Т^ = 1/ М(

фильтра, построенного методом ВЭФ, позволяет приближенно оценить коэффициент Косп величиной К « 10^(1/ Т\). Для фильтра (14) Т= = I— , что дает практически тот же результат К « Косп что и представленный на рис. 5.

Таблица 2. Соотношения параметров плотности (13)

№ кривой 1 2 3 4 5

a 0,05 0,1 0,2 1 5

Dt (a,p) 21,03 11,07 6,13 2,47 2,03

Выводы

1. Полученные выражения функционалов (6)-(8) имеют достаточно общий вид и могут быть использованы для решения задач обнаружения, распознавания, оценки и собственно фильтрации сигналов. В записи функционала (7) раскрыта его внутренняя структура, содержащая коэффициенты корреляционных связей (р^ (0), р., (0)) шумового процесса заданного СДУ

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

чиной T, = 1/ M (

).

2. Семейство интегральных уравнений фильтров максимального правдоподобия (табл. 1) отличается от известных наличием множителя , который «участвует» в образовании дискриминационной характеристики фильтра (/ (Ъ )/ '(Ъ )/ Ь) и приводит к появлению в уравнениях слагаемого А/ - текущей динамической разности.

3. Эффективность исследуемого фильтра, оцениваемая формулой (16), заключается в комбинации амплитудно-частотного подавления шума. Причем, такой эффект увеличивается с ростом корреляции помехи.

4. Результаты моделирования подтверждают теоретические расчеты, полученные в работе. Так синтезированный фильтр (14) остается не только работоспособным при малых входных ОСПвх, но и значительно выигрывает у фильтра на основе ВЭФ (рис. 5) при уменьшении ОСПвх. Дискриминационная характеристика фильтра (14) (рис. 3) отвечает практическим соображениям, имеет возможность адаптивной настройки (параметры а, в) и может быть реализована на базе дифференциальных управляемых каскадов усиления.

FUNCTIONALS AND FILTERS THE MAXIMUM LIKELIHOOD BASED ON STOCHASTIC DIFFERENTIAL EQUALIONS

A.V. AUSIANNIKAU

Abstract

The technique of filter design maximum likelihood based on SDE, qualitatively different from the classical optimal filters. Reviewed and analyzed the internal structure of functional likelihood cor-

responding to the SDE. A specific example is slowly fluctuating noise conditions shows the efficiency obtained in filters. The simulation results confirm the the theoretical calculations.

Литература

1. Стратонович Р.Л. // Изв. Вузов. Радиофизика. 1959. Т. 2, №6, С. 892-901.

2. Стратонович Р.Л. Избранные вопросы теории флюктуации в радиотехнике. М., 1961.

3. Амиантов И.Н. Избранные вопросы статистической теории связи. М., 1971.

4. Тихонов В.И., Кульман Н.К. Нелинейная фильтрация м квазикогерентный прием сигналов. М., 1975.

5. Фомин А.Ф. Аналоговые и цифровые синхронно-фазовые измерители и демодуляторы. М., 1987.

6. Овсянников А.В. // Докл. БГУИР. 2009. №8. С 13-21.

7. Чердынцев В.А., Овсянников А.В., Козел В.М. // Адаптивное устройство подавления помех Опубл. 15.05.1992. Бюл. №18.

8. Фомин В.Н. Рекуррентное оценивание и адаптивная фильтрация. М., 1984.