Учет сезонности в моделях авторегрессии

Ранее мы уже упоминали о том, что для моделирования сезонности в моделях авторегрессии существует специальная модификация. Она называется "SARIMA" ("Seasonal AR1MA") и все так же была подробно изучена в книге Дж. Бокса и Г. Дженкинса^[1]. В наши дни она иногда считается просто частью модели ARIMA. Рассмотрим ее форму записи, основные свойства и то, как можно с ее помощью получать прогнозы.

Идея моделирования сезонности с помощью моделей авторегрессии опирается на все те же лаговый оператор и зависимость текущего значения от предыдущих. Модель сезонной авторегрессии SAR(P)s. в первом приближении может быть записана так:

(8.64)

где – ошибка модели на наблюдении t; ФP – функция от оператора сдвига с лагом s, такая, что

(8.65)

где Ai – коэффициент при значении ряда на i-м лаге.

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

Функция (8.65) позволяет моделировать зависимость текущего значения не только от значения в прошлом сезоне, но и несколько сезонов назад.

Модель SAR(P)s во многом похожа на модель AR(р). Поэтому и различные свойства стационарности в ней в целом соответствуют моделям авторегрессии. Останавливаться подробно на этой модели не имеет смысла.

Определяется порядок сезонной авторегрессии с помощью коррелограммы по тем же принципам, что и простая авторегрессия, только при анализе нужно уже учитывать периодичность s. Так, при наличии значимого частного коэффициента автокорреляции на лаге Ps и дальнейшем обрыве ЧАКФ, имеет смысл построить модель SAR(P)s.

На рис. 8.25 приведен пример сезонного ряда данных (ряд № 1100) и его коррелограммы.

Как видим по графику и коррелограммам, в ряде данных наблюдается сезонность, причем квартальная. По автокорреляционной функции можно сказать, что она имеет затухающий вид, что сигнализирует о возможности построения SAR(P)s. При этом по частной автокорреляционной функции видно, что значимыми оказались коэффициенты на лагах 2 и 4. Скорее всего, в данных условиях достаточно будет построить модель SAR(1)4. Возможно, построение такой модели уберет автокорреляцию второго порядка. Если же она останется, то в дополнение к SAR(1)4 нужно будет построить модель с простой авторегрессионной компонентой.

Рис. 8.25. Ряд № 1100 и его коррелограммы

Если бы по коррелограмме на рис. 8.25 в частной автокорреляционной функции значимыми оказались коэффициенты на лагах 4 и 8 с последующим обрывом ЧАКФ, то имело бы смысл построить модель SAR(2)4.

Модель SMA((Q)s – сезонной скользящей средней – записывается следующим образом:

(8.66)

где

(8.67)

Сi – коэффициент при значении ряда на i-м лаге.

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

Как видим, эта модель также позволяет моделировать зависимость текущего значения от значений в нескольких сезонах назад. У этой модели те же условия обратимости, что и у модели МА(q), и тот же принцип прогнозирования: когда модель начинает опираться на расчетные значения, элемент скользящей средней становится равным нулю. Так, модель позволяет дать прогноз на s шагов вперед, опираясь на имеющиеся фактические значения.

Определение порядка модели сезонного скользящего среднего происходит по аналогии с моделью простого скользящего среднего, также "зеркально" по отношению к модели SAR(P)s, поэтому отдельно на этом вопросе мы останавливаться не будем.

Объединяя (8.64) и (8.66), мы получим модель SARMA(P,Q)s.:

Продолжая аналогии с оригинальной моделью ARIMA, построение модели можно осуществлять лишь в том случае, если ряд стационарен в сезонах, т.е. из сезона в сезон не наблюдается постоянной тенденции к росту либо какие- нибудь еще виды нестационарности. Для определения сезонной стационарности ряда существуют различные тесты. Одним из таких применяемых на практике является тест Канова – Хансен (Canova – Hansen test, СН test), который проверяет наличие тренда в сезонных компонентах^[2]. Он во многом напоминает KPSS-тест. Модель, на основе которой происходит тестирование, имеет вид

где хt – вектор регрессионных элементов, который может включать в себя любые факторы, в том числе и авторегрессионные; β – вектор коэффициентов при регрессорах; dt – вектор фиктивных переменных, характеризующих

(8.68)

(8.69)

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

(8.70)

где γ0 фиксирована и .

Нулевая и альтернативная гипотезы СН-теста формулируются идентично тому, как это было в KPSS-тесте:

Для проверки этой гипотезы так же рассчитывается LM-статистика па основе ковариационной матрицы коэффициентов при фиктивных переменных.

Если в результате проведения теста исследователь приходит к тому, что в ряде данных наблюдается сезонная нестационарность (что, например, может выражаться ростом амплитуды сезонных колебаний, то есть наличием мультипликативной сезонности), то для оценки модели SARMA требуется взятие сезонных разностей. В модель (8.68) они так же включаются по аналогии с моделью ARIMA:

Как видим, модель (8.71) очень похожа на модель ARIMA. Она называется чистой сезонной моделью в том случае, если ошибки модели ηt не автокоррелированы и распределены нормально. Однако на практике чистые сезонные процессы встречаются крайне редко, поэтому из-за того, что в такой модели не учтены элементы простой авторегрессии и скользящей средней, ошибки ηt будут автокоррелированы. Для исключения этой автокорреляции остатки ηt требуется описать моделью ARIMA:

(8.71)

откуда можно выразить ηt :

и подставить в (8.70):

чтобы получить окончательную модель SARIMA:

(8.72)

Если в модель нужно включить константу, то это делается путем добавления ее в правую часть по аналогии с простой ARIMA.

Модель (8.72) обычно обозначается как SARIMA(p,d,q) × (P,D,Q)S и чисто внешне состоит из элементов:

• 1) AR(р);
• 2) I(d);
• 3) МА("7);
• 4) SAR(P);
• 5) SI(D);
• 6) SMA(Q);
• 7) лаг сезонности s.

Такая форма записи, как (8.72), помимо прочего позволяет моделировать сезонность более сложного характера и "экономить" на коэффициентах. Покажем, например, что будет представлять собой модель SARIMA( 1,0,1) × (1,0,1)4, подставив параметры в (8.72):

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

(8.73)

Затем раскроем скобки в (8.73):

Далее выразим уt через остальные элементы, чтобы прийти к финальной записи модели:

.(8.74)

Как видим, в модели (8.73) помимо включенных элементов па наблюдениях t – 1 и t – 4 за счет перемножения появился еще и элемент t – 5, который в случае, если модель стационарна, будет меньше, чем c1 и С1 по отдельности. Включение его в модель позволяет более точно описать зависимость текущего значения от предыдущих и убрать из ошибок автокорреляцию. При этом, если бы мы оценивали модель (8.74) напрямую, нам нужно было бы оценить шесть коэффициентов вместо четырех. В этом как раз и заключается "экономия" на коэффициентах.

Заметим, что для построения SARIMA без сезонных разностей в распоряжении исследователя должно быть минимум 2s наблюдений. И то в таком случае фактически построение модели будет происходить лишь по второму сезону, так как данные из первого сезона становятся своеобразными "стартовыми значениями" и выпадают из аппроксимации. Взятие сезонных разностей в таком случае в принципе невозможно, да и выявить порядок модели становится крайне сложно. В таких условиях говорить о точности модели не приходится. Поэтому с практической точки зрения строить модель имеет смысл при наличии хотя бы 3s наблюдений (что, вообще говоря, так же мало). При прогнозировании обратимых процессов с помощью модели SARIMA сохраняется стандартный статистический принцип: чем больше данных, тем лучше.

В параграфе 8.3 мы определили, что по дифференцированному ряду № 2568 можно построить модель ARMA(3,3). Построим теперь по ряду № 2568 модель SARIMA. Ранее мы уже идентифицировали ряд элементов модели, запишем их:

• 1) s= 12;
• 2) D= 1;
• 3) р = 3;
• 4) q = 3;
• 5) d = 0.

Обратимся к коррелограммам по остаткам модели (8.61). Для удобства приведем их еще раз на рис. 8.26.

На этом рисунке обращают на себя внимание выбивающиеся за пределы доверительного интервала коэффициенты корреляции на 12-м лаге. После, на 24-м и 36-м лагах, затухающего характера АКФ и ЧАКФ не наблюдается. Поэтому для начата мы попробуем оценить модель с Р= 1. Если полученная модель не будет нас устраивать, мы еще раз обратимся к коррелограмме и попробуем вновь идентифицировать модель.

Рис. 8.26. Коррелограммы по ряду остатков по модели (8.61)

В итоге пошаговый подбор, рассмотренный нами в нескольких параграфах, приводит нас к тому, что по ряду в логарифмах нужно построить модель SARIMA(3,0,3) × (1,1,0)12:

Оценка этой модели дает значения коэффициентов, приведенные в табл. 8.4.

Таблица 8.4

Результаты оценивания модели SARIMA(3,0,3) × (1,1,0)12

Примечание: Coefficients – строка со значениями коэффициентов; s.e. строка со значениями стандартных ошибок коэффициентов; p-value – строка со значениями остаточных вероятностей по значимости коэффициентов.

Запишем эту модель в виде уравнения, округляя значения коэффициентов до сотых (в скобках указаны стандартные ошибки):

(8.75)

Как видим, в полученной модели незначимым оказался только один коэффициент – мри скользящей средней второго порядка. Однако исключать его из модели некорректно, так как из-за этого нарушится вся структура SARIMA. На рис. 8.27 приведены ряд остатков модели и его коррелограммы.

Рис. 8.27. Ряд остатков но модели (8.75) и его коррелограммы

Как видим, автокорреляции в остатках не наблюдается, в результате чего можно сделать вывод, что добавление каких-либо еще элементов вряд ли улучшит модель. Однако остатки все так же нельзя признать нормально распределенными из-за сильно выбивающегося остатка в апреле 1987 г. (наблюдение № 40 на рис. 8.27). Если обратиться к исходному ряду данных, то можно заметить, что на этот месяц пришелся пик, который действительно выбивается из общей картины. Если признать его случайным "выбросом" (значением, которое не может быть объяснено моделью по каким-то фундаментальным причинам), убрать из рассмотрения, переоценить модель, то тест Шапиро – Уилка дает остаточную вероятность в 0,2862, что уже больше 5% и не дает оснований отклонить нулевую гипотезу о том, что остатки распределены нормально. При прогнозировании на практике следует отдельно выяснить, что же произошло в апреле 1987 г., что именно привело к такой ситуации. Информация о причинах такого значения позволит принять правильное решение для дальнейшего моделирования и прогнозирования.

Дадим прогноз по модели (8.75). Графически он представлен на рис. 8.28.

Рис. 8.28. Фактические значения ряда № 2568, расчетные по модели (8.75) и прогноз

Заметно, что на последних трех значениях прогноз оказался менее точным, чем но предыдущим моделям. Возможно, из-за этого, а возможно, и из-за неточного прогноза и по нескольким другим точкам, но ошибка прогноза в модели (8.75) составила sMAPE = 5,38%, что больше, чем в случае с построением модели ARMA по ряду сезонных разностей логарифмов. Это, однако, не говорит о том, что модель SARIMA неприменима вообще. Возможно, что на других рядах данных она даст более точные прогнозы.

Здесь нельзя не вспомнить один из выводов, сделанных по результатам испытаний М3 – Competition: сложные, статистически обоснованные модели необязательно дают более точные прогнозы, чем простые^[3].

Авторы: Светуньков Иван, Светуньков Сергей

• [1] Box George Е. Р., Jenkins Gwilym M. Time series analysis, forecasting and control. Holden-day, Inc., 1976. P. 300.
• [2] Canova Fabio, Hansen Вruсе Е. Are Seasonal Patterns Constant over Time? A Test for Seasonal Stability // Journal of Business & Economic Statistics. 1995. Vol. 13. > 3. P. 237-252.
• [3] Makridakis S., Hibon, М. The М3 – competition: Results, conclusions and implications // International Journal of Forecasting. 2000. № 16. P. 451–476.