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

УПРАВЛЕНИЕ

------\ В СОЦИАЛЬНО-ЭКОНОМИЧЕСКИХ СИСТЕМАХ

УДК 519.258

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

Л. А. Осипов,

доктор техн. наук, профессор А. М. Кричевский,

аспирант

Санкт-Петербургский государственный университет аэрокосмического приборостроения

Рассмотрены модели временных рядов, характеризующихся наличием долговременной зависимости (долгой памяти). Для идентификации таких рядов предложено использовать модели класса ARIMA (p, d, q) с дробным показателем d. Показаны пути оценки параметра «памяти» временного ряда, решения задачи прогнозирования в таких рядах.

Time series models with long memory are considered. For the identification of such time series the ARIMA (p, d, q) model with fractional parameter d is proposed. Ways of estimation of the «memory» parameter and methods of forecasting in such time series are also studied.

Изучение временной структуры рядов различной природы, например экономика, телекоммуникация, астрономия и т. п., играет ключевую роль в моделировании и получении прогнозных оценок в различных областях науки и техники. В последнее время значительный интерес проявился к временным рядам (ВР), которые можно охарактеризовать термином «временные ряды с долгой памятью — time series with long memory» [1]. Существующие синонимы для этих рядов: долговременная зависимость (long-range dependence), сильная зависимость (strong dependence) или персистентность (persistence). Под рядами с долгой памятью понимаются не только стационарные ряды, но также и нестационарные, в которых зависимость от времени спадает очень медленно.

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

Долгая память обычно описывается в виде автоковариаций или спектральной плотности. Положим, что Xt, t = 0, ± 1,... является временным рядом. Если ряд стационарный, то E(Xt) = ц, еоу (Xt, Xt+j) = р(/) не зависят от t. В случае, если Xt имеет непрерывную функцию распределения, то его спектральная плотность выражается следующим образом:

1 ~

S(f) = -^ Е Р(/)е, -п<f <п. (1)

2п j=~

Здесь S(f) — неотрицательная четная функция с периодом 2п при ее продолжении за диапазон Найквиста [-п, +п].

Принято считать, что ряд Xt имеет долгую память, если

^0) = -1 Е Р0) =“, (2)

2п j=-.

т. е. S(f) имеет полюс на нулевой частоте.

Противоположная ситуация с нулевым значением спектральной плотности на нулевой частоте:

^0)=1 Ер( ]) = 0, (3)

2п j

определяет отрицательную зависимость или антиперсистентность.

Учитывая выражения (2), (3), можно сказать, что ряд Xt имеет короткую память, если

0 < 5(0) < го .

Кратко поясним возникновение новой модели, отталкиваясь от методологии Бокса—Дженкинса [2]. Наиболее распространенными моделями для стационарных ВР являются модели авторегрессии и скользящего среднего. Авторегрессионную модель порядкар, которая сокращенно обозначается AR(p) (autoregressiveprocess), можно записать в виде

Xt = <&1Xt-1 +Ф2Xt-2 + ... + ФрХ^p + at, (4)

где Ф1, Ф2,..., Фp — весовые коэффициенты.

Выражение (4) определяет процесс авторегрессии порядкар, в котором текущее значение ряда в момент t выражается через конечное число прошлых значений и величину возмущения at, не зависящую от прошлого. С помощью оператора сдвига B = Xt-1 /Xt модель (4) можно записать в эквивалентной форме

(1-Ф1В-Ф2B2 -...-ФpBp)yt =at,

которая после введения оператора авторегрессии Ф(Б) принимает вид

Ф (B)Xt = at. (5)

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

Xt = at -01 at -1 - ... -0qat-q , (6)

где символы 01, . . . , 0q используются для обозначения конечного набора весовых параметров.

Соотношение (6) определяет процесс скользящего среднего порядка q, или сокращенно МА (q), который представляет собой линейную комбинацию текущего и прошлых значений шума at. Используя оператор сдвига, можно записать для процесса (6) эквивалентное выражение

Xt = (1-01B -... 0qBq)at - 0(B)at. (7)

Моделями AR (р) и МА (q) за счет выбора их порядков р и q можно удовлетворительно описывать многие реальные процессы. Однако на практике для достижения большей гибкости в подгонке моделей к наблюдаемым ВР иногда целесообразно объединить в одной модели и авторегрессию, и скользящее среднее; при этом цель должна состоять в построении наиболее экономных моделей, дающих хорошую аппроксимацию с помощью небольшого числа параметров. Достижению этого помогает рассмотрение смешанных моделей авторегрессии — скользящего среднего, т. е. моделей ARMA (р, q), которые имеют вид

или в другой эквивалентной форме через операторы Ф(В), 0(B)

Ф(В)у =0(B)at. (9)

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

Существует неограниченное число различных проявлений нестационарности. Однако можно выделить обширный класс встречающихся в приложениях ВР со специфической однородной неста-ционарностью, которая удовлетворительно описывается стохастической моделью, являющейся модифицированной формой моделиARMA. Условие стационарности модели (9) означает, что корни полинома Ф(.В) лежат вне единичного круга [2]. Естественный путь получения нестационарного процесса, описываемого таким же уравнением, заключается в ослаблении этого ограничения. В частности, оказывается, что во многих случаях наблюдаемые в реальности процессы хорошо описываются моделями типа (9), у которых один или несколько корней Ф(.В) равны единице. Такой класс моделей называется процессами авторегрессии — проинтегрированного скользящего среднего. В английской аббревиатуре такой процесс запишется как ARIMA (с добавлением к уже известному сокращению ARMA слова integrated).

Рассмотрим модель

<p(B)Xt =0(B)at, (10)

где в отличие от равенства (9) ф(-В) — нестационарный оператор авторегрессии порядкаp + d, такой, что d корней уравнения ф(-В) = 0 равны единице, а остальные р корней лежат вне единичного круга; оператор же скользящего среднего 0(В) по-прежнему обладает порядком q и является обратимым (все его корни лежат вне единичного круга).

Тогда можно записать, что

p(B) = Ф(В)(1 - B)d,

где Ф(.В) — уже стационарный порядкар оператор авторегрессии (т. е. с корнями вне единичного круга). Если ввести оператор разности А = 1 - B, AXt = = Xt - Xt-1, то ф(.В) запишется как p(B)AdXt = 0(B)at и модель (10 ) можно представить в виде

Ф^А^ =0(B)at. (11)

Здесь d-я разность ряда Xt вычисляется по формуле

wt =AdXt =(1 -B) dXt и, следовательно, удовлетворяет уравнению

Xt = ФА-і + ... + Ф pXt-p + at - 6!at-L - .

.. 6qat-q,

(В)

Ф(B)Wt =6(B)at,

(і2)

т. е. является уже стационарным обратимым процессом ARMA (p, q).

Таким образом, процесс Xt, описываемый уравнением (11), можно получить d-кратным суммированием (или интегрированием) процесса {w}, являющегося в соответствии с (12) процессом ARMA. Вследствие этого процесс, задаваемый моделью (11), называют процессом ARIMA. Если в формуле (11) оператор авторегрессии Ф(В) имеет порядок р, а оператор скользящего среднего 0 (В) — порядок q, то кратко модель (11) записывается как ARIMA (p, d, q). В частности, при d = 0 эта более общая модель сводится к смешанной модели ARMA (p, q). Тем самым модель ARIMA (p, d, q) охватывает широкий класс как стационарных (при d = 0), так и нестационарных (при d > 1) процессов.

В работах [3, 4] было впервые предложено рассмотреть дробные значения d из интервала d е (-1/2, 1/2), что привело к дробной (fractional) авторегрессионной модели скользящего среднего порядков p, d, q (ARFIMA (p, d, q) или FARIMA (p, d, q)).

Примем, что Xt удовлетворяет следующему дифференциальному уравнению:

(1 - B)dXt = at, at - N (0, o2a). (13)

При расширении показателя d в (13) до нецелых степеней результатом является ряд с преобразованием, которое включает в себя разложение члена (1 - B)d по биномиальной теореме для нецелых показателей:

(і - B)d = j (-і)k і d 1 Bk,

k=Q

где

' d ^= d(d-1)Сd-2) •••(d - k-і)

k k!

V J

Применяя это разложение к Xt, получаем

(і - B)d Xt = j (-і)k {d ] BkXt = j AkXt-k =

k=Q

V

k 4 _ _

"t,

k=Q

где коэффициенты авторегрессии Лk выражаются через гамма-функцию:

Лі, = (-«* і ^ 1- r(k+ d)

k

V .

r(-d)r(k +1)

Грейнджер и Хоскинг показали [3, 4], что характеристики таких временных рядов обладают важными свойствами: например, Хг является стационарным и обратимым для й є (-1/2, 1/2). Кроме того, оказывается, что положительная или отрицательная зависимости определяются знаком при параметре й, т. е. автокорреляционные коэффициенты процесса Х{ имеют тот же знак, что и й.

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

Простейшей реалистической моделью для стационарного ряда является параметрическая модель, которая выражает p(j) для всех j и S(f) для всех f как параметрическую функцию только двух параметров: d и неизвестного масштабного фактора. Возможно, что самой ранней такой моделью явился фрактальный шум, концепция которого возникла из рассмотрения самоподобия.

Определим самоподобный процесс следующим образом [1]: непрерывный стохастический процесс {y(t); -го < t < го} является самоподобным с параметром самоподобия Hе (0; 1), если для любого a > 0 процесс {y(at); -ro<t<го} имеет то же распределение, что и процесс {aHy(t); -ro<t <го}.

У каждого реального самоподобного процесса должен быть наибольший и наименьший масштаб: нельзя бесконечно увеличивать или уменьшать масштаб. Например, в случае броуновского движения (БД), представляющего образец самоподобного процесса, диапазон масштабов, в пределах которого сохраняется самоподобие, охватывает много порядков величины: от размеров сосуда с жидкостью (допустим, 0,1 м) до длины свободного пробега молекул между столкновениями, которая для малых частиц может достичь 10-9 м. Во многих случаях объект считается самоподобным, если его можно масштабировать с коэффициентом подобия 10 и меньше (до трех дискретных шагов) [5]. Оценим спектральную плотность броуновской функции S(f), которую определим как проекцию БД на одно пространственное направление в зависимости от времени. БД порождается независимыми приращениями и имеет плоский спектр. Следовательно, сумма (интеграл) приращений обладает спектральной плотностью, пропорциональной f ~2. Отметим, что в общем случае зависимость спектральной плотности от частоты характеризуется степенным законом вида f ~в. Среди шумов большой известностью пользуется белый шум со спектральным показателем в = 0. Иначе говоря, спектр белого шума не зависит от частоты. Проинтегрировав белый шум один раз по времени, получаем коричневый шум (проекцию БД на одно пространственное измерение), который имеет спектральную плотность, пропорциональную f ~2. Но белый и коричневый шумы далеко не исчерпывают все спектральные возможности: между ними располагается розовый шум со спектром f _1, а за коричневым — черный, пропорциональный f ~в, где в > 2.

В качестве примера самоподобного процесса приведем результаты моделирования коричневого и розового шумов. Генерирование коричневого шума сводится к суммированию независимых случайных чисел и реализуется сравнительно легко посредством табличного расчета в Excel. На рис. 1 приведена реа-

Значения

Спектральная плотность 3000

2000

1000

0

реализации

0,20

0,40 Частота

Рис. 1. Реализация процесса коричневого шума вместе с автокорреляционной функцией и спектральной плотностью

1

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

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

Релаксационный процесс с дискретными значениями времени хп можно задавать с помощью генератора случайных чисел, который позволяет получать независимые случайные числа гп, подставляемые затем в рекуррентное соотношение [5]

Хп+1 =РХп +41 -Р2 гп, х0 =0, (14)

где Р — требуемый коэффициент корреляции между соседними случайными значениями. Со временем релаксации т этот коэффициент связан соотношением р = ехр(-1/т). Таким образом, для набора значений времени релаксации, каждое из которых в 10 раз превосходит предыдущее (т =1, 10, 100, ...), коэффициенты корреляции получаются вычислением последовательных корней десятой степени (т. е. р = 0,37; 0,90; 0,99, ...). Результаты

расчета, выполненные по (14), приведены на рис. 2 вместе с автокорреляционной функцией и спектральной плотностью.

При использовании модели класса ЛЯ¥1ША (р, й, q) важно правильно определить параметры этой модели. Неверное определение параметров р и q приводит к несогласованной оценке коэффициентов АВ- и МА-моделей, но ошибка в оценке й дает неверную интерпретацию обеих моделей из-за потери идентификации. Асимптотическое поведение спектра указывает на то, что при моделировании рядов с короткой памятью несущественно влияние очень низких частот и очень длинных лагов, т. е. в ситуациях доминирования й. Вследствие этого становится понятно, что оценки й должны основываться на информации о низких частотах или длинных лагах. С точки зрения требований устойчивости, оценки должны быть основаны на очень малой части данных при увеличении объема выборки данных, поэтому естественно ожидать более медленной скорости сходимости, чем для оценок, полученных по модели с целыми значениями параметров. Однако в очень длинных временных рядах (в экономике, метеонаблюдениях) доступное число степеней свободы может быть достаточным для обеспечения адекватной точности. Такие оценки обычно называются полупараметрически-ми (semiparametric) [1].

Значения

АКФ

1

0,5

141

281

421

561

Спектральная плотность

400

300

200

100

10 Смещение

701

841 Номер реализации

0 0,05 0,15 0,25 0,35 0,45 Частота

0

5

■ Рис. 2. Реализация процесса розового шума вместе с автокорреляционной функцией и спектральной плотностью

Существует несколько методов для оценки параметра й, являющихся в то же время и тестами для обнаружения долгой памяти во временных рядах. Параметр й можно оценить во временной или частотной областях [6].

Во временной области для получения оценки й используется асимптотическое представление для ковариаций [1]

РО) ~С1 ]М-1 при ] ^го,

где с1 > 0 — постоянная; ] — временной лаг.

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

• нелинейная регрессия выборочных автоковариаций;

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

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

Однако распределения этих оценок достаточно

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

В частотной области одна из первых оценок й была предложена в работе [7]. Сущность метода заключа-

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

Для экспериментальной проверки изложенного в отношении ВР с долгой памятью здесь воспользуемся данными из экономики за 19531992 гг. Изменение выбранного параметра во времени показано на рис. 3. Характер автокорреляционной функции и периодограммы, являющейся оценкой спектральной плотности, указывает на то, что рассматриваемый временной ряд, скорее всего, обладает долгой памятью и характер его изменения можно описать модельюАВ¥1МА (р, й, q).

Для анализа ВР, показанного на рис. 3, воспользуемся модулем динамического моделирования PcGive из программного пакета GiveWin2.20. Этот модуль позволяет идентифицировать и тестировать модель и оценивать ее параметры.

Значения

АКФ

1,G

G,5

10 Сдвиг

Спектральная плотность 15GG

1GGG

5GG

G

'Л-

J

1,G Частота

Рис. 3. Временной ряд, автокорреляционная функция, периодограмма

5

5

Для оценки параметров модели в данном модуле используются точный метод максимального правдоподобия (Exact Maximum Likelihood — EML) и нелинейный метод наименьших квадратов (Nonlinear Least Squares — NLS). Решение задачи выполняется в форме последовательного диалога за несколько шагов, в частности: выбор данных для анализа, установка параметров модели, оценка параметров модели. Для данного ряда была выбрана модель ARFIMA (1, d, 2). Параметр d оценивался с помощью нелинейного метода наименьших квадратов. Результаты расчетов параметров модели приведены в таблице.

Значения AR-1 и MA-1 не приводятся, так как при редактировании модели были выбраны фик-

■ Результаты расчетов

b Results uwm\

k Non-linear least squares estimation of ARFIMA(1 , d, 2) model л

The estimation sample is: 1960 (2) - 1990 (1)

The dependent variable is: CONS (data. in7)

Coefficient Std.Error t-value t -prob

d parameter 0.219975 0.1320 1.67 0. 098

MA-2 0.19 8 9 08 0.09466 2.10 0. 038

COHS_l 0.935435 0.05361 17.4 0. 000

Constant 56.2894 46.90 1.20 0. 233

no. of observations 120 no. of parameters 5

AIC.T 529.027308 AIC 4 .4085609

mean (COWS) 874.17 var(CONS) 177.849

sigma 2.12141 sigmaA2 4.50036

BFGS using numerical derivatives (epsl=0.0001; eps2= 0.005): 1

Strong convergence ■

Used starting values:

0.00000 0.26390 0.98673 11.433 -

< >

а)

Уровень

б)

Плотность

вероятности

-3 -2 -1 G 1

2 3 Интервал

Рис. 4. Исходный (-----) и подобранный (--) ряды (а) и гистограмма остатков и кривая нормального распре-

деления (б)

Уровень

■ Рис. 5. Прогноз ряда на 4 интервала времени вперед:-------— исходный ряд;---— подобранный ряд

сированные лаги для AR-1 и MA-1, равные единице. На рис. 4, а приведены исходный ВР и подобранный, который описывается моделью ARFIMA (1, d, 2) при d = 0,22. На рис. 4, б приведена гистограмма остатков для проверки адекватности модели. Подобранную модель можно использовать и для прогноза. Например, на рис. 5 видно, что предсказанные значения сохраняют тенденцию ряда. Таким образом, в работе показана возможность использования моделей ВР с долгой памятью для анализа и прогнозирования.

Литература /

1. Time Series with Long Memory / Ed. P. M. Robinson. Oxford University Press, 2003. 2. Бокс Дж., Дженкинс Г. Анализ временных рядов. Прогноз и управление. Вып. 1. М.: Мир, 1974. 3. Granger C. W. , Joyeux R. E. An introduction to long-memory series models and fractional differencing // Journal of Time Series Analysis. 1980. Vol. 1. P. 15-29. 4. Hosking J. R. M. Fractional differencing // Biometrika. 1981. Vol. 68. P. 165-176. 5. Шредер М. Фракталы, хаос, степенные законы / НИЦ «Регулярная и хаотическая динамика». Ижевск, 2005. 6. Breidt F. J., Crato N., Lima P. The detection and estimation of long memory in stochastic Volatility // Journal of Econometrics. 1998. Vol. 73. P. 325-348. 7. Geweke J., Porter-Hudak S. The estimation of long memory time series models // Journal of Time Series Analysis. 1983. Vol. 4. P. 221-238.

/