CC BY
78
ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2012 Управление, вычислительная техника и информатика № 1(18)
УДК 519.2
Ю.Б. Буркатовская, С.Э. Воробейчиков, Е.Е. Сергеева ОЦЕНИВАНИЕ ПАРАМЕТРОВ И ОБНАРУЖЕНИЕ МОМЕНТА ИХ ИЗМЕНЕНИЯ ДЛЯ ОБОБЩЕННОГО АВТОРЕГРЕССИОННОГО ПРОЦЕССА С УСЛОВНОЙ НЕОДНОРОДНОСТЬЮ
Рассматривается задача обнаружения момента изменения параметров GARCH(p,q)-процесса. Авторегрессионные параметры процесса предполагаются неизвестными до и после момента разладки. Предлагается последовательная процедура оценивания параметров, основанная на взвешенном методе наименьших квадратов. Выбор весовых коэффициентов и момента остановки позволяет строить оценки с ограниченным среднеквадратическим отклонением, зависящим от параметра процедуры H. Процедура определения момента разладки основана на сравнении оценок неизвестных параметров процесса на различных интервалах наблюдения. Получены верхние границы для расчета вероятностных характеристик предложенной процедуры: вероятности ложной тревоги и ложного спокойствия.
Ключевые слова: GARCH(p,q), момент разладки, метод наименьших квадратов, среднеквадратическое отклонение, гарантированное оценивание.
В 1986 г. T. Bollerslev [1] впервые предложил использовать модель GARCH -обобщенную авторегрессионную модель гетероскедастичности для анализа временных рядов. При описании процессов типа GARCH предполагается, что на текущую изменчивость дисперсии влияют как предыдущие изменения показателей, так и предыдущие значения дисперсии. Модели типа GARCH часто используются при обработке информации в задачах последовательного анализа данных, имеющих эконометрическую направленность, а именно, при управлении финансовыми рисками, так как пренебрежение определением структурных изменений может приводить к финансовым потерям.
В настоящее время интерес к данной модели не снижается, о чем свидетельствует большое количество работ в этой области [2-5] и др. Так, в работе [2] E. Hillebrand предложил алгоритм оценивания, основанный на функции логарифмического правдоподобия, с использованием ненаблюдаемого процесса условной вариации. Davies и др. [3] для определения изменения применяют обобщенное отношение правдоподобия, которое приводит к квадратичной форме. Gombey и Serban в [4] используют эффективный вектор вклада в последовательной процедуре, когда необходимо определить изменение параметра, если начальное значение задано, а остальные компоненты являются мешающими параметрами. Работа Е. Gombey [5] посвящена апостериорному методу обнаружения разладки, когда полностью доступна последовательность наблюдений и не определены начальные параметры, однако необходимо оценивать все параметры модели. Для определения разладки автор использует функцию log-правдоподобия и вектор эффективного вклада.
Во многих практических приложениях задача обнаружения момента разладки случайных процессов оказывается тесно связана с задачей оценивания параметров этих процессов. Для оценивания параметров GARCH-модели часто используется
оценка квазимаксимального правдоподобия. Таким оценкам, в частности, посвящены статьи [6-8]. Pan и др. [9] изучали вероятностные и асимптотические свойства оценки квазимаксимального правдоподобия параметров пороговой модели GARCH. Для стандартной модели GARCH асимптотические свойства, в частности асимптотическую нормальность, оценок такого типа рассматривали Berkes [7] и Francq и Zakoian [8]. Straumann и Mikosch [10] установили, что оценки квазимак-симального правдоподобия для общего класса моделей с условной гетероскеда-стичностью являются асимптотически нормальными. В работе [11] рассматривался класс пороговых GARCH-моделей и строились оценки по методу квазимаксимального правдоподобия и по методу наименьших квадратов. Авторы показали, что асимптотически МНК-оценка параметров является более точной, чем оценка максимального правдоподобия. Робастные оценки параметров модели GARCH рассматривались в работах [12-14]. Baillie and Chung [15] предложили оценку с минимальным расстоянием для модели GARCH(1,1), которая основывается на автокорреляционной функции квадратов наблюдений. В статье [16] предлагается метод оценивания, основанный на автоковариационной функции квадратов наблюдений, не требующий знаний о функции распределения.
В данной работе предлагается последовательная процедура обнаружения момента разладки процесса GARCH(p,q) с неизвестными авторегрессионными параметрами. Метод обнаружения разладки для случая с известными параметрами рассмотрен в [17]. Предложен метод оценки неизвестных параметров процесса, основанный на модифицированном взвешенном методе наименьших квадратов [18], позволяющий получить оценки с гарантированной точностью.
1. Постановка задачи
Рассматривается устойчивый случайный процесс GARCH(p,q)
Xn = CTnSn , n = 0,1... (1)
где {en} - последовательность независимых одинаково распределенных случайных величин с нулевым средним, единичной дисперсией и известным распределением. Плотность распределения шумов положительна на всей числовой прямой. Условная вариация процесса xn представляет собой случайный процесс вида
ст2 = a + V Xix2 , +V ц,-ст2 ,.
П I П—1^^^ >1 n—i
i =1 i=1
Параметры {a, X,} предполагаются неизвестными, а параметры ц, - известными. Параметры процесса удовлетворяют условиям
p q
a > 0, X, > 0, ц j > 0,1 < i < p,1 < j < q,0 < ^ X, + ^ ц j <1. (2)
i=1 j=1
Тогда процесс {ст;} является стационарным [1]. Вектор параметров Л = [a, Xj,..., X p ] меняет свое значение с Л0 на Л1 в неизвестный момент времени
0. Начальное и конечное значения параметров удовлетворяют условию
||Л 0 — М2 > Д >0,
где Д является известным значением, определяющим минимальное расстояние между значениями параметров до и после момента разладки. Требуется по наблюдениям за процессом {xn} определить момент разладки.
2. Построение оценок параметров
Так как значения параметров до и после момента разладки являются неизвестными, то необходимо получить их оценки. Преобразуем сначала процесс (1), записав его в матричной форме
= XnЛт + MSn-1, Xn : qх(р +1), M : q х q;
Sn = [стП—l,•••,a2-q]T , Л = 1 1 , V ., v p ];
Л 2 2 1 Xn—1 . Xn—p "ц . • М q—1 М q
О О О , M = 1 . 0 0
0 0 ... 0 _ 0 . 1 0
X„ =
Используя это представление, получим
Sn = Xn Лт + MS„_1 = Xn Лт + M (X„-1AT + MSn_ 2 ) = ( + MXn_!)T + M 2 Sn_ 2 = ...
... = ( + (Xn-1 +...+MkXn-k ) Лт + Mk+1Sn_k _1. (3)
Из ограничений на параметры (2) следует, что для любого /-мерного вектора S с положительными коэффициентами выполнено неравенство
(q Л
max (MqS), <1£ц/ I max (S)..
1<]<q' 11
■ i 1 1< 1< q' '1 V1=1 / 1 4
Следовательно, тах(Мк+'£) ^ 0 при к , и последним слагаемым в выра-
1< ]< д' ']
жении (3) можно пренебречь.
Запишем аппроксимацию для вектора 8п более компактно:
Бп = О (п, ц, х )ЛТ.
Здесь О (п, ц, х) - случайная матрица, зависящая от реализации процесса {хк }1<к <п _1 и подчиняющаяся следующим рекуррентным соотношениям:
О (п, ц, х) = Хп + МО (п _ 1, ц, х);
О (no, Ц, х) = Xn0, по = тах { Р}.
Условная дисперсия процесса (1) ст^ - это первая компонента вектора 8п. Обозначив первую строку матрицы О (п, ц, х) через ¥ (п, ц, х), получаем
= F (n, ц, х )Л
т
¥(п,ц,х) = [¥0 (п, ц,х),¥х (п,ц,х),...,¥р (п, ц,х)].
Рекуррентные уравнения для компонент вектора ¥ (п, ц, х) принимают вид
д д
¥о (пцх) =1+Xц}¥о (п_j,цх), ¥ (пцх) = х2_ + Xц¥ (п_ 1,ц,х), 1 = 1,Р
1=1 1=1
¥ (по _ i, ц, х) = [1, хпо _1_ 1, • • •, хпо _ р_ ] ], 1 =0, д _1.
(4)
Чтобы привести процесс (1) к удобному для исследования авторегрессионному процессу, используем подход, который был предложен в [19]. Учитывая пред-ставлиние (4), запишем процесс (1) в виде
х2 = ст2 + стп(е2 _1) = ¥ (п цх )лТ + ¥ (п цх )лТ ( _1).
Введем следующие обозначения:
хп _ ¥ (п, ц, х)
", «п.г _
, 0 < i < p,
Un = [^n,0 ,
J, Л = [а,V-,Xp], Zn = ^A B2 = E( -l) bn = ВЛПп.
,
п, р J у |_ у 1у у р ^ ^ п в
Теперь представим процесс (1) в векторном виде:
уп _ Л^ + ЬСп. (5)
Построим оценку величины дисперсии шума Ьп в (3). Поскольку 0 < «п,. < 1 ,
(
bn = БЛи„ < B
\
a + ХЛ V i=i У
(6)
Отсюда следует, что дисперсия шума в процессе (5) ограничена сверху, но неизвестна.
Для компенсации неизвестной дисперсии помех вычисляется статистика Г^ . Учитывая (4), имеем
( р 1
х2 ^ т1п{ (n,ц,х)-,¥р (n,ц,х) « + Хх. 822
V г=1 /
Введем обозначение
2
min
Множитель ГN можно выбрать в одной из двух форм:
а) Г N = CN
X z
V1=N0 У
Cn = B2 E
V2
X^n
Vl=N0 У
(7)
b) rN = Cn X z2, Cn = B2E X
i=N0
Vl=N0 У
Плотность распределения шумов {el} должна быть такова, чтобы существовала константа Cn в (7). Типы таких плотностей рассмотрены в [19]. Промежуток [ N0, N ] выбирается из следующих соображений: значение N0 должно быть достаточно велико, чтобы использование аппроксимации (4) было корректным; кроме того, значения функций min {F0 (n, ц, х),..., Fp (n, ц, х)} не должны быть близки
к нулю. Поскольку функция F0 (n, ц, х) возрастает с ростом n , а плотность распределения шумов процесса положительна на всей числовой прямой, то такой промежуток может быть выбран за конечное время.
Zn =
Из соотношений (7) следует, что
( р 12
Е— < —, С = В2 Гп С
а + ХЛ;
V г=1 У
(8)
Для нахождения оценки неизвестного вектора параметров Л используем последовательный метод оценивания, предложенный в [18]. Оценка имеет вид
( N 1
Л (N1) =
X Упу(п х)итп А 1(N1),
V n_N+1 у
N1
А(^)_ X Кп,х)Ппитп .
п_ N+1
Для фиксированного значения Н>0 определим момент остановки т _ т(Н
т _ т(Н) _ > N +1: (N1) > н}, (9)
где итт(^) - минимальное собственное значение матрицы А( N1).
Положительные весовые функции у(п, х) на интервале [ N +1, N + ст_ 1] задаются следующим образом:
у(п, х) _—:= 1 , (10)
7г N<0,,
где с - наименьшее значение N , для которого матрица А(^) не вырождена. Веса у(п, х) на интервале [N + ст, т_1] находятся из условий
^Г^ = X ^(п,х)иПип. (11)
N n_N+ст
Последняя весовая функция у(т, х) находится из условий
^Г^> 'X ^(п,х)иПип, От1П(т)_Н. (12)
Г N п _ N+ст
Из соотношенй (10) - (12) следует, что весовые функции удовлетворяют условию
1
0 < у(п, х) < - --
а/г^
Отсюда и из (6), (8) получаем, что
Еу1 (п, х)Ь2 < 1.
Оценка параметров Л*(Н) в момент времени т имеет следующий вид:
Л*(Н)_ (X у(п, х)упиП | А_1(т), А(т) _ X у(п, х)и„иП. (13)
V п_ N+1 У п_ N+1
Свойства предложенной оценки определяет теорема 1.
Теорема 1. Для любого значения параметра процедуры Н>0 момент прекращения наблюдений т(Н) конечен с вероятностью единица и средний квадрат
*
нормы отклонения оценки Л (Н) от истинного значения вектора параметров Л удовлетворяет неравенству
е\|л*(Н)-Л2 <Н±Р.
II II Н2
Доказательство. Момент прекращения наблюдений т(Н) является конечным тогда и только тогда, когда расходится почти наверное ряд
ГN Ё у2(П’ X)UTUn = «•
n=N+1
Для сходимости почти наверное данного ряда необходимо и достаточно, чтобы выполнялось условие [20]
Ve> 0: P ІГ n Ё v2(n, x)UTnUn >є
n=k
k
->0.
(14)
Так как UTnUn > 1, это условие может выполняться только при Г NV2 (n, x) ^ 0 по вероятности. Перепишем уравнение для определения весовых коэффициентов (11) в виде [21]
Umm (n) = Г- min Sn (a (n -1) + v(n, x)UnUTn ) S )
Umrn(« - + v2(n, x)UTnUn •
Г„ Г^х|И|=1 V — " г„
Отсюда получаем, что для произвольного вектора 5 : 1|£|| = 1 верно
V2 (n, х)и Тпипг N - х) (итп8 )2 - (а (п -1) 5 - от1П (и -1)) < 0.
Приравнивая левую часть к нулю и решая квадратное уравнение, получаем, что оно имеет два корня, один из которых не положителен, а второй - не отрицателен. Коэффициент v(и, х) удовлетворяет условию
v(n, x) <
КS)2 +yj(UUS)4 + 4rNUUUn (SnA (n - 1)S - omin (n -1)j
2Г NUUUn
Равенство здесь достигается, когда £ является собственным вектором матрицы А (п), соответствующим ее минимальному собственному значению. Отсюда получаем оценку члена ряда в (14)
ГNv (n, x)Un Un > min
N n n Sil S||=1
(uUs )4
2Г NUUUn
+ (SnA (n - 1)S — Umin (n -1)). (15)
Правая часть этого неравенства сходится к нулю, если и только если оба ее слагаемых сходятся к нулю. Для этого вектор 5 должен сходиться к вектору, соответствующему минимальному собственному значению матрицы А(и-1), и одновременно к вектору, ортогональному ип. Поскольку вектор ип зависит от случайной величины ^п, не зависящей от А(и-1) , правая часть (15) будет больше некоторой константы с положительной вероятностью.
Используя неравенство Коши - Буняковского и соотношение (4), получаем
ЕЛ*( H) -Л < E
Ё bnZnv(n x)UU I A~1(t)
n=N+1
Ё bnZnv(^x)UU
n=N+1
(16)
Введем усеченный момент остановки t(N) = min{т,N1}, тогда t(N)^т при N1 ^ да . Выражение (16) преобразуется к виду, отличающемуся от исходного только пределом суммирования. Обозначим через Fn = ст {,..., х С«0+1, —, Zn } ° - алгебру, порожденную случайными величинами {,...,х ,Z^ +1,---,Zn}. Используя свойства условного математического ожидания получаем
-E
H2
£ Ь2пv2(n,x)UTnUnZХхт)
т( лу
£ ЬпZnV(n, x)Un = -^ EE
п= N+1
N1 n-1
£ £ bAv(k>x)v(n x)UlUnQkCnxт)
n-1
+A ee H2
=N+2 k = N+1 N,
n-1
H
N
1 ^V1
-E £ bn2v2(n, x)UTUnX(näX)E[Zn|¥n-1 ]-
n=N+1 n-1
+—7Е £ £ x)v(n, с^Х(и<т)Е[с„1¥„-1 ] =
7 п =Ы + 2 к =Ы+1
=777Е £ ^2<п.х)иГипХ(п„) —^-ТтЕ £ Ъ£.
7 п = N+1 7 п = N+1
Учитывая оценки (6), (8) и условия для выбора весов (10)-(12), получаем
1 т — N+а-1 — т
—Е £ Ъ1 У2(п,х^Ц2 < —Е £ У2(п,х)|и||2 + —Е £ у2(п,х)^ <
7 п=N+1 7 п=N+1 7 п=N+а
< с (7+р )
н2
'-E — , H+Z.
Tn H2
Теорема доказана.
3. Построение процедуры обнаружения разладки
Построим процедуру определения момента разладки. Пусть матрица
Т2
А(71,Г2) = £ у(п, х)ипитп , а ит1п П2) - ее минимальное собственное значение.
п =П1
Весовые функции у(п, х) определяются аналогично (10)-(12). Построим последовательность моментов остановки
Т0 = N; Т = т1п { > Т-1 : итіп (Т-1 + 1 п) ^ Н}, І ^ 1.
На каждом интервале [тг-1 +1, тІ ] найдем оценку параметров ДІ (Н) процесса (1), построенную по формуле, аналогичной (13). Далее сравним оценки параметров, полученные на интервалах, отстоящих друг от друга на т шагов. Если интервал [т;-1 +1, тІ ] не содержит момент разладки 0, то вектор параметров Л на этом интервале является постоянным и его значение равно или начальному значению Л0, или конечному Л1. Если, для определенного І, разница между значениями параметров на интервалах [хг-т-1 +1, хг-т ] и [хг--1 +1, тІ ] не меньше, чем заданная величина Д, то
тІ-т < 9 < хг-1. Свяжем статистику Ji с интервалом [хг--1 +1, тІ ] для І>т
•г = (Л* -Л*-т)Т (Л* -Л*-т). (17)
Эта статистика характеризует квадрат нормы различия оценок с номерами г и г—т.
Обозначим отклонение оценки Д* (Н) от истинного значения вектора параметров
через о,-. Если выполняется условие 0 > тг , то до момента тг значения параметров остаются неизменными и статистика (17) имеет вид
•г = 1К -СТ,-т|Р .
Если т,-т < 9 < тг-1, то есть изменение значений параметров произошло на интервале [тг-т,тг-1 ], то статистика • примет следующий вид:
•г =||Л1 -Л0 + аг -аг-т IГ
Алгоритм определения точки изменения параметров процесса заключается в следующем: значение заданной статистики (17) сравнивается с пороговым значением 5. Решение о наличии разладки принимается при превышении значения статистики Ы ^ значения 5.
Важными характеристиками любой процедуры обнаружения разладки являются вероятности ложной тревоги и ложного спокойствия. Благодаря использованию взвешенного метода наименьших квадратов для построения оценок, в каждом цикле наблюдений можно обеспечить заданную вероятность ложной тревоги и ложного спокойствия, выбирая параметр процедуры 7 соответствующим образом.
Теорема 2. Пусть 0 <5<Д, тогда вероятность ложной тревоги Р0 и вероятность ложного спокойствия Р1 на любом интервале наблюдений [тг-1+1, тг] являются ограниченными
Ро(7,8) < ^ *7.8) < -4^.
78 72 (Д-л/б)
Доказательство. Рассмотрим вероятность ложной тревоги. В этом случае значение статистики •] превышает порог 5 до момента разладки 0. Используя свойства нормы вектора и неравенство Чебышева, а также утверждение теоремы 1, получаем
Ро (7, 8) = РЫ > 8 | т, < 9} = Р {а,- - ст^|2 >8} 2Е{СТг 11+N1 } < 17+^.
1 ’ 8 7 8
Для нахождения вероятности ложного спокойствия рассматриваются случаи, когда момент разладки уже наступил, а значение статистики (17) не превысило пороговое значение 5. Вероятность Р1 имеет вид
Р(7, 8) = Р { <8 |т,-т <9 <Тг-1} = Р { -Л о + а( -аг-т||2 < 8} =
= Р {Л1 -Л0 + аг -аг-Л < ^} ■•
Учитывая, что ||Л1-Л0||2 >Д >0 и используя свойства нормы, неравенство Чебышева и утверждение теоремы 1, получаем
Р (7, 8) < Р {Д - ||а - аг-л ^л/5} = Р {а - а-Л о/Д-^/б} ци+р) .
7 2 (Д-^б)
Теорема доказана.
Заключение
В работе построена и исследована последовательная процедура обнаружения момента изменения значений параметров процесса GARCH(/>,q), параметры которого предполагаются неизвестными до и после момента разладки. Для построения процедуры обнаружения разладки использовались оценки неизвестных параметров процесса, вычисленные на различных интервалах наблюдений. Предлагаемый метод построения оценок основан на использовании модифицированного взве-шанного метода наименьших квадратов с гарантированным среднеквадратическим отклонением. Таким образом, точность получаемых оценок неизвестных параметров зависит от выбора заданного параметра процедуры. Найдены характеристики процедуры обнаружения разладки.
ЛИТЕРАТУРА
1. Bollerslev T. Generalized autoregressive conditional heteroskedasticity // J. Econometrics. 1986. V. 86. P. 307-327.
2. Hillebrand E. Negleting parameter changes in GARCH models// J. Econometrics. 2005. V. 129. P. 121-138.
3. Davies R.A., Huang D., Yao Y.-C. Testing for change in the parameter value and order of autoregressive model // Ann. Statist. 1995. V. 23. P. 282-304.
4. Gombey E., Serban D. Monitoring parameter change in AR(p) time series models // Statistics Centre Technical Reports 05.04, The University of Alberta, Edmonton, Canada, 2005.
5. Gombey E. Change detection in autoregressive time series // J. Multivariate Analysis. 2008. V. 99. P. 451-464.
6. Berkes I., Horvath L. The efficiency of the estimators of the parameters in GARCH processes // Annals of Statistics. 2004. V. 32. P. 633-655.
7. Berkes I., Horvath L., Kokoszka P.S. GARCH processes: Structure and estimation // Bernoulli. 2003. V. 9. P. 201-227.
8. Francq C., Zacoian J.M. Maximum likelihood estimation of pure GARCH and ARMA -GARCH processes // Bernoulli. 2004. V. 10. P. 605-637.
9. Pan J., WangH., TongH. Estimation and power-transformed and threshold GARCH models // J. Econometrics. 2008. V. 142. P. 352-378.
10. Bai J. Least squares estimation of a shift in linear process // J. Time Series Analysis. 1994. V. 15. N.5. P. 453-472.
11. Hamadeh T., Zakoian J.-M. Asymptotic properties of LS and QML estimators for a class of nonlinear GARCH processe // J. Statist. Plann.Inference. 2011. V. 141. P. 488-507.
12. Boudt K., Croux C. Robust M-estimation of multivariate GARCH models // Computational Statistics and Data Analysis. 2010. V. 54. P. 2459-2469.
13. Muler N., Yohai V.J. Robust estimates for GARCH models // J. Statistical Planning and Inference. 2008. V. 138. N. 10. P. 2918-2940.
14. Peng I., Yao Q. Least absolute deviations estimation for ARCH and GARCH models // Biometrica. 2003. V. 90. N. 4. P. 967-997.
15. Baillie R.T., Chung H. Estimation of GARCH models from the autocorrelation of the squares of a process // J. Time Ser. Anal. 2001. V. 22. No. 6. P. 631-650.
16. Storti G. Minimum distance of GARCH(1,1) models // Computattional Statistics and Data Analysis. 2006. V. 51. P. 1803-1821.
17. Буркатовская Ю. Б., Воробейников С.Э. Гарантированное обнаружение момента разладки GARCH-процесса // Автоматика и телемеханика. 2006. № 12. С. 56-70.
18. Meder N., Vorobejchikov S. On guaranteed estimation of parameters of random processes by the weighted least square method // Proc. 15th Triennial World Congress of the International Federation of Automatic Control. Barcelona, 2002. No. 1200.
19. Дмитренко А.А., Конев В.В. О последовательной классификации процессов авторегрессии с неизвестной дисперсией помех // Проблемы передачи информации. 1995. Т. 31. Вып. 4. С. 51-62.
Буркатовская Юлия Борисовна Воробейчиков Сергей Эрикович Сергеева Екатерина Евгеньевна Томский государственный университет,
E-mail: tracey@tpu.ru; sev@mail.tsu.ru; sergeeva_e_e@mail.ru Поступила в редакцию 5 ноября 2011 г.
Burkatovskaya Yulia B., Vorobeychikov Sergey E., Sergeeva Ekaterina E. (Tomsk State University, Tomsk Polytechnic University). Parameter estimation and their change-point detection for generalized autoregressive process with conditional heteroscedasticity.
Keywords: GARCH(p,g), change-point, least squares method, mean square error, guaranteed estimation.
The problem of change-point detection of the parameters of GARCH(p,q) process is considered. The utoregressive parameters of the process before and after the change point are supposed to be unknown. A sequential procedure for estimating the parameters based on the weighted least squares method is developed. The choice of the weights and the stopping rule allows one to construct an estimator with a preassigned mean square error depending on parameter H of the procedure. The procedure of change-point detection is based on comparison of parameter estimators on different observation intervals. The upper bounds for probabilities of the false alarm and the delay are found.
