Медведев Г.А., Морозов В.А. Практикум на ЭВМ по анализу временных рядов

Подождите немного. Документ загружается.

(5.59) удобно определять критич еск ие значения F

1−α

(P, N), учитывая, что

lim

N → ∞

1−α

(P, N) = F

1−α

(P, ∞) = χ

1−α

(P )/P .

Задание 2. Получить реали зац ию временного ряда с элементами вида

(5.60), (5.61), содержащими негармоническую периодическ ую составляющую с

известным целочисленным периодом. Рассмотреть два варианта, когда объем

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

это не выполняется. Исследовать несколько различных функций (5.60).

Задание 3. Рассмотреть применение тестов Фишера и Колмогорова — Смир-

нова к реализациям, которые были получены при выполнении предыдущих за-

даний. Используя критерий Колмогорова — Смирнова, нужно учитывать, что

при N → ∞ пределы имеют значения K

0,01

(N) → 1,628, K

0,02

(N) → 1,517,

0,05

(N) → 1,358, K

0,1

(N) → 1,224, K

0,2

(N) → 1,073.

Задание 4. Получить реализацию временного ряда с известной спектраль-

ной плотностью f(λ), например, при помощи имитации процесса АРСС (p, q).

Вычислить периодограмму временного ряда I(ω

). Использова ть тесты Фише-

ра и Колмогорова — Смирнова для решения воп ро са о соответствии временного

ряда спектральной плотности f (λ). При вычислении периодограммы применить

алгоритм быстрого преобразования Фурье (БПФ).

Глава 6

МНОГОМЕРНЫЕ ВРЕМЕННЫЕ РЯДЫ

6.1. Оценивание числовых характеристик

многомерных числовых рядо в

Встречающиеся на практике временные ряды, как правило, могут (а в неко-

торых случаях и дол жны) изучаться как компоненты опред еле нного много-

мерного СП. Так, например, количество осадков и урожайно сть сельскохозяй-

ственных культур, рассматриваемые как функции времени, составляют двух-

мерную векторную функцию времени, компоненты которой изменяются слу-

чайным образом по соответствующей стохастической модели и стохастически

зависимы между собой. В общем случае имеет смысл рассматривать временной

ряд {x

, t ∈ T } с элементами, являющимися n-мерными векторами x

= (x

, . . . , x

)

∗

. Для многомерных временных рядов второго порядка требуется,

чтобы M{x

} < ∞ для всех t и i. Математическим ожиданием такого времен-

ного ряда является вектор

= M{x

} = (M{x

}, M{x

}, . . . , M{x

})

∗

= (m

, m

, . . . , m

)

∗

. (6.1)

Корреляционные свойства многомерного временного ряда определяются мат-

ричной функцией, каждый элемент которой является взаимной корреляци-

онной функцией соответствующих компонент временного ряда: R(t, t + τ) =

140

= M{x

∗

t+τ

} = (M{x

t+τ,j

}) = (R

(t, t + τ)), 1 ≤ i ≤ n, 1 ≤ j ≤ n. Чаще

пользуются ковариационными матрицами (КМ) временных рядов

(t, t + τ) = M{(x

− m

)(x

t+τ

− m

t+τ

)

∗

} = (R

(t, t + τ)). (6.2)

Если элементы временного ряда распределены по нормальному закону, то со-

вокупность (m

, R

(t, t + τ )) полностью характеризует статистические свойства

временного ряда.

Для стационарных временных рядов будем иметь вместо (

6.1) и (6.2) равен-

ства

= m = (m

, m

, . . . , m

), t = 1, 2, . . . ;

(t, t + τ) = r(τ ) = (r

(τ)), t = 1, 2, . . . , (6.3)

т.е. в этом случае вектор m

не зависит от временного параметра t и является

постоянным, а матрица R

(t, t + τ) также не зависит от параметра t и опреде-

ляется только временным сдвигом τ.

Обычно удобнее пользоваться КМ, составленной из нормированных вз аим-

ных корреляционных функций компонент временного ряда

ρ(τ) = (ρ

(τ)), 1 ≤ i ≤ n, 1 ≤ j ≤ n, (6.4)

где элементы матрицы ρ(τ) задаются соотношением

(τ) = r

(τ)

(0)r

(0). (6.5)

Свойства функций r

(τ) и ρ

(τ) рассматривались в гл. 4.

Многомерный белый шум описывается как векторный процесс W

, для ко-

торого

M{W

} = 0, M{W

∗

} = Q,

M{W

∗

t+τ

} = 0, τ > 0, t = 1, 2, . . . . (6.6)

Как прав ил о, временные ряды являются реализациями СП, представляемых

в виде разложения

= m

+∞

j=−∞

t−j

, t = 1, 2, . . . , (6.7)

где C

— по след ова тел ьно сть матриц с абсолютно суммируемыми элементами.

Это значит, что если C

— матрица с элементами C

(k, l), то для любой пары

k, l выполняется неравенство

+∞

j=−∞

(k, l)| < ∞. При C

= 0 для всех j < 0

данное разложение известно под назван ием “разложение Волда”. Для процесса

(

6.7) имеем

M{x

} = m

, r(τ ) =

+∞

j=−∞

∗

j+τ

. (6.8)

141

Рассмотрим стационарный временной ряд {x

, 1 ≤ t ≤ N } с неизвес тными m

и r(τ ), где x

= (x

, x

, . . . , x

)

∗

. Несмещенной оценкой m является выборочное

среднее

m =

t=1

= ( ˆm

, ˆm

, . . . , ˆm

)

∗

, (6.9)

где ˆm

t=1

, 1 ≤ i ≤ n.

Для стационарных временных рядов сходимость оце нок (

6.9) к истинным

значениям характеризуется следующими результатами.

Если r

(τ) → 0 при τ → ∞, 1 ≤ i ≤ n, то

lim

N → ∞

M{(

m − m)

∗

(

m − m)} = 0. (6.10)

Если

+∞

τ=−∞

(τ)| < ∞, 1 ≤ i, j ≤ n, то

lim

N → ∞

NM{(

m − m)

∗

(

m − m)} =

i=1

+∞

τ=−∞

(τ) = r

∗

. (6.11)

Оценка

m в (

6.10), (6.11) находится из формулы (6.9) и является асимптотически

нормальной, если выполняются некоторые дополнительные предположения.

Если x

представляется в виде (6.7), где W

— НОР СВ со свойствами (6.6), а

— последовательность (n ×n)-матриц с абсолютно суммируемыми элемента-

ми, то

m, вычисляемое по формуле (6.9), является асимптотически нормальным

с параметрами

m} = m, M{(

m − m)

∗

(

m − m)} =

+∞

k=−∞

+∞

l=−∞

∗

= G

. (6.12)

Это свойство может быть использовано при построении доверительных обла-

стей для m. Например, если КМ G

в (

6.12) является невырожденной и извест-

ной, то асимптотическая доверительная область для m на уровне значимости α

определяется как множество

{m |(m −

∗

−1

(m −

m) ≤ χ

1−α

(n)}, (6.13)

где χ

1−α

(n) — соответствующий квантиль χ

-распределения с n степенями сво-

боды. Такое построение доверительной области практически не используется,

так как редко матрица G

оказывается известной в случае, когда m неизвестно.

Если может быть найдена состоятельная оценка

матрицы G

, то в (6.13)

подставляют ее вместо G

и получают необходимую доверительную область.

Однако в общем случае трудно находить состоятельную оценку матрицы G

142

Более простой п одход — построение доверительных интервалов для каждой ком-

поненты m

вектора m, основанное на выборках {x

, x

, . . . , x

}, 1 ≤ i ≤ n.

Комбинируя найденные таким образом доверительные интервалы, можно пол у-

чить доверительную область для вектора m.

Если условия справедливости (6.12) выполняются, то компонента ˆm

выбо-

рочного среднего (6.9) является а симпт оти ческ и нормальной с математическим

ожиданием m

и дисперсией

+∞

τ=−∞

(τ), 1 ≤ i ≤ n. (6.14)

Пусть {τ

} — последовательность чисел, удовлетворяющих условию

(τ

/N) → 0 при N → ∞. Состоятельной оценкой является

|τ|≤τ

1 −

|τ|

ˆr

(τ), 1 ≤ i ≤ n, (6.15)

где ˆr

(τ) — оценка автоковариационной функции i-й компоненты x

век-

тора x

Состоятельной оценкой автоковариационной функции ˆr

(τ) является следу-

ющая выборочная ковариационная функция:

ˆr

(τ) =

N − τ

N−τ

t=1

− ˆm

)(x

t+τ,i

− ˆm

), 1 ≤ i ≤ n. (6.16)

Таким образом, на уровне значимости α доверительный интервал дл я i-й ком-

поненты вектора m при достаточно больших N задается соотно шени ем

− ˆm

| ≤ Φ

−1

1 −

, 1 ≤ i ≤ n, (6.17)

где Φ

−1

(β) — квантиль нормального распределения на уровне β;

вычисляется

по формулам (6.15), (6.16).

Так как

P {|m

− ˆm

| ≤ Φ

−1

(β)

, 1 ≤ i ≤ n} ≥ 1 −

i=1

P {|m

− ˆm

| ≤ Φ

−1

(β)

то доверительная область для вектора m, в которой он находится с вероятно-

стью, не меньшей чем (1 − α), определяется следующим образом:

− ˆm

| ≤ Φ

−1

1 −

, 1 ≤ i ≤ n

. (6.18)

Для больших значений n эта доверительная область будет значительно шире

области, в которую попадает m с вероятностью (1 −α). Тем не менее, поскольку

143

(6.18) строить немного легче, чем (6.13), для небольших n удобнее применять

доверительную область (6.18).

В качестве оценки КМ r(τ ) естественно пользоваться выборочной КМ

(ВКМ), которая задается соотношением

ˆr(τ ) = (ˆr

(τ)) =

N−τ

t=1

t+τ

−

m)(x

−

∗

, 0 ≤ τ ≤ N − 1. (6.19)

Через элементы матрицы определим оцен ки элементов КМ, составленной из

нормированных взаимных ковариационных функций

ˆρ

(τ) = ˆr

(τ)/

ˆr

(0)ˆrjj(0), 1 ≤ i, j ≤ n. (6.20)

В условиях справедливости формул (6.12) оценки (6.19) и (6.20) для двух-

мерных временных рядов являются сос тоятельными по вероятности, т.е. п ри

N → ∞ для всякого фиксированного τ ≥ 0

ˆr

(τ) → r

(τ), ˆρ

(τ) → ρ

(τ), 1 ≤ i, j ≤ 2. (6.21)

Скорость сходимости в этом случае характеризуется сле дующими соотно-

шениями. Если матрица Q в (6.6) является диагональной (т.е. компоненты W

независимы между собой), то для всякого τ ≥ 0 распределение ˆρ

(τ) асимпто-

тически нормально с нулевым средним и дисперсией

+∞

j=−∞

(j)ρ

(j). (6.22)

Если τ

6= τ

, то вектор (ˆρ

(τ

), ˆρ

(τ

))

∗

асимптотически нормален с нулевым

средним компонент, дисперсиями компонент, определяемыми (

6.22), и ковариа-

цией компонент, вычисляемой по формуле

+∞

j=−∞

(j)ρ

(j + τ

− τ). (6.23)

Это дает возможность конструировать различные тесты относительно ком-

понент при анализе их корреляционных свойств. Например, если одна из ком-

понент является белым шумом, то ˆρ

(τ) имеет дисперсию d

= N

−1

, какими

бы свойствами ни обладала другая компонента.

Рассмотрим задачу проверки на независимость двух стационарных компо-

нент временного ряда. Как следует из (

6.22), асимптотическое распределение

ˆρ

(τ) зависит как от ρ

(·), так и от ρ

(·). Поэтому тест независимости двух

компонент ряда не может быть основан только на оценива нии значений ρ

(τ)

без учета других свойств компонент ряда. Эту трудность можно преодолеть

144

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

борочную взаимную ковариацию ρ

(τ). В данном случае декорреляция дости-

гается соответствующим линейным преобразованием. В частности, если {x

} и

} — обратимые процессы АРСС (p, q), таким преобразовани ем является

∞

j=0

(i)

t−j,i

, i = 1, 2, (6.24)

где коэффициенты θ

(i)

определяются из разложения

∞

j=0

(i)

= a

(i)

(u)/b

(i)

(u), |u| < 1, i = 1, 2,

в котором a

(i)

(u) и b

(i)

(u) — полиномы, задающие АР и СС соответственно.

Так как на практике истинная модель обычно неизвестна и наблюдений {x

t ≤ 0} не имеется в налич ии , удобно за менит ь последовательности {z

}, i = 1, 2,

остатками {

, 1 ≤ t ≤ N}, которые рассматривались в § 4.3. Если подобранная

модель АРСС (p, q) является истинной, то последовательности остатков будут

последовательностями белого шума, а с учетом асимптотической нормально-

сти — последовательностями НОР СВ. Таким образом, при справедливой ги-

потезе о независимости компонент временного ряда выборочные автокорреля-

ции ˆρ

(τ

) и ˆρ(τ

), τ

6= τ

, последовательностей {z

} и {z

} асимптотически

нормальны, независимы, имеют нулевые сред ни е и дисперсии, равные N

−1

. По-

этому тест на независимость компонент на уровне значимости α заключается в

проверке неравенства

|ρ

(τ)| ≤ Φ

−1

1 −

−1/2

, τ = 0, 1, 2, . . . , (6.25)

где Φ

−1

(β) — квантиль стандартного нормального распределения на уровне β.

Если декоррелировать только одну из двух исходных компонент (например,

}), то при справедливой гипотезе о независимости компонент из (6.23) сле-

дует, что выборочные автокорреляции ˆρ

(τ

) и ˆρ

(τ

), τ

6= τ

, последователь-

ностей {z

} и {x

} асимптотически нормальны с нулевым средним, дисперсией

−1

и ковариацией N

−1

(τ

−τ

), где ρ

(·) — АКФ последовательности {x

Следовательно, и в этом случае тест на независимость компонент на уров не зна-

чимости α заключается в проверке выполнения неравенс тва (6.25) для любого

значения τ.

При получении (6.22) и (6.23), на которых основан описанный тест, пр едп ола -

галось, что компоненты временного ряда удовлетворяют условиям разложения

(6.6). Для нормальных временных рядов этого можно и не требовать, поскольку

для них справедлива формула Бартлетта. Если {x

} — двухмерный гауссовский

процесс и его АКФ удовлетворяет неравенству

+∞

τ=−∞

| r

(τ)| < ∞, 1 ≤ j, i ≤ 2,

145

то

lim

N → ∞

Ncov(ˆρ

(τ

), ˆρ

(τ

)) =

+∞

j=−∞

(j)ρ

(j + τ

− τ

+ρ

(j + τ

)ρ

(j − τ

) − ρ

(τ

)

(j)ρ

(j + τ

) + ρ

(j)ρ

(j − τ

)

−

−ρ

(j)

(j)ρ

(j + τ

) + ρ

(j)ρ

(j − τ

)

+ ρ

(τ

)ρ

(τ

)

(j) + ρ

(j) +

(j)

´´

. (6.26)

Отсюда, в частности, следует, что если компонентами {x

} являются белые

шумы и ρ

(τ) = 0, τ ∈ [a, b], то

lim

N → ∞

ND{ˆρ

(τ)} = 1, τ ∈ [a, b], (6.27)

т.е. асимптотическая дисперсия ˆρ

(τ) в этом случае тоже равна N

−1

Среди мно гомерн ых СП своей практической распространенн остью выделя-

ются про цесс ы АРСС, которые в многомерном варианте определ яю тся как про-

цессы {x

}, удовлетворяющие многомерному разностному уравнению

i=1

t−i

+ W

j=1

t−j

, t = 1, 2, . . . . (6.28)

Здесь x

= (x

, x

, . . . , x

)

∗

— n-вектор; A

, . . . , A

, B

, . . . , B

—

(n × n)-матрицы; {W

} — случайные векторы с компонентами: W

= (W

, . . . , W

)

∗

, для которых выполняются условия (

6.6). Последовательность

векторов x

, заданных уравнением (6.28), называется многомерным (при n > 1)

процесом АРСС (p, q).

Для рассмотрения свойства многомерного процесса АРСС (p, q) опред ели м

матричные полиномы

A(z) = I − A

z − . . . − A

, B(z) = I + B

z + . . . + B

, (6.29)

где I — единичная (n × n)-матрица.

Критерий каузальности процесса x

. Есл и detA(z) 6= 0 для всех |z| ≤ 1,

то (6.28) имеет единственное стационарное решение x

∞

j=0

t−j

, где (n× n)-

матрицы Ψ

удовлетворяют равенствам

∞

j=0

= Ψ(z) = A

−1

(z)B(z), |z| ≤ 1. (6.30)

Критерий обратимости процесса x

. Если detB(z) 6= 0 д ля всех |z| ≤ 1

и {x

} — стационарное решение (

6.28), то W

∞

j=0

t−j

, где матрицы Φ

146

определяются равенством

∞

j=0

= Φ(z) = B

−1

(z)A(z), |z| ≤ 1. (6.31)

Матрицы Ψ

и Φ

находятся из следующих рекуррентных соотношений:

= I, Ψ

i=1

j−i

+ B

, j = 1, 2, . . . ; (6.32)

= I, Φ

= −

i=1

j−i

− A

, j = 1, 2, . . . , (6.33)

где Ψ

= 0 для j < 0; A

= 0 для j > p; B

= 0 для j > q.

Ковариационная матр ичн ая функция каузального процесса АРСС (p, q)

определяется равенством

r(τ) = (r

(τ)) = M{x

∗

t+τ

} =

∞

k=0

QΨ

∗

k+τ

, τ ≥ 0, (6.34)

где Q — КМ белого шума W

со свойствами (

6.6). Существует γ ∈ (0, 1) и

константа K, такие, что элементы матрицы r(τ) удовлетворяют неравенствам

(τ)| < Kγ

|τ|

, 1 ≤ i, j ≤ n, τ = 0, ±1, ±2, . . . . (6.35)

КМ r(τ) может быть также определена из уравнения Юла — Уолкера

r(τ) −

i=1

r(τ − i) =

i=τ

QΨ

∗

i−τ

, τ = 0, 1, 2, . . . . (6.36)

Первые (p + 1) уравне ний (

6.36) могут быть р азр ешен ы относительно матриц

r(0), r(1), . . . , r(p), если использовать тот факт, что r(−τ ) = r

∗

(τ). КМ r(τ ) для

τ > p получаются рекурр ент но из (6.3 6 ).

Производящая функция КМ, представленная

G(z) =

+∞

τ=−∞

r(τ)z

выражается в виде

G(z) = Ψ(z)QΨ

∗

−1

) = A

−1

(z)B(z)QB

∗

−1

)(A

−1

))

∗

. (6.37)

Таким образом, если набор матриц {A

} и {B

} известен , то корреляционные

характеристики процесса (

6.28) могут быть найдены из (6.32), (6.34) или (6.36).

В том случае, когда {A

}, {B

} неизвестн ы, необходимо строить процедуру их

оценивания.

147

Рассмотрим случай гауссовского процесса АРСС (p, q). При этом для оце-

нивания параметров процесса может быть использован ММП. Чтобы по-

строить гаус совс кую функцию правдоподобия (ГФП) для выборки {x

1 ≤ t ≤ N}, составим вектор с nN компонентами x = (x

, . . . , x

)

∗

= (x

∗

, x

∗

, . . . , x

∗

)

∗

. Матрица ковариации такого

вектора имеет размер (nN ×nN ) и может б ыть представлена в блочной форме:

R = M{xx

∗

} =







r(0) r(1) . . . r(N − 1)

∗

(1) r(0) . . . r(N − 2)

. . . . . . . . . . . .

∗

(N − 1) r

∗

(N − 2) . . . r(0)







. (6.38)

В (

6.38) блоками матрицы R являются матрицы r(τ ), 0 ≤ τ ≤ N −1, определен-

ные в (6.34).

ГФП с учетом введенных обозначений имеет вид

L(R) = (2π)

−nN/2

(detR)

−1/2

exp

−

∗

−1

. (6.39)

Так как матрица R — симметри чес кая, она задается nN (nN + 1)/2 элементами.

Для оценивания этих пара метров необходимо найти такую совокупность элемен-

тов, которая бы доставляла L(R) максимум. Эта задача чрезвычайно сложная

в вычислительном отношении и вряд ли удобная в с мысле практического во-

площения. Вместе с тем элементы матрицы R, точнее ее блоки, находятся по

формулам (

6.32) — (6.34), выражаются через {A

}, {B

}, Q и содержат в общей

сложности около n

(

+ p + q) независимых элемен тов , поэтому вычислить зна-

чения ГФП (

6.39) можно, если известны n

(

+p+q) чисел, при помощи (

6.32) —

(6.34) и (6.38). Поиск максимума ГФП (6.39) в пространстве этих переменных —

задача менее трудная, чем указанная выше, так как N обычно существенно боль-

ше, чем n. Тем не менее она остается еще достаточно сложной. Более удоб ным

является использование при построении ГФП оценок предсказания, о которых

речь пойдет в § 6.2, где будет обсуждена и проблема оценивания параметров.

В случае, когда процесс {x

} не является гауссовским, получать оценки неиз-

вестных корреляционных характеристик или параметров разностного уравнения

(6.28) можно точно таким же путем, максимизируя ГФП (6.39). Но только в

этом случае полученные оценки не будут оценками ММП, а будут М-оценками

неизвестных параметров, идея построения которых рассматривалась в § 2.4.

Лабораторная работа 22. Оценивание числовых характеристик

многомерных временных рядов

Цель работы. Освоить технику статистической обработки временных рядов.

Провести сравнительный анализ различных подходов к построен ию многомер-

ных доверительных областей.

Порядок выполнения работы

148

Построить модель стационарного многомерного временного ряда с n > 2 на

основе разложения (6.6). Для этого задать конечную последовательность матриц

, α ≤ j ≤ β}, такую, что α < 0 < β; задать положительно определенную

матрицу Q, у которой отдельные, но не все, недиагональные элементы были

бы нулевыми; задать ве кто р m математических ожиданий компонент m

ряда

}. Провести имитацию этого многомерного временного ряда и получить его

реализацию {x

, 1 ≤ t ≤ N }.

Задание 1. Основываясь на (6.8), (6.5), вычислить матрицы r(τ ) = (r

(τ)),

(ρ

(τ)), а также величину r

∗

i=1

+∞

τ=−∞

(τ) и матрицу G

по формуле (

6.12).

Поскольку матриц ы C

= 0, j 6∈ [α, β], неограниченные суммы в формулах для

∗

и G

окажутся ограниченными.

Задание 2. Построить оценки

m, ˆr(τ) и ˆρ(τ ) по формулам (6.9), (6.19) и

(6.20) соответственно. Выяснить, насколько сильно отличается от нуля (

m−

−m)

∗

(

m −m) и от r

∗

величина N(

m −m)

∗

(

m −m), что подтверждает свойства

(6.10) и (6.11). Сравнить матрицы (

m − m)(

m − m)

∗

и G

, проиллюстр ир овав

свойство (6.12). Провести сравнительный анализ матриц ˆr(τ) и r(τ), а также

ˆρ(τ) и ρ(τ ) для некоторых значений τ = 0, 1, 2, . . . .

Задание 3. Для нескольких значений α постр оит ь доверительные области

(6.13) и (6.18). Найти объемы этих областей и сравнить их. Объем доверительной

области (6.18) вычисляется как объем n-мерного параллелепипеда

(18)

2Φ

−1

1 −

´´

i=1

Доверительная област ь (

6.13) — это n-мерный эллипсоид. Его объем вычис-

ляется по формуле

n/2

|detT |

1 +

(χ

1−α

(n))

n/2

i=1

где T — матрица, с тол бца ми которой являются собственные векторы, соответ-

ствующие указанным собственным значениям; Γ(·) — гамма-функция; {λ

} —

собственные значения матриц ы G

−1

. Иначе гово ря, для вычисления объема n-

мерного эллипсоида, задаваемого матрицей G

−1

, необходимо диагонализировать

эту матрицу, т.е. представить ее в виде G

−1

= T diag{λ

∗

Задание 4. Выбрать две некоррелированные компоненты исследуемого вре-

менного ряда и проверит ь их на независимость, используя декорреляцию (6.24)

и критерий (6.25). Проделать то же самое для двух коррелированных компонент

рассматриваемого временного ряда.

Задание 5. Задать две разли чные модели каузальных обратимых процессов

АРСС (p, q). Провести имитацию двухмер ног о времен ног о ряда с компонентами,

149