Прогнозируемость и интервальное прогнозирование нестационарных процессов на основе стохастических дифференциальных уравнений Текст научной статьи по специальности «Математика»

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

2013 № 7 (77)

УДК 519.2:005

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

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

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

Поступила в редакцию 30 мая 2013

Рассмотрена задача интервального прогнозирования нестационарных процессов описываемых моделями стохастических дифференциальных уравнений. Определена прогнозируемость таких процессов. Получены алгоритмы интервального прогнозирования в дискретном и непрерывном времени.

Ключевые слова: нестационарный процесс, интервальный прогноз, алгоритм.

Введение

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

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

Информационная прогнозируемость и задача ИП

Рассмотрим вначале задачу ИП в непрерывной форме. Модель наблюдения имеет вид y(t) = f (t) + ^(t) , t e[t0 ,T], где y(t) рассматривается как сумма детерминированной f (t) = f (t,A) (Л - вектор неизвестных постоянных параметров тренда) и случайной составляющей ^(t) - «шума» наблюдения описываемого СДУ:

m+a(t^)=g(tm,m=^, (1)

где a(t, g (t) - функции, удовлетворяющие условию Липшица; C,(t) - нормальный белый шум с нулевым средним (С(0) = 0 и дельтообразной корреляционной функцией

(C(t)C(t -т)) = N5(i)/2 , b(t) = Hm(N /2At)£+^g2(x)dx - диффузия (в случае g = const

диффуззия Ь = ^2/2), N - односторонняя спектральная плотность. Решая уравнение (1) известными аналитическими способами, может быть получена нестационарная плотность (НП) Р = Р(?,Е, 12), (Р(?0,Е, 12) = ро, Р(?, ±до | 2)=0), где 2 - набор неизвестных параметров НП.

Верхний а+ (?) и нижний а- (?) доверительные пределы процесса Е(?) определяются из уравнения р (а—) -ф(?) < м <а+(?) -ф(0) = у, где м = у -ф(?, /, т ), ф(?) - систематическая составляющая наблюдаемой переменной у(?), т (?) - в общем случае не равное нулю математическое ожидание Е(?), У - заданная доверительная вероятность. Следовательно, выражение для верхнего и нижнего предела а± (?) = ф(?) ± Е— 2,БЕ (?, 2,0) ±у / 2J, которое при симметричной НП преобразуется к виду а± (?) = ф(?) ± Е 1 [?, 2, У /2], где Б— - функция, обратная функции распределения Е (?, 2, а± (?) -ф(?)) = У /2 . Таким образом, формируя некоторым известным способом оценки X* — [Л*, 2* ] на интервале наблюдения ? е[?0,Т], оценка функции а± (?') для момента времени > Т принимает вид

а±*(?') = q*(t') ± Е-1 [t2*, у/2] , qT (t') = q(t/(Л ), ^ (2; )). (2)

Сформулируем следующее утверждение. Если система оценок X* обеспечивает их состоятельность, асимптотическую эффективность, X* ^ X по вероятности, то и оценки непрерывных по t функций сходятся: а± (t',X*) ——>а±(t',X) . Таким образом, Р(а* (t') < y(t') < а+* (t')) = у ^ Р(а- (t') < y(t') < а+ (t')) = у, что соответствует постановке задачи ИП в непрерывной форме. Ширина интервала прогноза Ду(t') = 2F^1 [t',2*,у /2j .

Следовательно, формируя НП, определяющую однозначное интегральное преобразование Е 1,

решение задачи ИП сводится к вычислению (2) в моменты времени t' > T.

Будем называть информационной прогнозируемостью (в дальнейшем прогнозируемостью) стохастического процесса ^(t) по Фишеру неотрицательную функцию

Е(t) lnр /= J(ölnр /, t > 0, представляющую собой зависимость от

времени, при выполнении условий регулярности для р . Физический смысл функции Е (t) следующий. Для чисто диффузионного процесса (a(t, = 0, g(t) = g ) прогнозируемость L (t) = (btи lim L (t) = 0, а для процесса, описываемого линейным СДУ (a(t, £) = ц£,, g(t) = g ), прогнозируемость L (t) = а-2 (1 - е~2ц)-1 и limL (t) = а-2 = 2ц / b = const.

Покажем, что для СДУ (1) с параметрами a(t, £) = ц^(0, g(t) = gtv, g, t, ц> 0, v > 0 информационная прогнозируемость стохастического процесса определяется функцией

E (t) = 2v (-ц^ ^e^ (Г(1 + v, -2ц) - vr(v))]-1, ^ = b/2ц, b = #g2/2, (3)

i1 да

ta- e~ dt - неполная гамма-функция.

X

Действительно, коэффициенты сноса и диффузии равны: K(t,Q = -|j^(t), K2(t) = btv. Тогда уравнение Фоккера-Планка-Колмогорова принимает вид

( Y

5Р /5х = Р + £йР /5^ + а? - 52Р /, т = ц, а2 = b/2ц .

Решение этого уравнения с начальными и граничными условиями Р(0, Е) = 5(, —Ео). Р(х, ±да) = 0 имеет гауссовскую форму нестационарной плотности

Г-—-exp

ржЩ

(I ~ )2 2D (t)

D (t) = 2-v (-цГ [Г(1 + v, -2ц) - vr(v)], (4)

p (t) = lnp / ^ = (t) 1. Таким образом, формула (3) доказана.

следовательно р(t) = ^(ölnp / у = D"4 1

Из выражения (3) следует: во-первых, чисто диффузионный нестационарный процесс с параметрами ц = 0, g (t) = gtv обладает информационной прогнозируемостью L (t) = limD. (t)-1 =(1 + v) t-(1+v) / b; во-вторых, алгоритм ИП на фоне шумового процесса, описываемого СДУ (1) согласно уравнению (2), и с учетом (4) имеет вид

+* аг

(t')=Ф; (t') ^V2D(t)Erf - Чу), Ф; (t')=f(t ■ л; ). (5)

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

использовать конструктивный подход - последовательное формирование НП из одномерных одношаговых плотностей перехода (ОПП), определяемых оптимальными разностными схемами [6] СДУ (1) в форме Ито:

= ^ -ДлгЧ +AgiC,i, (6)

где Д = ^ - интервал времени, на котором слева и справа интегрируется СДУ (1),

ам = а(р1, ^ ^, й,. = (1/ А)| ' , = -Да^ - математическое ожидание и дисперсия

^ = Ь(^ )Д. Глобальная среднеквадратическая погрешность разностной схемы (6) определяется

/ мТ

\2

величиной [6] g<a(JqM (a^g(t))

1/2

dt) , T = kA.

Iт

.10

Таким образом, в дискретной форме рассматриваемая задача ИП на г шагов вперед в условиях регулярного статистического эксперимента состоит в определении доверительных оценок а±*г = ю± (у) (ю - борелевская функция) таких, что F% (ак+г е|а^+г, <х+к+г}) = у /2,

У = У ,•••, Ук]Т, У = Я , г = 1,к . Уравнения (2) преобразуются к виду

а^) = фТ(tk+r)±F? \tk+r,2,у/2], ф!(tk+r) = ф^г,f(Л*),m(2k)). (7)

Если система оценок X* обеспечивает их состоятельность, асимптотическую эффективность, асимптотическую нормальность, то по обобщениям теоремы Слуцкого для момента времени t' > Ц рациональная функция а± (I',X*) сходится к функции а±(t',X), т.е. а±*^',X*) ——^а± ^',X) и Р(а*^') < у^') <а+*^')) = у^ Р(а- (t') < y(t') <а+ (t')) = у, что соответствует постановке задачи ИП.

Формирование НП и последовательные алгоритмы идентификации

Представим многомерную НП относительно у-Г в виде

Л © = С,1 ехр Bhl-l ^, л, © = р0, (8)

где Б.I-1 = Б^, ^г-1, Ег. ] - семейство параметрических функций, = ) - значение параметра НП в 7-ый момент времени, Б10 = Б[^, , ^ ]. В дальнейшем, для упрощения будем рассматривать временную зависимость параметра Ег в виде Ег- = . Если выбрать функции Вг,г-1 =- 1п(с ,г-1лг,г-1), где лг,г-1 - ОПП, еи-1 - константы, связь между выражением (8) и

к

выборочным эмпирическим функционалом (ВЭФ) определяется зависимостью Wk = ^к .

,=1

В этом случае функция В; ^ представляет собой частную информационную функцию потерь (ФП), выбор которой определяет обратное преобразование в (6) и качество оценок

параметров плотности X*. Положив Ск = (^^^Г1^ ) , получим щ = ^^ Щ м .

В случае последовательного поступления информации применимы рекуррентные оценки, так, например, для ВЭФ Щ имеем Wi = Щ_1 (Д.;_1 - Щ_1), где п = 1, 0<п, < 1,

Пк ^ 0 и ^ п = ю при к ^ да. В случае ^ = 1/, этот алгоритм определяет формулу ВЭФ Щ приведенную выше. Для обобщенной ФП В =_ 1п(СЩ), используя (8) получаем рекуррентную формулу В. = В._1 + В,.,_1.

Целесообразность выбора того или иного метода оценки вектора параметров плотности X* = [Л*, Е*]г определяется постановкой задачи и объемом априорной информации. Так, например. применение стохастической аппроксимации к (8) приводит к последовательным оценкам X* = X0 — КВг|х_х* , где Ух = д / ЭХ, - найденное одним из известных

способов предварительное значение оценки; Кг - матрица, определяющая конкретный вид стохастической аппроксимации; начальные условия для уравнения X10, К выбираются исходя из имеющейся в наличии априорной информации. Приведенные выше уравнения для переменных Вг, X* образуют замкнутую систему идентификации НП (8) и оценки ее параметров. Конкретизация уравнений обеспечивается выбором вида функции тренда / и параметров масштаба Е исходя из практических целей задачи прогнозирования.

Алгоритм ИП на основе конструктивного подхода формирования ОПП

Выбор частной ФП связан с выбором ОПП. Воспользуемся конструктивным подходом формирования ОПП, формально заменяя в одномерной плотности Р(Е | Е) переменную

^ ^Е,- _и , требуя при этом сохранения условий нормировки. Строго говоря,

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

В качестве такой конструкции, совпадающей по форме с обобщенно-нормальным распределением, рассмотрим ОПП следующего вида:

щ,._1 = В(т)ехр(_Л(т)(|Е, _Му|/^Ц.)"), т е{1;2} , / = 2Л (9)

где

Л(т) = (Г(3 / т)Г_1 (1 / т))"'2 , В(т) = Л(т)1/т / (^/Ц~Г(1 +1/ т) и Ц = Е, = .

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

к

(г,)=((д 1п я. .._1 / эе, )2=Л(д 1п Щ, _1/ эе, )2р(е,, Е^а е, а .

Можно показать, что при использовании конструкции ОПП (9) значение функции „ ±

распределения случайной величины у , зависящее от параметра аг , имеет вид

Е (а*) = 0,5

1 - Г(1 / т)1 / т, А(т) (|а± + ср, | / ^

(10)

Алгоритмы, параметры и характеристики ИП, оценки параметров НП имеют вид выражений, приведенных в таблице. Обозначено: Л^. =*0 -M*i, = у - f0,

M*\ = **-1-Ла,*-!, С = у-! -f-1, f0 = f(A?), = Ж.,), Vf = f/5Л, К - матрица

коэффициентов усиления, Erf-1 - функция обратная интегралу вероятности ф* = f + M* ..

Формула (10) получается интегрированием (9) по у . Частная ФП BiJ-x =- ln(ciji-1wiji-1),

ее производная V*iBiJ-1^; и одношаговая информационная прогнозируемость р (t)

определяются аналитическими вычислениями по определенным выше формулам. Явная выборочная оценка параметра D, по результатам k наблюдений, получается применением МПП -5 ln щ (£)/5D = 0, где щ определяется выражением (8). Полученная таким образом МПП-оценка состоятельна, асимптотически эффективна, асимптотически нормальна и позволяет использовать итерационные процедуры ее вычисления. Алгоритм ИП приведенный в строке 5 таблицы непосредственно следует из (1 0) при подстановке конкретного значения параметра m = {1;2}.

Параметры и характеристики алгоритма ИП

№ Характеристика Выражение

1 Частная ФП В , г,г—1 A(m)-M*)| /JD. )m

2 A(m)mD -m/2 m-1 sign(^ )

3 Информационная прогнозируемость (3 - m)/Dt

4 Выборочная оценка параметра О / \2/m D; = \mA(m)k-Yi=i(\ * - M*. )| / fif J

5 Верхний и нижний интервальный пределы а++г 1Ф*+r ln(1-y)/V2, m = 1, [ф+, ±^12D*;tk+rErf"'(y), m = 2

6 Верхний и нижний интервальный пределы а+*(?'), Г> Т jV(tln(1 -y)/ -Д, m = 1 [ф* (t) ± ^ 2D*t 'Erf-1 (y), m = 2

7 т*(г) JO, m = 1 J* exp(-bt /2D*), m = 2

Следует учесть, что при вычислении функции ф*+г для моментов времени ^+г > ^, когда данные наблюдения временного ряда отсутствуют, следует использовать статистическое

* * I * \ * *

усреднение фк+г = /к+г + = fк+r + т*(к + г) . Рассматривая уравнения относительно

а,к+г для двух последовательных моментов времени , и переходя к пределу при

параметре Д^О, получаем обобщенные уравнения оценки верхнего и нижнего доверительного предела в непрерывной форме (строка 6 таблицы). Величина я^(г')

определится решением однородного дифференциального уравнения + = 0 (строка 7 таблицы).

Алгоритм ИП на основе марковской модели ОПП

Форма ОПП марковского процесса имеет вид

= ехр{—(2О0— "М^]2}/фЩ . (11)

Одношаговая информационная прогнозируемость для ОПП (11) по определению имеет вид Lq (t) = (Dti )—1. Дискретный алгоритм ИП для процесса с ОПП (11) следующий:

af+r = .Cr + m(k + r) ±V2D*tk+rErf—), D* = Ь*А, (12)

где по результатам обработки к наблюдений определяется величина Ь* = (1/ kА^2 / Aty .

Для доказательства (12) рассмотрим СДУ 4(0 + а(Е)= S'V^CCO и разностную схему ЕУ = ЕУ—1 — Ааг_1 + g [' \[tC,(t)dt, для которой определяем параметры сноса K = —аг-1, диффузии

J ty —А

K2 = btt, t = (t + t-i) / 2, среднего значения M^ { = Е,-^ — Aa^j и дисперсии D = Dty , D = ЬА.

Значение функции распределения случайной величины y ОПП (11), зависящее от параметра af, имеет вид (10) при m = 2, откуда следует = Ф*+Г ±^/2D*^~Erf—1(у). Функция фк+r состоит из суммы двух составляющих: fr+r и M*k+r. Поскольку будущие значения шума в к+r момент времени неизвестны, то вместо величины M*k+r используем ее математическое ожидание m*(k + r).

МПП-оценка — дln щ (£)/ CD = 0 параметра D, получающаяся из (11), с учетом представления (по результатам обработки к наблюдений) следующая:

D* = (1 / кА^2 / ty . Тогда Ь* = D* / А = (1 / k)^*=1 А^2 / Aty. Таким образом, алгоритм, представленный формулой (15), доказан.

Устойчивые к виду ОПП шума алгоритмы ИП

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

соотношением а(г, Е,-) = _(Ь / 2)( д 1п я,.^ / д,,) то, как показано в [5], выбирая «наихудшую»

*

в некотором заданном классе ОПП Щ ,_1 (одинаковую V £■), для которой справедливо неравенство 3 [В , я

]>3[В , я] (или 1,[В , я ]<1,[В , я]) можно обеспечить построение

робастного (по отношению к виду используемой ОПП) интервального прогноза. В частности, приведенные в таблице 1 алгоритмы, являются робастными. Так, для т = 1 это класс всех невырожденных ОПП (9), когда известно только то, что она существует ^ = |я: я(0) > 8 > 0|.

Класс ОПП (9) с ограниченной дисперсией в сечении =|я )< Ц | соответствует

параметру т = 2.

Количественно робастность для приведенных классов может быть выражена величиной связанной с одношаговой информационной прогнозируемостью следующим образом: 1,[В* ,_1, я, ,_1] > I,(г,) = (3 _т) /Цг,. Это неравенство говорит о том, что при любой другой

плотности из заданного класса, кроме «наихудшей», величина одношаговой информационной прогнозируемости процесса будет больше.

Результаты моделирования и их обсуждение

В качестве краткой иллюстрации рассматриваемого подхода на рисунке представлены результаты моделирование алгоритмов (у = 0,9) приведенных в таблице (строка 5; т = 1,2). Для алгоритма А1 - т = 2, А2 - т = 1. Для моделирования была выбрана модель тренда

/(г) = 100 +102 г и шума заданного уравнением (1) с параметрами а (г,,) = 103 ,(г), g(г) = •*/?,

^ = 1. Зависимости верхнего и нижнего доверительного предела ^), полученные

алгоритмами А1-А2, строились от точки .

Анализ результатов моделирования показывает, что учет «тонкой» структуры шума позволяет: во-первых, более точно настроить диапазон ИП в сравнении с ИП

полученным алгоритмом А2; во-вторых, алгоритм ИП чувствителен к модели шума, поэтому при отсутствии сведений о структуре его модели целесообразен выбор устойчивых (робастных) алгоритмов ИП.

600 500

Г Прогноз

А2 . -

АО

>"4 ^щГ"

ju*-

1000 2000 3000 4000 5000 6000 7000 8000 9000 10000

Результаты моделирования алгоритмов ИП

400

300

200

100

0

0

t

Заключение

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

PREDICTABLE AND NON-STATIONARY PROCESSES OF INTERVAL PREDICTION BASED ON STOCHASTIC DIFFERENTIAL EQUATIONS

A.V. AUSIANNIKAU

Abstract

The task of interval prediction of non-stationary processes of stochastic differential equations described by models is considered. Predictability of such processes is defined. Algorithms of interval prediction in the discrete and continuous time are received.

Список литературы

1. Орлов Ю.Н., Осминин К.П. Нестационарные временные ряды: Методы прогнозирования с примерами анализа финансовых и сырьевых рынков. М., 2011.

2. Brahim-Belhouari S., BermakA. // Computational Statistics & Data Analysis. 2004. № 47. P. 705-712.

3. Лукашин Ю.П. Адаптивные методы краткосрочного прогнозирования временных рядов. М., 2003.

4. Давнис В.В., Тинякова В.И. Адаптивные модели: анализ и прогноз в экономических системах. Воронеж, 2006.

5. Овсянников А.В. // Изв. НАНБ. Сер. физ.-мат. наук. 2010. № 4. С. 21-28.

6. Никитин Н.Н., Разевиг В.Д. // Журн. выч. матем. и математич. физики. 1978. Т. 18, № 1. С. 106-117.