CC BY-NC-ND
358
ЛЕКЦИОННЫЕ И МЕТОДИЧЕСКИЕ МАТЕРИАЛЫ
Анализ временных рядов Канторович Г.Г.
В этом номере продолжается публикация курса лекций «Анализ временных рядов». В прошлом номере были рассмотрены основные определения понятия стационарных случайных процессов, их характеристики и свойства, а также класс стационарных случайных процессов типа ARMA и нестационарных - типа ARIMA. Проанализирован подход Бокса-Дженкинса к идентификации временных рядов. В настоящем номере продолжается рассмотрение подхода Бокса-Дженкинса, в частности, оценивание параметров типа ARIMA и прогнозирование с помощью этих моделей, приведены примеры применения подхода Бокса-Дженкинса.
Рассмотрены особенности поведения некоторых нестационарных временных рядов, нестационарные ряды типа TS и DS, тест Дикки-Фуллера для проверки гипотезы о типе ряда.
Лекция 5
Оценивание коэффициентов моделей типа ARMA
Для оценки коэффициентов модели AR(p) Xt = alXt_1 +... + apXt_p + et можно применить обычный метод наименьших квадратов (МНК). Поскольку регрессоры относятся к предыдущим моментам времени, а st - белый шум, то корреляция
регрессоров со случайным возмущением et отсутствует. МНК не дает несмещенных оценок, но дает состоятельные оценки, несмотря на то, что присутствует стохастический регрессор. Разумеется, мы должны тщательно проверить по остаткам, действительно ли наши предположения о том, что st - белый шум, выполнены. Если дополнительно белый шум является гауссовым, то значения Xt распределены нормально, а оценки коэффициентов, произведенные МНК, состоятельны и асимптотически нормальны. Подчеркнем, что нормальность оценок достигается только асимптотически.
Если данные имеют ненулевое выборочное среднее, то можно либо вычесть это среднее из данных и строить регрессию без свободного члена, либо просто строить регрессию со свободным членом в. Другими словами, для оценки математического ожидания процесса AR(p) можно использовать две статистики: уже упомя-
Канторович Г.Г. - профессор, к. физ.-мат. н., проректор ГУ-ВШЭ, зав. кафедрой математической экономики и эконометрики.
— 1 т в
нутую в прошлой лекции X = ^ X Xt и j = i—^-^— , в которой использованы
i=i 1 р
оценки МНК. Для гауссова процесса st обе оценки состоятельны и асимптотически нормальны. Более того, обе они асимптотически независимы от оценок параметров модели, полученных МНК.
Для моделей скользящего среднего невозможно аналитически выразить остаточную сумму квадратов X через значения реализации Xt и параметры модели, а следовательно, применить МНК. Для оценки параметров существуют 2 возможности, по-разному реализованные сегодня в специализированных компьютерных программах.
Первая возможность состоит в применении метода максимального правдоподобия. В предположении нормальности ошибки st выражаем ковариационную матрицу ошибок (ее элементы - значения автокорреляционной функции) через параметры , Р2, ■■■, bq обратимой модели MA(q). Функция правдоподобия для
нормально распределенного вектора X имеет вид:
r 2 I r 1 I _1 X S X х r
L(p,s |X) = —, — ( /detSX ) exp(--^—), где х - отклонения от среднего.
(J2PJ)T
Здесь через SX обозначена ковариационная матрица процесса X. Элементы этой
матрицы выражаются, как мы видели, через параметры модели , Д2, ■.., Pq,
поэтому процедура численной оптимизации позволяет найти оценки метода максимального правдоподобия, которые будут обладать обычными свойствами состоятельности и асимптотической нормальности. Кроме того, оценки параметров модели MA(q) и оценка дисперсии случайного возмущения асимптотически независимы.
Второй подход, позволяющий облегчить вычислительные затраты, применили Бокс и Дженкинс [2]. Они предложили использовать процедуру нелинейной оптимизации: процедуру поиска на сетке (grid-search procedure). Рассмотрим идею этого подхода на примере процесса MA(2). Найдем сначала оценку математического ожидания процесса Х в предположении о его стационарности и эргодичности. Затем выберем некоторые значения (, Р2), например (0,5; 0,5). Исходя из уравнения процесса, Х\ выражается через предыдущие ошибки, но они не известны. Поэтому будем считать, что X1 = X + e1. Это дает нам e1 = X1 _ X . Затем полагаем e2 = X2 _ X _ Д1е1, e3 = X3 _ X _ Д1е2 _ Д2е1 и так далее.
Теперь можно вычислить сумму квадратов отклонений X е2 = S(Р1 ,Д2 ) . Поскольку она посчитана для выбранных значений параметров (, Д2 ) , мы можем рассматривать ее как функцию этих переменных. Далее какой-либо численной процедурой, например поиском на сетке, можно перебирать комбинации p1 и p2
или численно
искать минимум этой функции: min Xe* ■ Коэффициенты (ß1*,ß*),
ß ß 12
которые обеспечивают минимум выражения ^ e2 , и будут оценками коэффициентов модели MA(2).
В общем случае модели MA(q) мы фиксируем начальные значения JJ0
и определяем остатки через наблюдения: e1 = Xj - X , e2 = X2 - X - p{°ex и так далее. Процедура получения остатков достаточно проста, все остатки с отрицательными и нулевым индексами в используемых выражениях заменяем нулями. Затем находим min ^ e2 . Численным поиском находим оценки параметров. Вопросы
J
сходимости этой процедуры и влияния выбора начальных значений на результат мы оставим за рамками рассмотрения. Упомянем только, что показана асимптотическая эквивалентность этой процедуры оцениванию методом максимального правдоподобия.
Аналогичная процедура может быть использована для оценки параметров модели ARMA(p,q). Для их оценивания может применяться комбинация метода наименьших квадратов с поиском на сетке. Чтобы проиллюстрировать, как это делается, рассмотрим модель ARMA(2,2). Пусть уравнение модели имеет вид
(1 - a1L - a2L2)Xt = (1 + JXL + J2L2 )et.
(1 + P,L + J L2 )s.
Его можно переписать в виде Xt =---2———. Введем вспомогательный слу-
1 - a1L - a2 L
чайный процесс Zt =-1-- et. Тогда процессы Xt и Zt связаны соотноше-
1 -a1L -a2L
нием Xt = (1 + PXL + J2L2)Zt, которое напоминает уравнение МА(2), только вместо st стоит Zt. Определим «наблюдения» Zt через наблюдения Xt, сконструировав Zt так же, как ранее остатки модели МА, т.е. значения Zt , которые еще не определены (с нулевыми и отрицательными индексами), полагаем равными нулю. Тогда:
Z\ = X
Z2 = X2 - J 1 Z1
Z 3 = X3 - J1Z 2 - J2 Z1
Например, из выражения Х2 = Z2 + + Ь2%0 следует второе из соотношений.
Далее значения Zt связаны с остатками в( исходной модели следующими соотношениями: 2 = ег. Относительно процесса модель стала АИ(2). Зная «реализацию» Zt для выбранных значений коэффициентов Д >Р2, можно оценить коэффициенты а1 и а2 с помощью обыкновенного МНК. В результате находим оценки коэффициентов (а1 ,а2, Д , Р2) , но получены они по-разному.
Оценки (Д ,Д2 ) заданы как начальные значения. А потом, исходя из них, построены оптимальные оценки (аг,а2). Применяя численный метод оптимизации (например поиск на сетке) оцениваем значения параметров, обеспечивающие минимум остаточной суммы квадратов Xе2 .
Если белый шум гауссов, для оценивания модели ARMA(p,q) можно также применять метод максимального правдоподобия. Общая схема его применения такова. Выражаем значения автокорреляционной функции процесса через параметры модели, формируем ковариационную матрицу порядка T, записываем функцию правдоподобия для имеющейся выборки X1 ,X2,...,XT , решаем (как правило, численно) систему уравнений правдоподобия относительно оценок коэффициентов. Реальные алгоритмы метода максимального правдоподобия, реализованные в специализированных компьютерных пакетах, отличаются от этой схемы, но принципиально она осуществима и отражает логику получения оценок метода максимального правдоподобия.
Напомним процедуру, которую предложили Бокс и Дженкинс. На первом шаге мы определяем, является ли ряд стационарным. Для этого смотрим на поведение выборочной автокорреляционной функции и выборочной частной автокорреляционной функции. Если их поведение не свидетельствует о стационарности исследуемого процесса, преобразуем исходный ряд взятием последовательных разностей и смотрим на поведение тех же функций для преобразованного ряда. Добившись стационарности, мы тем самым определили параметр интеграции d временного ряда. После этого, анализируя для стационарного ряда поведение выборочной автокорреляционной и частной автокорреляционной функций, подбираем параметры p и q. Все это вместе принято называть идентификацией модели. Это некая смесь точных знаний, интуиции и искусства. Затем, после определения (p,d,q), мы тем или иным способом оцениваем коэффициенты соответствующей ARIMA модели.
Прежде чем рассмотреть примеры применения этого подхода, обсудим, как оценить качество полученной модели. Конечно, надо проверить, хорошо ли легла модель на данные, т.е. хорошо ли мы подогнали модель. А также - не нарушены ли нами некоторые предположения о характере случайного возмущения et, что заставит нас сменить, возможно, методы оценивания.
Диагностика модели ARMA
Ясно, что исходной информацией для диагностики служат остатки модели. Мы предполагали, что случайное возмущение является белым шумом, поэтому, прежде всего, надо проверять некоррелированность остатков. Итак, проверке подлежат, в первую очередь, 2 момента. Первый - это качество модели. Нам нужен индикатор типа критерия множественной детерминации R2. И второй - нам надо обязательно проверить некоррелированность остатков, иначе оценивать AR^-часть методом наименьших квадратов нельзя, так как мы получаем несостоятельные оценки.
Как правило, при построении моделей временных рядов критерии качества подгонки моделей применяются для сравнения моделей между собой. Поскольку оценки коэффициентов проводятся путем оптимизации, фактически речь идет о
выборе порядка модели, т.е. о сравнении моделей с различным числом параметров. Абсолютные критерии, типа стандартного коэффициента множественной детерминации R2, не применяются. Из значительного числа критериев познакомимся с двумя наиболее часто применяемыми. Более подробно о критериях подгонки моделей можно прочитать в известной монографии [5].
Наиболее распространенными в настоящее время является предложенный Акаике в 1974 г. [1] критерий AIC (Akaike information criterion) и предложенный Шварцем в 1978 г. [14] BIC (Bayesian information criterion). Оба эти критерия построены примерно одинаковым способом. Информационный критерий AIC для модели ARMA(p,q) выглядит следующим образом:
AIC(p,q) = lns2 + 2,
где Т - число наблюдений. А критерий Шварца имеет несколько другой вид:
BIC(p,q) = lnâ2 + lnT(P+q) . Разумеется, â2 =—^^— . Обратите внимание, как
Т T - p - q
построены эти критерии: логарифмы остаточной суммы квадратов плюс штраф за уменьшение числа степеней свободы. Нам нужно так подобрать значения параметров p и q, чтобы получить минимальное значение каждого из критериев, в то время как по критерию R2 мы отдавали предпочтение модели с большим его значением. Обратите внимание, что в отличие от коэффициента детерминации ничего нельзя сказать ни о диапазоне изменения, ни о знаке критериев Акаике и Шварца, поскольку численная величина оценки ст2 зависит от единиц измерения.
Критерий Акаике базируется на обобщении принципа максимального правдоподобия. Приведенное выражение подразумевает, что случайное возмущение является гауссовым. Шибата [13] показал, что для процессов AR критерий AIC переоценивает порядок модели, следовательно, оценка порядка модели на основании этого критерия несостоятельна.
Критерий Шварца BIC основан на байесовском подходе и имеет более фундаментальное теоретическое обоснование. Показано, что оценка порядка модели по этому критерию является состоятельной. Тем не менее, на практике чаще используется информационный критерий AIC. Существует масса работ, в которых сравнивается применение этих критериев по отношению к разным моделям, но окончательного вывода пока не сделано. На практике разные критерии могут привести к выбору различных моделей. Есть некоторый накопленный опыт по поводу того, к каким типам моделей приводит один критерий и к каким приводит другой. Если вы публикуете какие-то исследования, считается вполне нормальным просто указать, какой критерий вы применяете без специального обоснования вашего выбора. Популярный пакет для анализа временных рядов Econometric views рассчитывает оба эти критерия.
Теперь перейдем к проверке автокорреляции остатков. В этом вопросе существует расхождение между тем, как следовало бы поступать с теоретической точки зрения и тем, что делается практически. Еще в 1970 г. Бокс и Пирс предложили статистический критерий для проверки автокорреляции временного ряда, который сегодня принято называть 0*-статистикой. Тест Бокса-Пирса проверяет гипотезу о совместном равенстве нулю всех автокорреляций временного ряда до порядка m включительно, т.е. гипотезу Но: рх = р2 = ... = pm против альтернативной
гипотезы Н^: X Pi > 0. Бокс и Пирс показали, что при увеличении длины выбор-
i=1
m
ки статистика Q* = T X rk имеет асимптотическое распределение С • Этот тест
k=1
было предложено применять для проверки наличия автокорреляции остатков модели типа ARMA(p,q), при этом число степеней свободы уменьшается на (p+q).
Работа Бокса и Пирса была опубликована в 1970 г. В том же году Дарбин опубликовал в другом журнале свою работу, в которой показал, что статистика Дарбина-Уотсона не применима для проверки автокорреляции остатков при наличии в качестве регрессора объясняемой переменной с лагом. Одновременно Дарбин предложил так называемую h-статистику и альтернативную процедуру Дарби-на для проверки наличия автокорреляции в таких моделях. Эти работы появились в разных журналах примерно в одно и то же время, и каждая породила последователей. Не сразу было замечено, что в моделях с лаговой объясняемой переменной в качестве регрессора статистика Бокса-Пирса не применима по тем же причинам, что и статистика Дарбина-Уотсона. Поскольку в моделях временных рядов проверка наличия автокорреляции критически важна, статистика Бокса-Пирса нашла широкое применение и до сих пор входит в состав большинства специализированных эконометрических пакетов.
Поскольку обнаружилось, что статистика Бокса-Пирса имеет малую мощность, Бокс и Льюнг [6] в 1978 г. предложили использовать для тех же целей
m
улучшенную Q-статистику Q = (T + 2)TX (T _ i)_ r2 , которую, как и статистику
i =1
Бокса-Пирса, часто называют портманто-статистикой (portmanteau-statistics). По сравнению со статистикой Бокса-Пирса различным слагаемым приданы разные веса. Авторы показали, что эта статистика имеет то же асимптотическое распределение %m , но лучше им аппроксимируется при конечном числе наблюдений. К
сожалению, и статистика Бокса-Льюнга теоретически не применима для тестирования автокорреляции остатков в моделях ARMA по тем же причинам, что и статистики Дарбина-Уотсона и Бокса-Пирса. Тем не менее, статистика Бокса-Льюнга входит во все специализированные пакеты, множество исследователей ею пользуются, хотя она теоретически не состоятельна.
Для проверки наличия автокорреляции в моделях ARMA лучше воспользоваться более мощным и универсальным способом, а именно методом множителей Лагранжа (Lagrange multiplier - LM), применительно к проверке автокорреляции остатков его еще называют тестом Бройша-Годфри (Breusch-Godfrey). Он входит в «триаду» классических асимптотических тестов: отношения правдоподобий, Вальда, множителей Лагранжа и применим для широкого класса задач проверки ограничений на коэффициенты модели. Мы знакомились с этим тестом еще в курсе эконометрики, сформулируем его применительно к анализу автокорреляции остатков.
Пусть рассматривается модель множественной регрессии
Y = bo + PX1 +... + PkXk + et, где X1 ,...Xk - разные регрессоры, в том числе, возможно, и лаговые значения как объясняющих, так и объясняемой переменных. Проверим предположение, что et
подчиняется авторегрессионой схеме порядка р, т.е. задается уравнением е = у1£(-1 +... + 7р£г-Р + иг, где щ - белый шум. В терминах коэффициентов модели основная и альтернативная гипотезы принимают вид:
но • 71= ••• = 7_р = 0
Н : 712 + ••• + 7р > 0 .
Метод множителей Лагранжа для проверки этой гипотезы заключается в следующем. Методом МНК строится обычная регрессия вида
Ъ = До + Р X1 +•••+РкХк + е .
Обозначим ее остатки через е^
1. Строится регрессия либо той же объясняемой переменной Уь либо остатков на старые регрессоры и остатки с лагом до р включительно (то есть в качестве дополнительных объясняющих переменных используем et-1 ,et-2,...,е(- ).
2. Проверяется гипотеза о том, что группа дополнительных переменных является излишней. Если в качестве объясняемой переменной используются остатки, то статистика ТЯ2, где Т - число наблюдений, а Я2 - коэффициент множественной детерминации регрессии из пункта 1, имеет асимптотическое (при увеличении числа наблюдений) распределение с2 с р степенями свободы.
Для проверки гипотезы о равенстве нулю группы переменных мы привыкли использовать Е-статистику, проверяющую, фактически, статистическую значимость уменьшения остаточной суммы квадратов от включения дополнительных переменных. Разумеется, Е-статистика применима и в этом случае. Но, напомним, она применима только при нормальном распределении случайного члена. Применение же теста множителей Лагранжа не требует нормальности распределения, но «работает» только асимптотически. Современные эконометрические пакеты, в частности Еу1е'Б, сообщают пользователю при применении ЪМ-теста значения и критические вероятности (р-уа1иеБ) для обеих статистик, так что вы можете выбирать: использовать ли Е-отношение, верное для конечных выборок, но в предположении нормальности, либо не требовать нормальности и использовать статистику ТЯ2 , верную лишь асимптотически.
Для проверки нормальности существует множество тестов. Мы рассмотрим тест Харке-Бера (1ащие-Бета). Этот тест выглядит следующим образом. Он вычисляет выборочные значения для коэффициентов асимметрии £ = —^ V (е( - X)3
Та
и эксцесса К =—^ V(е1 - ц)4 , где X - выборочное среднее, а а - выборочное Та
среднеквадратичное остатков модели. При условии нормальности остатков, стати-
(Т - Р - Ц - 1) 2 1 2 2
стика Харке-Бера - 6-[£ + 4(К - 3) ] имеет с распределение с двумя
степенями свободы.
Лекция 6
Рассмотрим несколько примеров применения подхода Бокса-Дженкинса к временным рядам, представляющим экономические переменные.
Одним из наиболее часто анализируемых финансовых временных рядов является индекс реальной годовой доходности, рассчитываемый компанией Standard & Poors по 500 акциям для США (S&P-500). Рассмотрим ряд для периода 1872-1995 гг. [9]. Среднее значение временного ряда составляет 3,08 американских цента за год.
Выборочные значения автокорреляционной функции rk и соответствующие значения среднеквадратичных ошибок этих оценок приведены в табл. 6.1. В третьем столбце приведены значения статистики Бокса-Льюнга и соответствующие им критические значения вероятностей.
Таблица 6.1.
Номер лага rk s.e.(rk) Q(k)
1 0,043 0,093 0,24 (0,62)
2 -0,169 0,093 3,89 (0,14)
3 0,108 0,093 5,40 (0,14)
4 -0,057 0,094 5,83 (0,21)
5 -0,117 0,094 7,61 (0,18)
6 0,030 0,094 7,73 (0,26)
7 0,096 0,094 8,96 (0,25)
8 -0,076 0,096 9,74 (0,28)
9 -0,000 0,097 9,74 (0,37)
10 0,086 0,097 10,76 (0,38)
11 -0,038 0,099 10,96 (0,45)
12 -0,148 0,099 14,00 (0,30)
Ряд содержит 123 точки. Это довольно большое число, поэтому для проверки значимости можно пользоваться нормальным распределением. Совершенно очевидно, что ни одно из значений rk не значимо на традиционных уровнях значимости. Это видно просто по обычной t-статистике. Отношение оценки к ее сред-неквадратическому отклонению, как правило, даже меньше единицы.
Тот же вывод можно сделать по статистикам Q, все их значения не позволяют отвергнуть гипотезу о нулевых автокорреляциях. Наименьшее из значений критической вероятности p-value составило 0,14. Это значит, что при любом уровне значимости меньше чем 0,14 мы не можем отвергнуть гипотезу об отсутствии автокорреляции, т.е. о совместном равенстве нулю всех рк. Отсюда мы делаем первый вывод, что временной ряд годовых значений индекса S&P-500 является стационарным. Во-вторых, можно идентифицировать модель ARMA. Если бы у ряда была МА-часть, то мы знаем, что рк = 0, начиная с некоторого номера. У
нас уже pj = 0, т.е. МА-части нет. Если бы была AR-часть, то значения рк должны были бы экспоненциально убывать, что не наблюдается. Итак, мы установили,
что наш процесс имеет вид Xt = /л + et, другими словами - это просто белый шум,
или ARMA(0,0). Это давно известный результат. Обычно его трактуют как справедливость игры на бирже, т.е. полную непредсказуемость изменения доходности акций.
Рассмотрим, следуя Т. Миллсу [8], пример другого временного ряда. Объектом рассмотрения является так называемый interest rate spread, т.е. разница между ставками по коротким и длинным бумагам. В качестве длинных бумаг для Великобритании выбраны 20-летние бумаги (UK gilts), а в качестве коротких -бумаги на 91 день (91 day Treasury Bills). Данные являются квартальными, период наблюдения с I квартала 1952 г. по IV квартал 1988 г.
Для этого ряда были рассчитаны 12 значений выборочных автокорреляционной и частной автокорреляционной функций (см. табл. 6.2).
Таблица 6.2.
Номер лага rk s.e. (rk) F kk s-e- ( F kk )
1 0,829 0,082 0,829 0,082
2 0,672 0,126 -0,048 0,082
3 0,547 0,148 0,009 0,082
4 0,435 0,161 -0,034 0,082
5 0,346 0,169 0,005 0,082
6 0,279 0,174 0,012 0,082
7 0,189 0,176 -0,116 0,082
8 0,154 0,178 0,114 0,082
9 0,145 0,179 0,047 0,082
10 0,164 0,180 0,100 0,082
11 0,185 0,181 0,028 0,082
12 0,207 0,182 0,038 0,082
Вплоть до 9-ого лага значения выборочной автокорреляционной функции г^ убывают по величине, а потом они начинают немного расти. Можно сделать заключение, что этот процесс абсолютно точно не является белым шумом. В то же время сомнений в его стационарности нет. Пять первых автокорреляций статистически значимы, остальные нет. Похоже, что ряд может иметь порядок авторегрессии 1 или больше, возможно, увеличение значений функции после 9-ого лага вызвано присутствием пары комплексных корней в характеристическом уравнении.
Для проверки наличия автокорреляций у исходного ряда посмотрим значения портманто-статистик: Q(12)=305,0, а Q*(12)=295,2. Поскольку известно, что математическое ожидание случайной величины, распределенной по ск , равно к, то
даже без таблиц понятно, что, используя Q статистику, мы должны отклонить гипотезу о том, что первые 12 теоретических автокорреляций совместно равны нулю.
Стоит посмотреть на поведение выборочной частной автокорреляционной функции. Единственное статистически значимое значение выборочной частной автокорреляционной функции - первое. Можно сделать вывод, что это достаточно веское свидетельство в пользу модели АИ(1).
Раз ряд стационарен, то d=0. Следовательно, модель, скорее всего, имеет вид ARIMA(1,0,0). Конечно, это не единственная возможность. Давайте оценим эту модель. Модель AR(1) оценивается просто методом наименьших квадратов, что дает следующие результаты (в скобках под оценками коэффициентов - их стандартные ошибки): Xt = 0,176 + 0,856 Xt_1 + et; SS = 0,870 . Мы видим, что interest
( 0,098 ) ( 0,045 )
rate spread уже имеет смысл прогнозировать, так как некоторая информация о его будущем содержится в его текущем значении. Модель содержит значимый свободный член, это означает, что у ряда есть ненулевое математическое ожидание. Его оценка равна ¿и = ^l6, = 1,227, стандартная ошибка этой оценки со-
1 _ 0,856
ставляет 0,506. Q-статистика для ряда остатков этой модели составила 13,21, что очевидно свидетельствует об отсутствии автокорреляций, хотя мы помним, что теоретически этот тест здесь неприменим. Впрочем, тест множителей Лагранжа приводит к тому же выводу.
Наш подход все время строился на том, что мы подбирали наилучшую модель. Есть альтернативный подход. Если мы можем один и тот же процесс записать моделями достаточно близкими по характеристикам, то будем формировать так называемый портфель моделей (набирать группу моделей). В нашем случае в первую очередь имело бы смысл рассмотреть модели ARMA(1,1) и AR(2). После оценивания получим следующие модели:
AR(2): Xt = 0,197+ 0,927 X_ _ 0,079 Xt_2 + et; S = 0,869 .
(0,101) (0,054) (0,084)
Обратите внимание, оценка остаточной суммы квадратов стала чуть меньше, но не очень существенно. Коэффициент при втором регрессоре незначимый.
ARMA(1,1): Xt = 0,213+ 0,831 Xt_1 + et + 0,092 et_l; S = 0,870 .
(0,104) (0,051) 0,095
И в этой модели уменьшения оценки дисперсии ошибки не произошло, хотя мы добавили параметр. Естественно, что мы отдаем предпочтение AR(1). Метод множителей Лагранжа, примененный к каждой из двух моделей, показал отсутствие автокорреляции остатков.
Приведем еще один пример, который показывает, что отбор лучшей модели по критериям AIC и BIC может привести к различным результатам. Рассмотрим ряд месячных значений индекса деловой активности для Великобритании FTA All Share (Financial Times-Actuaries) за период c 1965 по 1990 гг. [8]. Ряд насчитывает около 300 точек, т.е. достаточно длинный. Вначале рассчитаем значения выборочных автокорреляционной и частной автокорреляционной функций (см. табл. 6.3).
Кроме того, была посчитана Q-статистика, Q(12)=19,3. Аргумент 12 указывает количество слагаемых, включаемых в сумму, т.е. количество автокорреляций, совместное равенство которых нулю проверяется. Нулевая гипотеза для этой статистики заключается в том, что вплоть до лага с номером 12 истинные значения рк равны нулю. Эквивалент р-value для Q(12) составил 0,08. Следовательно, на
уровне значимости 5% мы не можем отвергнуть нулевую гипотезу об отсутствии ненулевых значений автокорреляционной функции случайного члена вплоть до порядка, равного 12.
Таблица 6.3.
Количество лагов гк Б.е. (тк) Ф кк 3-е- ( Фкк )
1 0,148 0,057 0,148 0,057
2 -0,061 0,059 0,085 0,057
3 0,117 0,060 0,143 0,057
4 0,067 0,061 0,020 0,057
5 -0,082 0,062 -0,079 0,057
6 0,013 0,062 0,034 0,057
7 0,041 0,063 0,008 0,057
8 -0,011 0,064 0,002 0,057
9 0,087 0,064 0,102 0,057
10 0,021 0,064 -0,030 0,057
11 -0,008 0,064 0,012 0,057
12 0,026 0,064 0,064 0,057
Из табл. 6.3 трудно сразу выявить какое-то определенное поведение гк или Фк, нет убывания значений по абсолютной величине начиная с некоторого
номера. Значения гк для лагов 1 и 3 значимы на 5-процентном уровне значимости, если мы применим критические значения нормального распределения {=1,96. Во всяком случае, трудно увидеть в этих данных свидетельства в пользу точных значений параметров р и д. Пожалуй, имеет смысл рассматривать р и д не больше 3, поскольку значения с большим лагом незначимы. Итак, р < 3,д < 3 . Мы не пытаемся определить точно р и д, а выбираем некие их максимально возможные значения.
Поскольку значения выборочных автокорреляционной и частной автокорреляционной функций быстро становятся незначимыми, есть основания считать ряд
стационарным. В принципе поведение гк и фкк соответствует схеме поведения
стационарного ряда. У нас есть некоторое количество значимых выборочных значений, правда, не видно, чтобы они экспоненциально убывали, но, по крайней мере, стационарные модели с р и д меньше или равным 3 вполне разумные «кандидаты» на более подробное исследование.
Один из популярных подходов заключается в том, что строятся все модели
для р<3 и д<3. Здесь их совсем немного. Поэтому составим сводную таблицу для значений критериев А1С и В1С в зависимости от параметров модели.
Таблица 6.4.
Критерий А1С
Р д
0 1 2 3
0 3,701 3,684 3,685 3,688
1 3,689 3,685 3,694 3,696
2 3,691 3,695 3,704 3,699
3 3,683 3,693 3,698 3,707
Продолжение таблицы
Критерий В1С
р 0 Ч 1 2 3
0 3,701 3,696 3,709 3,724
1 3,701 3,709 3,730 3,744
2 3,715 3,731 3,752 3,759
3 3,719 3,741 3,758 3,779
Как видно, по критерию А1С, минимальный показатель у модели ЛКМЛ(3,0). По критерию Шварца наилучшей оказалась другая модель АЕМА(0,1), т.е. скользящее среднее первого порядка.
Рассмотрим оцененные модели, отобранные по обоим критериям. Для модели АИМА(3,0) было получено следующее регрессионное уравнение (в скобках - среднеквадратические ошибки оценок соответствующих коэффициентов):
X, = 4,41 + 0,173 Xt_l - 0,108 Х,-2 + 0,144 Х,-3 + et; ст = 6,25.
^ (0,62) (0,057) (0,057) (0,057)
Для модели АИМА(0,1) было получено:
Х1 = 5,58 + 0,186 ем + е{; ст = 6,29 .
( 0,35 ) ( 0,056 )
Прежде всего отметим, что все коэффициенты в обеих моделях значимы. Во-вторых, полученные модели выглядят весьма разными. Давайте посмотрим, так ли это. Начнем с оценок математических ожиданий. Для модели АИМА(0,1) оценка математического ожидания, очевидно, равна ¿и = 5,58, а для модели АИМА(3,0) оценка математического ожидания будет равна
4 41
и, =-—-= 5,58.
И 1 - 0,173 + 0,108 - 0,144 '
Видно, что по оценке математического ожидания эти модели эквивалентны.
Так как АИМА(0,1) - это, по сути, есть модель МА(1), и коэффициент 0,186 меньше 1, эта модель является обратимой. Поэтому она эквивалентна следующему выражению:
е =в0 -0,186ХН + (0,186)2Х_2 -(0,186)Х,-3 +. . . = в0 -0,186ХН + 0,035Х,-2 -0,006Х,-3 +. . .
Теперь можно сравнить соответствующие коэффициенты: 0,173 и 0,186 (это первые коэффициенты соответствующих моделей); 0,035 и 0,108 и, наконец, последние - 0,006 и 0,144. Если первые коэффициенты более или менее соответствуют друг другу, то последующие уже не очень. Это означает, что краткосрочная динамика моделей различна. Но если сложить все коэффициенты, то эти суммы оказываются достаточно близкими. Одна равна 0,843, а вторая - 0,791. Сумма коэффициентов модели с распределенными лагами характеризует полное или долгосрочное поведение. Поэтому долгосрочное поведение этих двух моделей не очень отличается друг от друга.
Как уже упоминалось, была предложена другая идея, идея использования портфеля моделей. Т.е. вместо требования выбора одной единственной модели мы выбираем набор, портфель моделей, которыми можно описывать исследуемый временной ряд.
Чтобы сформировать этот портфель, был предложен эвристический подход, т.е. без строгого обоснования, просто из разумных соображений. РобЫИ и Тге-шаупе в 1987 г. [12] предложили рассматривать в качестве меры близости моде-
„ „ {-2т [ АС( рм)-А1С (р,9)]>
лей следующую величину: Я = ехр 2 1 1 . Здесь параметры (р1, д1)
характеризуют наилучшую по рассматриваемому критерию модель. А для всех остальных моделей с параметрами (р, д) мы рассчитываем величину Я. Никакого содержательного смысла в ней нет, это просто мера «расстояния» между моделями. РоБкШ; и Тгешаупе предложили включать модель в портфель, если я <7ш.
Применим этот подход к нашему примеру, посмотрим к каким моделям мы приходим. Критерий А1С отбирает в нашем случае 6 моделей:
{(3,0); (0,1); (0,2); (1,1); (0,3); (1,0)> ,
которые разумно рассматривать.
А вот критерий В1С выбрал всего 3 модели: {(0,1);(1,0);(0,0)>. Любопытно, что по критерию В1С даже модель (0,0) включается в портфель, т.е. даже модель чистого белого шума не так уж плоха. Все модели очень простые (с малым числом коэффициентов), как мы уже отмечали у критерия В1С более серьезные штрафы за потерю степеней свободы.
Отобранные модели можно использовать для прогнозирования и для поиска взаимосвязей между переменными. Это более или менее жесткий отбор «кандидатов».
Лекция 7
Прогнозирование с помощью ARMA моделей
Прогнозирование будущих значений экономической величины является одним из основных способов применения моделей временных рядов. Использование моделей типа ARMA обладает некоторыми особенностями по сравнению с прогнозированием по модели множественной регрессии. При прогнозировании по правильно специфированной модели, вообще говоря, существуют 2 источника ошибок прогноза:
- неопределенность будущих значений случайной величины e;
- отсутствие точных значений коэффициентов модели (у нас есть только их оценки, оцененные по имеющейся выборке).
При прогнозировании по модели множественной регрессии мы оцениваем значение зависимой переменной при заданных значениях независимых переменных - регрессоров. Это приводит к тому, что характеристики прогноза, как случайной величины, являются, по сути, условными характеристиками при условии
имеющейся выборки независимых переменных. Иная ситуация в моделях типа ARIMA. Значение переменной прогнозируется для некоторого будущего момента времени, при этом лаговые значения этой переменной, служащие регрессорами модели, можно рассматривать фиксированными на выборочных значениях, или случайными. Первая возможность приводит к условному прогнозу, как и для модели множественной регрессии, а вторая - к безусловному прогнозу. Таким образом, при прогнозировании по модели типа ARIMA можно рассматривать как условный, так и безусловный прогнозы. Из курса теории вероятностей известно, что условная дисперсия случайной величины не превышает ее безусловную дисперсию, поэтому точность условного прогноза всегда выше.
При прогнозировании по модели ARIMA от имеющейся выборки зависят как оценки коэффициентов модели, так и значения регрессоров, поэтому сложно аналитически выразить условную дисперсию ошибки прогноза через имеющиеся значения временного ряда. Общепринято ограничиваться не очень реалистичным предположением о том, что коэффициенты модели известны точно. Разумеется, это предположение уменьшает дисперсию ошибки прогноза, чем увеличивает кажущуюся точность как условного, так и безусловного прогнозов.
Показано, что если мы хотим достичь минимума среднеквадратической ошибки (MSE), т.е. не требовать несмещенности, то надо взять условное математической ожидание: E{XT+}. Такое условное математическое ожидание
будет гарантировать получение MSE, или иногда ее обозначают MMSE, т.е. Minimum mean square error.
Начнем с модели MA(q): Xt = в + et + ß1et+... + ßqet_q . Мы полагаем, что коэффициенты модели точно известны, и что имеются значения Xt для t е [1, T]. Очевидно, что безусловным точечным прогнозом для любого момента времени будет математическое ожидание процесса, т.е. в. Условным прогнозом для момента времени T+1 будет условное математическое ожидание
Xt+1 = Е{в +
eT+1 + ßieT + ... + ßqeT_q+1 X1,..., XT } .
Среди случайных величин e, которые стоят в левой части, есть такие, которые связаны с имеющимися наблюдениями. Ведь наблюдение складывается из «модельного» значения и ошибки, поэтому условные математические ожидания всех
слагаемых, кроме eT+1, не равны нулю. Рассмотрим сначала E{eT|X1,...XT}, потому что схема одна и та же. Это математическое ожидание - остаток между наблюдением и расчетом, прогнозом по модели, т.е. eT = XT _ XT . Поэтому условные математические ожидания от всех предыдущих значений случайной составляющей надо заменить соответствующими остатками. Точно так же конструируется прогноз не на 1, а на 2 и вообще на h шагов вперед. Все последующие e заменяются нулями, а предыдущие - заменяются реально наблюдаемыми остатками. Таким образом, для модели MA(q) прогноз зависит от того, какие ошибки были на предыдущих шагах. Начиная с шага (q+1) условный прогноз представляет собой просто математическое ожидание в, т.е. условный прогноз совпадает с безусловным.
Рассмотрим условную дисперсию ошибки прогноза на 1 шаг:
уаг0Т+1 - Хт
х„..хт)=Е{({в+8т+1 + Ре + •••+Регч+л -в-Рет + •••+Р/т-^) х^х )2}=<?Е.
Аналогично дисперсия прогноза на 2 шага равна уаг(Хт+2 - Хт+2 X1,•••XT) = (1 + Р12)ст^2,
а дисперсия прогноза на Н шагов составит (1 + Р12 + ••• + Р1)&1 при Н<д. При к > д дисперсия ошибки условного прогноза становится равной дисперсии ошибки безусловного прогноза, т. е. просто дисперсии случайного процесса Х(. Теперь рассмотрим модель стационарного процесса АИ(р):
Х = в+«1Х,-1 + «2 Х - 2 + •••+«рА - Р + е ■
Для прогноза на 1 шаг вперед можно записать:
Хт +1 = Е{Хт
Х1,„Хт } = Е{в+а1 Хт + ••• + арХт - + ет+1 Х1,^Хт } = в+а1 Хт + ••• +арХт-р •
Мы просто подставляем в уравнение модели р предыдущих значений реализации временного ряда. Для прогноза на 2 шага вперед получаем:
Е{Хт + 2 | Х1-Хт } = Е{в+«1Хт+1 + ••• +арХт - р + 2 + Вт + 21Х 1,-Хт } .
Разумеется, математическое ожидание от случайной ошибки е опять даст 0, условные математические ожидания от Х1, Х2 Хт равны самим этим значениям, но в это выражение входит условное математическое ожидание от Хт+1, полученное на предыдущем шаге. Можно подставить его выражение и получить развернутую формулу через значения реализации, но обычно это не нужно. Удобнее рассматривать рекуррентное соотношение, связывающее последовательные значения прогноза. Это соотношение является линейным разностным уравнением порядка р, и его решение стремится при увеличении t к величине
1-в-, т.е. опять к безусловному прогнозу.
Условную дисперсию ошибки прогноза можно рассчитать аналогично случаю модели скользящего среднего, но выкладки становятся весьма громоздкими, даже для моделей малого порядка. Рассмотрим, например, модель АИ(2) без свободного члена. Тогда Хт +1 — «1Хт +«2 Хт-1 , а Хт+, = а,Хт +а2Хт-1 + еТ+, . Очевидно, что дисперсия ошибки прогноза на 1 шаг равна а2е . Для прогноза на 2 шага соответственно получаем:
Хт+2 = а1 Хт+1 +а2Хт =а1(а1 Хт +а2Хт-1) + а2Хт ,
Дисперсия ошибки прогноза на 2 шага равна: (1 + Р12)сте?. Для прогноза на 3 шага получим:
Хт+3 = а1 Хт21 +а2Хт+1 = а1(а1(а1 Хт +а2Хт-1) + а2Хт) +а2(а1 Хт +а2Хт-1) ,
Хт+з =а1 Хт+2 + а2 Хт+1 + ет+3 =
= а1 а (а1 Хт + а2 Хт -1 + Вт+1) + а2 Хт
+ ет+2) + а2 (а1 Хт + а2 Хт-1 + ет+1) + ет+3 •
2 2 2 4 2
Дисперсия ошибки прогноза на 3 шага равна (1 + Р\ + Р2 + 2Д Р2 + Р\ )®е ■
Очевидно, что дисперсия ошибки увеличивается от шага к шагу.
Значительно более простые выражения для дисперсии ошибки прогноза получаются, если перейти от АИ(р) представления модели к эквивалентному МА представлению Х( = в + е( + уе_\ +... + _ч + •••, хотя и с бесконечным числом
слагаемых. Тогда дисперсия ошибки прогноза на к шагов выражается очевидной
к_\
формулой ХУ (Уо = 1) ■
г = 0
Для общей модели АИМА(р^) нужно просто объединить то, что мы сейчас получили. Если мы хотим получить ММБЕ прогноз, то мы рассчитываем прогнозные значения по нашей модели, подставляя туда для времени [1,7] - наблюденные значения Х и рассчитанные значения остатков, а для всех последующих моментов времени - заменяем остатки нулями, а для значений Х подставляем их прогнозные значения. Для получения дисперсии ошибки прогноза переходим к МА представлению и пользуемся только что полученной формулой.
В качестве примера рассмотрим модель ЛИМЛ(1,1) Х( = аХ_ + е( + Ре(_\.
Тогда Х{ = (1 _аЬ) _\(1 + рЬ)ег = (1 + а1 + а21? + ...)(1 + рЬ)ег = ; + _\ + 72е _2 + •••, где
г]\ = а + р , г]2 = а (а + Р) = ак _Ча + Р). Начиная со второго коэффициенты
убывают в геометрической прогрессии. Теперь легко посчитать дисперсию ошибки
2 к_1 2 2п к_1 2г, оч2ч 2П (а + Р)(1 _а2к )ч „ прогноза на к шагов: и; Хч = стг(1 + Ха (а + Р)2) = ст; (1 + -- --). Срав-
г =0 г=0 1 _а
нив с результатом, полученным в лекции 3, видим, что и для этой модели дисперсия ошибки прогноза асимптотически равна дисперсии временного ряда.
Во всех рассмотренных случаях условный точечный прогноз асимптотически приближался к математическому ожиданию ряда, а дисперсия ошибки прогноза - к дисперсии ряда. Это означает, что для стационарного процесса влияние имеющейся информации о реализации на прогноз и его точность асимптотически убывает до нуля. К тому же при увеличении горизонта прогноза дисперсия ошибки не превышает дисперсии временного ряда. Напомним, что этот результат является следствием нереалистического предположения о том, что коэффициенты модели известны точно.
Так же, как и в моделях множественной линейной регрессии, хорошее качество подгонки модели АШМА не гарантирует высокой точности прогноза, т.е. высокой прогнозной силы. Для оценки прогнозных свойств модели одним из общеупотребительных приемов является разбиение имеющейся реализации на 2 части. По первым п наблюдениям выбирается и оценивается модель, а по последним (Т-п) наблюдениям проводится сравнение наблюденных и рассчитанных по модели значений. Такая процедура иногда называется постпрогнозом. Сравнивая прогнозную силу моделей, отобранных по критерию Поскитта-Тремейна, мы выбираем «наилучшую».
Нестационарные временные ряды
Сейчас мы переходим к рассмотрению нестационарных временных рядов. Теперь ситуация изменилась. Из теоремы Вольда следует, что модели типа ARMA охватывают все стационарные процессы. А вот с нестационарными временными рядами ситуация иная, фактически мы будем рассматривать только частные виды нестационарных временных рядов. С одним из таких видов мы уже встречались. По определению модели ARIMA(p,d,q), d - это степень интеграции ряда, т.е. ряд становится стационарным после применения d раз операции взятия последовательной разности. Мы начнем рассмотрение именно с нестационарных рядов, которые могут быть приведены к стационарному виду с помощью взятия последовательных разностей. Естественно, мы начнем с простейшего вида ряда 1(1). Оказалось, что 2 разных по свойствам типа нестационарных рядов приводятся к стационарному виду с помощью взятия последовательных разностей. Мы их уже рассматривали.
1 тип: процесс с детерминированным полиномиальным трендом.
Xt = Pk (t), где Pk (t) - полином степени k от t, а st - стационарный процесс, не обязательно белый шум. Если ограничиться рассмотрением только линейного тренда Xt = a + pt + et, то можно записать: DXt = Xt - Xt-1 = (1 - L)Xt = p + (st - et_1 ).
Поскольку et - стационарный процесс, то его первая разность также стационарный процесс, хотя если st - белый шум, то появляется МА-часть. В случае полиномиального тренда для приведения к стационарному виду нужно взять последовательную разность несколько раз.
2 тип: процесс случайного блуждания: Xt = m + Xt-1 + st. В этом случае
DXt = m + st, и процесс Xt называется случайным блужданием с дрейфом. Мы можем записать решение этого разностного уравнения в следующем виде:
t-1
Xt = mt + Zet- j ■
j=0
Поэтому, если применить подход Бокса-Дженкинса и, имея некоторую реализацию, перейти к разностям, оценить модель ARMA(p,q), то чтобы вернуться назад к процессу Xt, нужно выбрать, по какой схеме возвращаться. Как мы уже рассматривали, эти процессы ведут себя по-разному. В чем-то они схожи: у обоих есть линейный тренд, но они отличаются случайной частью. В первом случайная часть - это текущий шок, текущие возмущения, а во втором - это накопленные возмущения от всех предыдущих шоков.
Если случайное блуждание можно привести к стационарному виду только взятием первой разности, то ряд первого типа можно привести к стационарному виду также выделением линейного тренда, например построив линейную регрессию на время и рассмотрев стационарный остаток. Применим ли такой подход к случайному блужданию? Что произойдет, если на самом деле процесс является случайным блужданием, пусть даже для простоты с нулевым математическим ожиданием, т.е. имеет вид DXt = st, а мы строим регрессию на время вида
Xt = a + pt + et ? Можно интуитивно догадаться, что мы получим значимый по
t-статистике коэффициент р. Метод наименьших квадратов как бы преобразует непостоянную дисперсию в значимый тренд, т.е. в непостоянное математическое ожидание. Таким образом, мы эти 2 типа процесса спутаем. Вопрос о том, насколько опасно путать эти 2 процесса между собой, привлек внимание относительно поздно. Для краткости введем следующие названия для этих типов нестационарных процессов.
1 тип: процесс, приводимый к стационарному путем выделения линейного тренда - TSP (trend stationary process). Это процесс вида Xt = a + р t + et, он приводится к стационарному процессу путем включения в регрессию линейного тренда. Это, в принципе, процесс, у которого есть детерминированный тренд. Иногда такой процесс называют TS.
2 тип: процесс, приводимый к стационарному путем взятия первой разности - DSP (diferencing stationary process). Вид этого процесса таков: Xt = Xt_ + st.
Иногда такой процесс называют DS.
Сведем же рассмотренные нами различия в свойствах этих процессов в табл. 7.1.
Таблица 7.1.
Процесс ТБР Процесс DSP
• Не стационарен из-за непостоянного тренда. • Конечная память о шоках, он забывает об ошибке на предыдущем шаге сразу. Если вместо белого шума будет стоять более общий процесс АИМА(р^), то, конечно, шоки сказываются некоторое время, но их влияние со временем ослабевает. • Не стационарен из-за непостоянной дисперсии. • Так как в явном решении стоит сумма всех предыдущих e, то шоки помнятся все время. Это процесс с бесконечной памятью. Экономически это не очень понятно, шоки не должны сказываться постоянно.
В 1984 г. была опубликована работа Нельсона и Канга [10], в которой они задались вопросом: что произойдет, если процесс является DSP, а мы будем выделять линейный тренд, как для процесса типа TSP? Их результаты были получены с помощью большого количества испытаний Монте-Карло и заключались в следующем.
1. Линейная регрессия случайного блуждания без дрейфа на время дает коэффициент множественной детерминации R2 » 0,44 независимо от размера выборки, т.е. R2 значим даже для малого количества точек наблюдения.
2. Если дрейф присутствует (m* 0 ), то R2 > 0,44 и зависит от размера выборки. Это интуитивно можно понять, потому что как только m* 0 , появляется
свой тренд - нестационарное среднее. Также оказалось, что R2 ® 1 при T ® ¥ ,
т.е. для процесса DSP с дрейфом R2 ® 1 при увеличении объема выборки. Это явление значимого коэффициента тренда, когда в истинном уравнении тренд отсутствует, было названо «кажущиеся тренды» (spurious trends).
3. Оценка дисперсии остатков составляет примерно 14% от истинной дисперсии случайного возмущения, т.е. оценка дисперсии сильно занижена. Процедура оценивания указывает на значимый тренд и малую дисперсию, хотя на самом
деле процесс случайного блуждания характеризуется растущей без ограничений дисперсией.
4. Остатки регрессии оказываются коррелированными с коэффициентом
1 10
корреляции примерно равным р1 = 1 _ т ■
5. Разумеется, в этом случае t-статистика не годится для проверки гипотезы о значимости коэффициента при времени, она смещена в сторону принятия гипотезы о наличии линейного тренда.
6. Независимые случайные блуждания демонстрируют высокую корреляционную зависимость.
Этот пункт требует пояснения. Если взять 2 независимых случайных блуждания:
Xt = Xt _1 + et'
где ut и st - независимые белые шумы, и построить регрессию Yt = a + bXt + ut, то
из общих соображений мы ожидаем, что коэффициент b будет незначим. Но вы должны чувствовать из первых четырех результатов, что этого здесь не произойдет. Два независимых случайных блуждания показывают высокую регрессионную зависимость, коэффициент b оказывается значимым. Как можно это пояснить? Если оба процесса имеют значимые тренды, то они оба зависят от t. Регрессионный анализ покажет, что процессы зависимы между собой. А это означает, что если каждая из величин, между которыми мы ищем регрессионную зависимость, является DS-процессом, то регрессия между ними является «кажущейся». Таким образом, нею-торые зависимости являются кажущимися. Нам представляется, что связана динамика денежной массы и инфляции, но мы не учли, что оба процесса - DSP, и сделанный экономический вывод неправомочен. Этот результат Нельсона и Канга показывает опасность прямого применения обычного регрессионного анализа в том случае, когда данные являются типа DS.
Значимость этого результата особо подчеркнута более ранней, знаменитой работой Нельсона и Плоссера [11], которая связана с исследованием исторических макроэкономических рядов в США. Результат этого исследования показал, что практически все ряды макроэкономических показателей США, за исключением ряда уровня безработицы, оказались рядами типа DS, т.е. типа случайного блуждания с дрейфом. Это произвело впечатление разорвавшейся бомбы. Известный экономист Саржент заметил, что все, что сделано до сих пор в области макроэкономической динамики, подлежит пересмотру. Правда, дальнейшие исследования показали, что ситуация не столь драматична.
Две эти работы, Нельсона-Канга и Нельсона-Плоссера, поставили вопрос ребром. Оказывается, надо проводить четкое разграничение между рядами двух видов: DS и TS. Если тренд в модели типа TS отсутствует, то вопрос сводится к различию между стационарными и нестационарными рядами. Раньше мы рассматривали для этого поведение выборочной автокорреляционной функции.
В прошлой лекции мы рассматривали ряд значений спрэда в Великобритании и пришли к выводу, что он описывается моделью AR(1). Если взять первые разности этого ряда, то он будет также хорошо описываться моделью ARIMA(0,1,0). А это очень разные ряды. Один из них (TS с нулевым коэффициентом при трен-
де) является стационарным, а второй - нестационарным. Визуальное же исследование поведения выборочной автокорреляционной функции, т.е. применение подхода Бокса-Дженкинса не позволяет надежно выбрать, к какому из типов относится исследуемый ряд.
Мы нуждаемся в методе, позволяющем формально проводить различие между рядами типа TSP и DSP.
Эти 2 ряда описываются следующими моделями:
TSP
DSP
Xt = a + p + et Xt = m + Xt _ + et
Попробуем рассмотреть модель: Х( = а + рХ(-1 + pt + е(. Она вобрала черты
обеих моделей, и гипотезы о характере ряда можно записать в виде простых гипотез о ее параметрах.
Г=1
Р = 0-
Но: ряд является DS ^
Н1: ряд является ТБ ^ |р| < 1 (тогда е будет не просто белый шум, а некоторый стационарный ряд. Вообще-то мы это имели в виду, когда определяли ТБ, предполагая, что он приводится к стационарному с помощью выделения тренда).
Нулевая гипотеза относится к классу общих линейных гипотез. При традиционном подходе для ее проверки нужно оценить 2 регрессии: Х: = а + pXt-1 + pt + и
Х: = а + Х:-1 + . И затем проверить значимость разности остаточных сумм квадратов, используя F-статистику■
Тест, к рассмотрению которого мы сейчас переходим, появился сначала в более простой версии этой модели: Х: = а + pXt-1 + , т.е. без включения линейного тренда. В этом случае гипотезы принимают вид:
Н0: р = 1 ^ ББ Н1: р< 1 ^ ТБ.
Теперь можно проверить гипотезу о том, что р = 1 при помощи {-статистики. Уравнение можно переписать в другом виде. После вычитания из обеих частей Х-1, получим АХ: = а + (р-1)Х:-1 + . Пусть р-1 = у , тогда проверяемые гипотезы примут вид:
Н0: у = 0 Н1: у < 0.
В классической линейной регрессии для проверки такой гипотезы применяется односторонняя {-статистика. Но в случае выполнения нулевой гипотезы, ряд Х: является случайным блужданием, его дисперсия стремится к бесконечности
у
при увеличении t, и распределение статистики —не является распределением
ем Стьюдента. Следовательно, и асимптотическое распределение этой статистики не является нормальным. Причиной является невыполнение условий Центральной предельной теоремы в этом случае. Аналитическое выражение для асимптотиче-
У
ского распределения статистики —можно выразить через стохастические ин-
se.(g)
тегралы от винеровского случайного процесса, но мы этого делать не будем. Критические точки этого распределения приходится рассчитывать численно, используя симуляционные процедуры Монте-Карло. Впервые это распределение было выведено и затабулировано в работе Дикки и Фуллера [3] и носит их имя. Тест, использующий для проверки типа нестационарности это распределение, при условии у = 0 , т.е. когда процесс принадлежит типу DS, называется тестом Дикки-Фуллера, и обозначается как DF-тест. При условии, что нулевая гипотеза о том, что у = 0, выполнена, мы имеем процесс типа случайного блуждания (DSP). Именно для этого случая не работает t-статистика.
Позже МасКтпоп [7] расширил эти таблицы, рассмотрел некоторые другие случаи и предложил аппроксимирующие формулы для быстрого расчета критических точек в компьютерных программах. С помощью моделирования Дикки и Фуллер также рассчитали аналог F-статистики для случая, когда процесс является случайным блужданием. Это распределение также иногда называют распределением Дикки-Фуллера. Обычно из контекста ясно, о каком из распределений идет речь.
Вернемся к общей модели: Xt = a + pt + (р_ 1)Xt_1 + st. В этом случае мы имеем следующие гипотезы:
Н>:р = 1; p = 0 Н :р< 1.
Сравним, как ведет себя DF-статистика и F-статистика. Нас интересует критическое значение при заданом числе наблюдений и для заданного уровня значимости. Если Т - количество наблюдений, то в зависимости от него получим следующую табл. 7.2 для уровня значимости 5%, в которой F-статистика - это фактически F2 тз
Таблица 7.2.
T
DF
F-статистика
25 50 100
¥
7.24 6,73 6,49
6.25
3,42 3,20 3,10
3,00
На рис. 7.1 схематически показаны графики плотностей распределения Дикки-Фулера и Фишера. Мы видим, что ЮЕ-распределение значительно правее, чем Е-распределение. Это означает, что в большом количестве случаев применение стандартной F-статистики ведет к тому, что мы считаем, что ряд относится к типу ТБ, в то время как он типа ВБ, а именно, если F-статистика попадает в область (а; в).
F-распределение
DF-распределение
Рис. 7.1.
Кроме того, оказалось, что аналитическое выражение, а следовательно, и форма БЕ-распределения зависят от того, включены ли в модель регрессоры at, Р1 или нет ни того, ни другого. Поэтому для этих случаев нужны, по сути, разные таблицы.
Прежде чем двигаться дальше, давайте рассмотрим пример оценки конкретного ряда: логарифма индекса производства по данным Федеральной Резервной Системы США [4]. Рассматриваются данные с I квартала 1950 г. по IV квартал 1977 г. Специфицируем модель в следующем виде:
У = А + Ь •t +а1 • Yt_1 +а2(^-1 - ^2) + Ъ ■
Естественно, мы предполагаем, что ~ ШЫ(0, ст ), т.е. белый шум.
Эта модель была оценена методом наименьших квадратов, но перед оцениванием из обеих частей вычли . Это стандартный прием, тогда можно сравнивать коэффициент с нулем, а не с единицей. А в левой части уравнения оказывается объясняемая переменная ДУ . Поэтому были оценены 2 следующие модели, соответствующие нулевой и альтернативной гипотезе (в скобках - стандартные ошибки оценок коэффициентов).
Первая модель: у - у-1 = ду = 0,52 + 0,0012 t - 0,119 у-1 + 0,498(У-1 - у-2);
з.е. (0,13) (0,00034) (0,033) (0,081)
ВББ = 0,056448 (ВББ ж ) .
Вторая модель: Пусть
У, - У,-1 = 0,0054 + 0,447 (у-1 - У,-2).
(0,0025) (0,083)
А = 0; а1 = 1. МНК дает следующую модель: Соответственно КББ = 0,063211 (ЛХХЛ). Считаем
(0,063211 - 0,056448)/
стандартное F-отношение: Ес
отношение
0,056448
■ = 6,34 . А теперь по-
106
смотрите на приведенную выше табл. 7.2. Ряд содержит немного более 100 точек, поэтому воспользуемся строкой для 100 наблюдений. Если бы мы сравнивали F-отношение с квантилями распределения Фишера, то Екритическое для 5-процентного уровня значимости было бы равно 3,1. Соответственно мы бы пришли к выводу, что исследуемый ряд - типа ТБ. Но мы должны сравнивать статистику с
квантилями распределения Дикки-Фуллера, т.е. с 6,49. В результате на уровне значимости 5% гипотеза о наличии единичного корня не может быть отвергнута, а для нас это означает, что ряд относится к типу DS. Поскольку проверка типа ряда - TS или DS - сводится к тестированию наличия или отсутствия единичного корня, то эта процедура носит стандартное название unit root test. Речь идет о проверке наличия единичных корней в характеристическом уравнении. Но наличие или отсутствие единичного корня - это то, что в подходе Бокса-Дженкинса называется порядком интеграции - единица или нуль. Т.е. unit root test, который появился как возможность различения между TS и DS рядами, дает статистическую процедуру для определения порядка интеграции ряда в подходе Бокса-Дженкинса.
Наибольшее распространение получило использование распределения Дик-
9
ки-Фуллера для статистики —. Мы подробно рассмотрим ее применение в
s.e.(r)
следующих лекциях.
* * *
СПИСОК ЛИТЕРАТУРЫ
1. Akaike H. A New Look at the Statistical Model Identification, IEEE Transactions on Automatic Control, AC-19. 1974. P. 716-723.
2. Box G. E. P. and Jenkins G. M. Time Series Analysis, Forecasting and Control, rev. Ed. San Francisco: Holden-Day, 1976.
3. Dickey D. A. and Fuller W. A. Distribution of the Estimators for Autoregressive Time-Series with a Unit Root // Journal of the American Statistical Assiciation. 1979. Vol. 74. P. 427-431.
4. Dickey D. A. and Fuller W. A. Likelihood Ratio Statistics for Autoregressive Time Series with a Unit Root // Econometrica. 1981. Vol. 49. P. 1057-1072.
5. Judge G. G., Griffits W. E., Hill R. C, Lutkepohl H, Lee Tsoung-Chao. The Theory and Practice of Econometrics. Second edition. NY: John Willey and Sons, 1985.
6. Ljung G. M. and Box G. E. P. On a Measure of Lack of Fit in Time Series Models // Biometrika. 1978. 65. P. 297-303.
7. MacKinnon J. C. Critical Values for Cointegration Tests // UC San Diego Discussion Paper. 1990. P. 90-4.
8. Mills T. C. The Econometric Modelling of Financial Time Series. Cambridge University Press, 1993.
9. Mills T. C. The Econometric Modelling of Financial Time Series. Second edition. Cambridge University Press, 1999.
10. Nelson C. R. and. Kang H. Pitfalls in the Use of Time as an Explanatory Variable in Regression // Journal of Business and Economic Statistics. January 1984. Vol. 2. P. 73-82.
11. Nelson C. R. and Plosser C. I. Trends and Random Walks in Macroeconomic Time Series: Some Evidence and Implication // Journal of Monetary Economics. 1982. Vol. 10. P. 139-162.
12. Poskitt D. S. and Tremayne A. R. Determining a Portfolio of Linear Time Series Models // Biometrica. 1987. V. 74. P. 125-37.
13. Schibata R. Selection of the Order on an Autoregressive Model by Akaike's Information Criterion // Biometrika. 1976. 63. P. 147-164
14. Schwarz G. Estimating the Dimension of a Model // The Annals of Statistics. 1978. 6. P. 461-464.
