Analiza Szeregów Czasowych/Procesy stochastyczne

Z Skrypty dla studentów Ekonofizyki UPGOW

(Różnice między wersjami)
(Wariancja i odchylenie standardowe)
Linia 97: Linia 97:
: <math> Var(aX) = a^2 Var(X) = a^2 \sigma^2_X </math>
: <math> Var(aX) = a^2 Var(X) = a^2 \sigma^2_X </math>
oraz dla dwóch zmiennych niezależnych ''X'' i ''Y''
oraz dla dwóch zmiennych niezależnych ''X'' i ''Y''
-
: <math> Var(X + Y) = Var(X) + Var(Y)= \sigma^2_X \sigma^2_Y.</math>
+
: <math> Var(X + Y) = Var(X) + Var(Y)= \sigma^2_X + \sigma^2_Y.</math>
'''Odchylenie standardowe''' <math> \sigma_X </math> zmiennej losowej ''X'' jest to dodatni pierwiastek kwadratowy z wariancji ''X''.
'''Odchylenie standardowe''' <math> \sigma_X </math> zmiennej losowej ''X'' jest to dodatni pierwiastek kwadratowy z wariancji ''X''.

Wersja z 02:42, 30 gru 2010

Analiza Szeregów Czasowych
<<< Wstęp |  Dekompozycja szeregu czasowego>>>

Spis treści

Elementy teorii prawdopodobieństwa

Teoria prawdopodobieństwa jest zasadniczym narzędziem służącym do analizy szeregów czasowych. Jest to narzędzie mówiące nam tyle, że nie potrafimy jednoznacznie przewidywać jak zachowa się nasz układ i czy przewidywane wartości z właśnie znalezionego modelu będą rzeczywiście zrealizowane.

Aby zrozumieć kolejne podrozdziały wymagana jest elementarna wiedza o zbiorach. Można ją znaleźć we wstępie do siostrzanego podręcznika traktującego o procesach losowych. Wiele innych informacji niezbędnych do zrozumienia tego rozdziału omówiono właśnie tam. Część informacji zawartych w tym skrypcie będzie się powtarzać. Moim celem jednak było zachowanie najważniejszych rzeczy w całości, a odwoływanie się do wiedzy raczej ogólnej. Ponadto tutaj skupimy się raczej na zmiennych dyskretnych.

Prawdopodobieństwo 
definiujemy w odniesieniu do przestrzeni zdarzeń S, będącej zbiorem, którego elementy nazywamy zdarzeniami elementarnymi. Każde zdarzenie elementarne może być traktowane jako wynik eksperymentu, bądź kolejne notowanie wielkości instrumentu rynku. W eksperymencie polegającym na rzucie dwiema monetami na przestrzeń zdarzeń możemy patrzeć jako na złożoną ze zbioru wszystkich możliwych 2-słów[1]. nad zbiorem {O,R}.
S = {OO, OR, RO, RR}.
Zdarzenie 
jest podzbiorem przestrzeni zdarzeń S. Np. przy rzucie dwiema monetami zdarzenie pozyskania jednego orła i jednej reszki to {OR, RO}. Zdarzenie S nazywamy pewnym, zdarzenie \(\varnothing\) nazywane jest zdarzeniem niemożliwym. Mówimy, że dwa zdarzenia są wzajemnie wykluczające się (rozłączne) gdy \( A \cap B = \varnothing \). Z definicji wszystkie zdarzenia elementarne wzajemnie się wykluczają.

Aksjomaty teorii prawdopodobieństwa

Rozkład prawdopodobieństwa \(P(\cdot)\) 
na przestrzeni zdarzeń S jest to przyporządkowanie zdarzeniom z S liczb rzeczywistych w taki sposób, że spełnione są następujące aksjomaty
  1. \(P(A) \ge 0 \) dla każdego zdarzenia A. P(A) nazywamy prawdopodobieństwem zdarzenia A.
  2. \(P(S) = 1 \ \) (wybór jedynki to w zasadzie warunek normalizacji).
  3. \(P(A \cup B) = P(A) + P(B) \) dla każdych wzajemnie wykluczających się zdarzeń A i B.

Z powyższych aksjomatów i elementarnej teorii zbiorów natychmiast wynikają wnioski

  • \(P(\varnothing) = 0, \)
  • jeżeli \(A \subseteq B \ \) to \( P(A) \le P(B), \ \)
  • jeżeli \( \bar{A} = S - A \ \) to \( P(\bar{A}) = 1 - P(A). \ \)
  • dla dowolnych dwóch zdarzeń A i B zachodzi \(P(A \cup B) = P(A) + P(B) - P(A \cap B) \le P(A) + P(B). \ \)
Przykład 
Rzut dwoma monetami. Każde z czterech zdarzeń elementarnych ma prawdopodobieństwo 1/4. Wówczas prawdopodobieństwo otrzymania conajmniej jednego orła wynosi
\( P(\{OO, OR, RR \}) = P(OO) + P(OR) + P(RR) = \frac{3}{4}. \)

Dyskretny rozkład prawdopodobieństwa

Rozkład prawdopodobieństwa jest dyskretny, jeżeli jest zdefiniowany nad skończoną lub przeliczalną przestrzenią zdarzeń. Jeżeli S to przestrzeń zdarzeń to dla każdego zdarzenia A zachodzi

\( P(A) = \sum_{s \in A} P(s) \ \)

ponieważ zdarzenia elementarne wzajemnie się wykluczają.

Jeżeli S jest skończona i prawdopodobieństwo każdego zdarzenia elementarnego jest równe i wynosi

\( P(s) = |S|^{-1} \ \)

to mamy do czynienia z jednostajnym rozkładem prawdopodobieństwa na S. W takim przypadku eksperyment może być opisywany jako losowy.

Prawdopodobieństwo warunkowe

Czasami zdarza się że dysponujemy wiedzą częściową na temat danego procesu. Jeżeli rozpatrujemy szeregi czasowe to taka sytuacja jest na porządku dziennym - znamy przebieg danej zmiennej w konkretnej realizacji w przeszłości i na podstawie tych danych możemy z pewnym prawdopodobieństwem (warunkowym właśnie) prognozować przyszłość.

Prosty przykład z monetą i dwukrotnym rzutem pozwoli wyjaśnić ideę warunkowego prawdopodobieństwa. Jeżeli ktoś rzuci dwukrotnie monetą i powie nam, że w jednym (nie wiemy którym) rzucie wypadł orzeł a następnie zapyta "jakie jest prawdopodobieństwo, że dwukrotnie wypadł orzeł?" to powinniśmy właśnie skorzystać z własności prawdopodobieństwa warunkowego, by na to pytanie odpowiedzieć poprawnie. Jeżeli mamy już raz O, to możliwości dla drugiego rzutu (o którym nic nie wiemy) są nadal {O,R}. Zatem na pewno nie będziemy mieli dwóch reszek - czyli możemy śmiało wykluczyć RR. Nasza przestrzeń zdarzeń kurczy się więc do 3 możliwości: OO, OR, RO. Zatem prawdopodobieństwo, że wyrzuciliśmy OO jest 1/3!

Prawdopodobieństwo warunkowe określa w sposób formalny pojęcie posiadania wiedzy na temat zjawiska. Prawdopodobieństwo warunkowe zdarzenia A danego tak, że jednocześnie zachodzi inne zdarzenie B jest zdefiniowane następująco

\( P(A|B) = \frac{P(A \cap B)}{P(B)}, \) oczywiście pod warunkiem \( P(B) \ne 0. \ \)

Oznaczenie P(A|B) czytamy "prawdopodobieństwo że A pod warunkiem B". Jest to stosunek prawdopodobieństwa \( P(A \cap B) \ \) do zdarzenia \(P(B) \ \). Wracając do powyższego przykładu z dwoma rzutami (lub dwoma monetami) mamy: prawdopodobieństwo, że na obu monetach wypadł orzeł \( P(A \cap B) = 1/4 \) a prawdopodobieństwo, że mamy choć jednego orła w obu rzutach jest \(P(B) = 3/4 \). Dzieląc pierwsze przez drugie dostajemy \( P(A|B) = 1/3. \)

Jeżeli dwa zdarzenia A i B są niezależne (np: dwa niezależne od siebie rzuty monetą) to

\( P(A|B) = P(A). \ \)

Zdarzenia należące do rodziny zdarzeń \(A_1, A_2, \dots, A_n \) są parami niezależne jeżeli

\( P(A_i \cap A_j) = P(A_i) P(A_j) \ \)

dla wszystkich \( 1 \le i < j \le n. \) Mówimy, że zdarzenia są wzajemnie niezależne (lub po prostu niezależne) jeżeli każdy k-ty podzbiór \(A_{i_1}, A_{i_2}, \dots, A_{i_k} \) rodziny, gdzie \( 2 \le k \le n \) i \( 1 \le i_1 < i_2 < \dots < i_k \le n \) spełnia

\( P(A_{i_1} \cap A_{i_2} \cap \dots \cap A_{i_k}) = P(A_{i_1}) P(A_{i_2})\dots P(A_{i_k}). \ \)

Twierdzenie Bayesa

Z definicji prawdopodobieństwa warunkowego wynika, iż dla dwóch zdarzeń A i B zachodzi

\( \frac{P(A \cap B)} = P(B) P(A|B) = P(A) P(B|A). \)

Przekształcając to równanie możemy dostać wzór na prawdopodobieństwo że A pod warunkiem B

\( P(A|B) = \frac{P(A) P(B|A)}{P(B)} \)

który znany jest jako wzór Bayesa. Wiemy (z teorii zbiorów), że \( B = (B \cap A) \cup (B \cap \bar{A})\), co daje nam wzór na \( P(B) = P(B \cap A) + P(B \cap \bar{A}) = P(A) P(B|A) + P(\bar{A}) P(B|\bar{A}) \) i jednocześnie alternatywny zapis wzoru Bayesa

\( P(A|B) = \frac{P(A) P(B|A)}{P(A) P(B|A) + P(\bar{A}) P(B|\bar{A})}. \)

Zdarza się, że powyższe wzory potrafią bardzo skutecznie uprościć obliczanie prawdopodobieństw warunkowych.

Dyskretne zmienne losowe

Zmienna losowa (dyskretna) jest funkcją ze skończonej lub przeliczalnej przestrzeni zdarzeń S w zbiór liczb rzeczywistych. Przyporządkowuje ona liczbę rzeczywistą każdemu możliwemu wynikowi eksperymentu, co z kolei daje nam możliwość operowania wynikowym rozkładem prawdopodobieństwa określonym właśnie na zbiorze liczb.

Dla zmiennej losowej X i liczby rzeczywistej x możemy zdefiniować zdarzenie X = x jako \( \{ s \in S: X(s) = x \} \), czyli

\( P(X = x) = \sum_{ s \in S: X(s) = x} P(s). \ \)

Funkcja

\( f(x) = P(X = x) \ \)

nazywana jest funkcją gęstości prawdopodobieństwa zmiennej losowej X.

Często na tej samej przestrzeni zdarzeń definiuje się wiele zmiennych losowych. Wiemy jednak, że będziemy zajmować się szeregami czasowymi i pojedynczymi zmiennymi losowymi, a w zasadzie przebiegiem czasowym konkretnych realizacjami procesu rządzącego daną zmienną.

Wartość oczekiwana

Najprostszą i najbardziej użyteczną charakterystyką rozkładu zmiennej losowej jest jej średnia, czyli wartość oczekiwana określona wzorem

\( EX = \sum_x x P(X=x). \ \)

Czasami stosuje się zapisy:

  • \( \langle X \rangle, \)
  • \( \mu_X \) (lub po prostu \( \mu \) gdy zmienna określona jest jednoznacznie).
Kilka własności 
  1. Wartość oczekiwana sumy dwóch zmiennych jest równa sumie dwóch wartości oczekiwanych
    E(X+Y) = EX + EY.
  2. Dla dowolnej funkcji g(x)
    \(E(g(x)) = \Sigma g(x) P(X=x). \)
  3. Dla dowolnej stałej a
    E(aX) = a EX.
  4. Dla dwóch dowolnych zmiennych losowych X, Y i stałej a zachodzi
    E(aX + Y) = sEX + EY.
  5. Jeżeli dwie zmienne losowe są niezależne, to
    E(XY) = EX EY.

Wariancja i odchylenie standardowe

Wariancja zmiennej losowej X o wartości oczekiwanej EX wynosi

\( \sigma^2_X = Var(X) = E[(X - EX)^2] = E(X^2) - E^2X \)

Podobnie jak wyżej

\( Var(aX) = a^2 Var(X) = a^2 \sigma^2_X \)

oraz dla dwóch zmiennych niezależnych X i Y

\( Var(X + Y) = Var(X) + Var(Y)= \sigma^2_X + \sigma^2_Y.\)

Odchylenie standardowe \( \sigma_X \) zmiennej losowej X jest to dodatni pierwiastek kwadratowy z wariancji X.

Procesy stochastyczne

Proces stochastyczny - rodzina zmiennych losowych określonych na pewnej przestrzeni probabilistycznej S o wartościach w pewnej przestrzeni mierzalnej, najczęściej w \(\R\).

Definicja 

Niech \(T\) będzie niepustym zbiorem (zbiór indeksów), S będzie przestrzenią probabilistyczną oraz M będzie przestrzenią mierzalną. Rodzinę zmiennych losowych

\(X=(X_t)_{t\in T}\),

nazywamy procesem stochastycznym. Przestrzeń M nazywamy przestrzenią fazową albo przestrzenią stanów procesu \(X\).

Dużo bardziej dokładna definicja, własności i inne bajery podane są w sposób bardzo przystępny i szczegółowy w sąsiednim podręczniku i autor nie do końca widzi potrzebę powtarzania tych samych informacji w tym miejscu. Na nasze potrzeby wystarczy pamiętać, że w zagadnieniach poruszanych w tym skrypcie dziedziną procesu losowego zawsze będzie czas, a w szeroko pojętej ekonometrii wartości swe ów proces przybierał będzie ze zbioru liczb rzeczywistych.

W dalszej części zajmiemy się kilkoma charakterystykami procesów losowych (stochastycznych) najbardziej użytecznymi do analizy owej specjalnej klasy tych procesów zwanych szeregami czasowymi.

Definicja i rola funkcji autokowariancji (autokorelacji)

Funkcja kowariancji

Do wyznaczania zależności pomiędzy zmiennymi losowymi użyteczna bywa funkcja kowariancji.

Definicja 3.1
Dla dwóch zmiennych losowych \( \{X_t, t \in T\}\ \) oraz \( \{Y_s, s \in T\}\ \) funkcja
\( \begin{align} ~cov(X(r),Y(s)) = &E[(X_r - EX_r)(Y_s - EY_s)] = \\ &E(X_tY_s) - EX_t EY_s ~~~\text{dla} ~~~ r,s \in T \end{align} \)

określa liniową zależność pomiędzy powyższymi zmiennymi losowymi. Stopień współzależności owych zmiennych losowych można podać za pomocą tzw. współczynnika korelacji Pearsona \( r_{XY}\ \)

\( cov (X, Y) = r_{XY} \sigma_{X} \sigma_{Y}. \)

Wartość współczynnika korelacji Pearsona mieści się w przedziale domkniętym [-1, 1]. Im większa jego wartość bezwzględna, tym silniejsza jest zależność zmiennych losowych między zmiennymi. \(r_{XY} = 0\) oznacza brak liniowej zależności między cechami, \(r_{XY} = 1\) oznacza dokładną dodatnią liniową zależność między cechami, natomiast \(r_{XY} = -1\) oznacza dokładną ujemną liniową zależność między cechami, tzn. jeżeli zmienna \(X\) rośnie, to \(Y\) maleje i na odwrót. Współczynnik korelacji liniowej można traktować jako znormalizowaną kowariancję. Korelacja przyjmuje zawsze wartości w zakresie [-1, 1], co pozwala uniezależnić analizę od dziedziny badanych zmiennych.

Funkcja autokowariancji

W przypadku gdy analizujemy szereg czasowy opisywany poprzez ewolucję jednej zmiennej losowej możemy mówić najwyżej o funkcji autokowariancji. Dla szeregu czasowego \( \{X_t, t \in T\}\ \) możemy taką funkcję zdefiniować następująco.

Definicja 3.2
Jeżeli \( \{X_t, t \in T\}\ \) jest procesem dla którego wariancja zmiennej losowej dla każdej chwili czasu \( \sigma_{X_t} \) jest skończona, wtedy funkcja autokowariancji procesu \( \{X_t\}\ \) zdefiniowana jest jako
\( \begin{align} ~\gamma_X(r,s) = &K_{XX}(r,s) = cov(X(r),X(s)) = cov(X_r,X_s) = \\ &E[(X_r - EX_r)(X_s - EX_s)] = E(X_tX_s) - EX_t EX_s ~~~\text{dla} ~~~ r,s \in T. \end{align} \)

Analogicznie do funkcji kowariancji, autokowariancja określa liniową zależność pomiędzy tą samą zmienną losową w dwóch chwilach czasu t i s.

Funkcja autokorelacji

Jeżeli dowolny proces losowy \( \{X_r, r \in T\}\ \) posiada wartość oczekiwaną \( EX_r \) oraz wariancję \( \sigma_{X} \) to możemy zdefiniować funkcję autokorelacji procesu (swego rodzaju unormowaną funkcję autokowariancji) jako

Definicja 3.4
\( \rho(r,s) = \frac{E[(X_t - EX_t )(X_s - EX_s )]}{\sigma_{X}^2}\, , \)
Jeżeli \(k=t-s\;\) to możemy funkcję autokorelacji zapisać jako funkcję jednej zmiennej:
\( \rho(k) = \frac{E[(X_t - EX_t )(X_{t-k} - EX_{t-k} )]}{\sigma_{X}^2}\, . \)

Stacjonarność procesu stochastycznego

Stacjonarność w sensie szerokim

Definicja 3.3
Szereg czasowy \( \{X_t, t \in \Z\}\ \), gdzie zbiór indeksów zdefiniowany jest jako \( \Z = \{0, \pm 1, \pm 2,\cdots \}\) nazywamy stacjonarnym (w sensie słabym) jeżeli spełnione są poniższe punkty
\( \begin{align} (i) &~E | X_t |^2 < \infty ~~~ \text{dla} ~~~ t \in \Z \\ (ii) &~E X_t = m ~~~ \text{dla} ~~~ t \in \Z \\ (iii)&~\gamma_X(r,s) = \gamma_X(r+t,s+t) ~~~ \text{dla} ~~~ t \in \Z \end{align} \)
Uwagi
  1. Powyższa definicja odnosi się do tak zwanej słabej stacjonarności, stacjonarności w szerszym sensie lub stacjonarności rzędu dwa. Ma ona zastosowanie najczęściej podczas analizy szeregów czasowych. Na tym kursie analizy szeregów czasowych będzie to podstawowa definicja jaką będziemy rozpatrywali.
  2. Jeżeli proces \( \{X_t, t \in \Z\}\ \) jest stacjonarny
\( \gamma_X(r+t,s+t) = \gamma_X(r-s,0) \!\)
często lepiej jest tak przedefiniować funkcję autokowariancji aby była funkcją tylko jednej zmiennej
\( \gamma_X(r-s,0) = \gamma_X(h,0) = cov(X_{t+h},X_h) = \gamma_X(h), \, \mbox{ gdzie } \, h = r - s, \, h,r,s \in \Z \)
Wtedy \( h \) możemy utożsamić z opóźnieniem w czasie dwóch zmiennych losowych \( X_{t+h}, X_t\ \). Wtedy jednoargumentowa funkcja autokorelacji procesu losowego \( \{X_t, t \in \Z\}\ \) zdefiniowana będzie jako
\( \rho_X(h) = \frac{\gamma_X(h)}{\gamma_X(0)}, \, h \in \Z. \)
  1. W powyższej definicji nie musimy ograniczać się do \( \Z \) jako dziedziny czasu. Jednak definicja ograniczająca do liczb całkowitych jest łatwiejsza, a w przypadku szeregów czasowych całkowicie wystarczająca, ponieważ na naszych zajęciach "czas" będzie zawsze indeksowany.

Stacjonarność w sensie ścisłym

Dla porównania podamy teraz definicję stacjonarności w sensie ścisłym.

Definicja
Proces stochastyczny (bądź szereg losowy) jest stacjonarny w sensie ścisłym, gdy zmienne losowe \(X(t)\ \) oraz \( X (t+\epsilon)\ \) mają te same rozkłady n-wymiarowe (rozkłady łączne)
\( f(x_1, t_1; x_2, t_2;...; x_n, t_n) = f(x_1, t_{1+\epsilon}; x_2, t_{2+\epsilon};...; x_n, t_{n+\epsilon})\ \)
dla dowolnych \(n\) i \(\epsilon\).
Uwagi
  1. \( f(x_1, t_1) = f(x_1, t_{1+\epsilon})\ \)
oznacza, że funkcja rozkładu nie zależy od czasu, t.j.
\( f(x_1, t_1) = f(x_1).\ \)
Wynika z tego natychmiast, że wartość oczekiwana jest również stałą funkcją czasu
\( \langle m(t) \rangle = m.\ \)
  1. Równość
\( f(x_1, t_1; x_2, t_2) = f(x_1, t_{1+\epsilon}; x_2, t_{2+\epsilon})\ \)
zachodzi tylko wtedy, gdy funkcja rozkładu zależy tylko od różnicy czasu
\( f(x_1, t_1; x_2, t_2) = f(x_1, x_2, t_1 - t_2).\ \)
Oznacza to, że w tym przypadku funkcja korelacyjna też zależy tylko od różnicy czasu
\( \langle x(t +\epsilon) x(t) \rangle = R(\epsilon).\ \)
  1. Definicja stacjonarności ścisłej obejmuje definicję stacjonarności w sensie szerokim. Odwrotne twierdzenie niekoniecznie jest prawdziwe.
Rysunek 3.1. Stacjonarny szereg ARMA(2,2). \({\color{red} EX_{n \in[2517,13232]} = 0.0018499, \sigma_{X_{n \in[2517,13232]}} = 0.33475}\) \({\color{red} \mu_{X_{n \in[2517,13232]}} = 0,000536, x^{max}_{X_{n \in[2517,13232]}} = 1.5311}\) \(EX_{n \in[12122,22803]} = 0.0044365, \sigma_{X_{n \in[12122,22803]}} = 0.33196\) \(\mu_{X_{n \in[12122,22803]}} = 0.0035154, x^{max}_{X_{n \in[12122,22803]}} = 1.3274\)

Jako, że kurs ten obejmuje zagadnienia związane również z wizualną analizą szeregów czasowych, to intuicyjnie można sobie ową stacjonarność silną wyobrazić następująco. Dwa dowolne fragmenty szeregu np. \(\{X_n, n \in(2517,13232)\}\ \) oraz \(\{X_n, n \in(12122,22803)\}\ \) będą miały takie same (z pewną dokładnością) charakterystyki statystyczne (średnie, wariancje, mediany, maksymalne amplitudy...).








Funkcja autokowariancji procesu stacjonarnego

Funkcja autokowariancji (lub opcjonalnie autokorelacji) stacjonarnego procesu stochastycznego, leżąca za każdym rozpatrywanym procesem stojącym za dowolnymi danymi wykazującymi losowość, jest podstawowym narzędziem służącym do identyfikacji procesu. Dzięki tej funkcji możemy znaleźć parametry opisujące typ i rodzaj procesu stojącego za daną realizacją procesu losowego. Techniki te poznamy w części traktującej o analizie szeregów czasowych.

Podstawowe własności

Jeżeli \( \gamma(\cdot) \ \) jest funkcją autokowariancji procesu stacjonarnego \( \{ X_t, t \in \Z \} \ \), to

(i) \( \gamma(0) \ge 0, \ \)
(ii) \( |\gamma(h)| \le \gamma(0), \ \) dla każdego \( h \in \Z \) oraz \(\gamma\) jest funkcją parzystą
(iii) \( \gamma(h) = \gamma(-h), \ \) dla każdego \( h \in \Z \)
własność (i) 
jest oczywista, ponieważ
\( \sigma_{X_t} \ge 0, ~~~~~(Var(X_t) \ge 0); \ \)
własność (ii) 
wynika z nierówności Cauchyego-Schwarza
\( |\gamma(h)| = |cov(X_{t+h}, X_t)| \le \sqrt{\sigma_{X_{t+h}}} \sqrt{\sigma_{X_{t}}}; \ \)
własność (iii) 
\( \gamma(-h) = \ \)
\( = E[(X_{t-h} - EX_{t-h})(X_{t} - EX_{t})] = \ \)
\( = E[X_{t-h}X_{t} - EX_{t-h}X_{t} - X_{t-h}EX_{t} + EX_{t-h}EX_{t}] = \ \)
\( = E(X_{t-h}X_{t}) - EX_{t-h}EX_{t} = \ \)
podstawiamy \( s = t -h \ \)
\( = E(X_{s}X_{s+h}) - EX_{s}EX_{s+h} = E(X_{s+h}X_{s}) - EX_{s+h}EX_{s} = \ \)
\( = \gamma(h). \ \)
Ilustracja własności (iii) \( \gamma(h) = \gamma(-h). \ \)

Funkcja określona nieujemnie

Definicja 
Funkcja przeprowadzająca liczby całkowite w rzeczywiste \(\kappa : \Z \to \R \ \) jest określona nieujemnie wtedy i tylko wtedy, gdy
\( \sum_{i,j=1}^{n} a_i \kappa (t_i - t_j) a_j \ge 0 \ \)
dla wszystkich dodatnich liczb całkowitych n oraz dla wszystkich wektorów \( \bar{a} = (a_1, a_2, \dots, a_n) \in \R^n \ \) oraz \( \bar{t} = (t_1, t_2, \dots, t_n) \in \Z^n \ \) wtedy i tylko wtedy, gdy
\( \sum_{i,j=1}^{n} a_i \kappa (i - j) a_j \ge 0 \ \)
dla wszystkich powyższych n i \(\bar{a}\).
Twierdzenie 
Funkcja rzeczywista określona na liczbach całkowitych jest autokowariancją stacjonarnego szeregu czasowego wtedy i tylko wtedy, gdy funkcja ta jest parzysta i nieujemnie określona. [bez dowodu]
Uwaga 1
Dla każdej funkcji autokowariancji \(\gamma (\cdot) \ \) istnieje stacjonarny Gaussowski szereg czasowy (proces losowy) który posiada funkcję autokowariancji.
Uwaga 2 
Czasem łatwiej jest zdefiniować proces stochastyczny oraz policzyć jego funkcję autokowariancji by sprawdzić czy jakaś funkcja jest określona nieujemnie.
Uwaga 3 
Funkcja autokorelacji spełnia wszystkie własności funkcji autokowaariancji oraz dodatkowo jedną trywialną \( \rho(0) = 1. \ \)

Funkcja autokowariancji próby

Zajmując się szeregami czasowymi zwykle dochodzimy do momentu, gdzie niezbędne staje się oszacowanie funkcji autokowariancji (autokorelacji) po to by uzyskać wiedzę na temat zmiennych procesu. Problem polega najczęściej na tym, że jedyne co jest dostępne to nie cały proces stochastyczny \( \{ x_t \} \) a pojedyncza jego realizacja \( \{ x_1, \dots, x_n \} \). Obliczenie tej funkcji to jeden z szeregu niezbędnych kroków w konstrukcji modelu posiadanych danych. Do oszacowania funkcji autokowariancji \( \gamma (\cdot) \ \) procesu stochastycznego stojącego za posiadanymi danymi losowymi używamy funkcji autokowariancji próby (z ang. sample autocovariance function - acf).

Definicja 
Funkcja autokowariancji próby dla realizacji \( \{ x_1, \dots, x_n \} \) zdefiniowana jest jako
\( \begin{align} (i) ~~~&\hat{\gamma}(h) := \frac{1}{n} \sum_{j=1}^{n-h}(x_{j+h} - \bar{x})(x_j - \bar{x}), 0 \le h \le n, -n < h < \le 0 \\ (ii) ~~~&\hat{\gamma}(h) = \hat{\gamma}(-h)\\ (iii)~~~&\bar{x} = \frac{1}{n} \sum_{j=1}^n x_j. \end{align} \)
Matlab / Octave
function [z] = asc_kow(x, stopien = 0)
 
  if ( isempty (x) || (!isvector(x)) ) 
    error ("kow: pierwszy agrument musi byc niepustym wektorem");
    return;
  endif
 
  %funkcja autokowariancji jest parzysta
  h = abs(stopien);
 
  %dlugosc wektora
  lx = length(x);
 
  if ( !isnumeric(stopien) || h >= lx)
    error ("kow: drugi argument musi byc liczba naturalna < dim(x)");
    return;
  endif
 
  %srednia
  srednia = mean(x);
 
  %suma
  korelacje = 0;
  for j=1:(lx-h)
    korelacje += x(j+h) * x(j);
  endfor
 
  z = (
    korelacje - 
    (sum(x(1:lx-h)) + sum(x(1+h:lx))) * srednia + 
    power(srednia,2) * (lx-h)
    ) / lx;
 
endfunction
Uwaga 1 
zastosowanie \((n)\) w mianowniku pierwszej równości gwarantuje, że funkcja
\( \Gamma_n := [\hat{\gamma} (i-j)]_{i,j=1}^n \ \)
jest określona nieujemnie.
Uwaga 2 
Funkcja autokorelacji próby określona jest przez funkcję autokowariancji próby poprzez
\( \hat{\rho}(h) = \frac{\hat{\gamma(h)}}{\hat{\gamma(0)}}, ~~~|h| < n. \ \)
Odpowiadająca macierz autokorelacji próby
\( \hat{R}_n := [\hat{\rho} (i-j)]_{i,j=1}^n \ \)
również określona jest nieujemnie.
Matlab / Octave
function [z] = asc_kor(x, stopien = 0)
 
  if ( isempty (x) || (!isvector(x)) ) 
    error ("kor: pierwszy agrument musi byc niepustym wektorem");
    return;
  endif
 
  if ( !isnumeric(stopien) )
    error ("kor: drugi argument musi byc liczba naturalna");
    return;
  endif
 
  %funkcja autorelacji jest parzysta
  h = abs(stopien);
 
  if ( h == 0)
    z =  1;
  else
    z = asc_kow(x,h) / asc_kow(x,0);
  endif
 
endfunction
Uwaga 3 
Z funkcji autokorelacji (właśnie ta funkcja wykorzystywana jest szerzej podczas badania szeregów czasowych) można wyczytać kilka informacji
a) Bardziej zmienne dane wejściowe dają w wyniku bardziej zmienną funkcję autokorelacji.

b) Dodatnia wartość funkcji autokorelacji dla przesunięcia 1 odpowiada takiej sytuacji, że kolejna obserwacja będzie leżeć po tej samej stronie co średnia.
c) Ujemna wartość funkcji autokorelacji dla przesunięcia 1 odpowiada takiej sytuacji, że kolejna obserwacja będzie leżeć po przeciwnej stronie co średnia.
d) "Periodyczność" funkcji autokorelacyjnej często odpowiada sezonowości w danych wejściowych.

Uwaga 4 
Funkcję autokorelacji (autokowariancji) możemy policzyć dla dowolnych danych (dowolnego szeregu czasowego), stacjonarnego ale i szeregu wykazującego trend czy sezonowość (bądź jedno i drugie). Funkcja autokorelacji policzona dla danych z trendem wykaże stosunkowo wolne zanikanie wraz z wzrostem h (odległości korelacji). Funkcja autokorelacji dla danych z czynnikiem periodycznym wykaże się również przebiegiem z okresowością. W pewnym sensie funkcja autokorelacji może nam służyć za test stacjonarności szeregu czasowego.

Funkcja autokorelacji częściowej

Funkcja autokorelacji częściowej (fac, partial autocorrelation function - pacf) podobnie jak zwykła funkcja autokorelacji zawiera bardzo istotne informacje dotyczące zależności pomiędzy zmiennymi w dla procesów stacjonarnych. Podobnie jak i funkcja autokorelacji, fac (pacf) zależy jedynie od własności statystycznych drugiego rzędu. Funkcję autokorelacji częściowej rzędu k można rozpatrywać jako korelację pomiędzy \( X_1 \ \) oraz \( X_{k+1} \ \) ale dopasowaną do danych pośrednich.

Definicja 
Funkcja autokorelacji częściowej \( \alpha (\cdot) \ \) procesu stacjonarnego zdefiniowana jest jako

\[ \begin{align} \alpha(1) &= corr(X_2, X_1) = \rho(1), \\ \alpha(k) &= corr(X_{k+1} - P_{\bar{sp}(1,X_2,\dots,X_k)}X_{k+1}, X_{1} - P_{\bar{sp}(1,X_2,\dots,X_k)}X_{1})),~k \ge 2, \end{align} \]

gdzie \(P_{\bar{sp}(1,X_2,\dots,X_k)}X_i \ \) to najlepszy liniowy predyktor X_i w oparciu o podzbiór \( \{ 1,X_2,\dots,X_k \} \ \).

Predyktor liniowy a średnia warunkowa

Nie do końca wgłębiając się w szczegóły teoretyczne \( P_\mathcal{M}X \ \) nazywamy operatorem rzutowym, rzutującym zmienną X na podprzestrzeń \(\mathcal{M}\).

Definicja 
Warunkowa wartość oczekiwana \(E_\mathcal{M}X \ \).

Jeżeli \(\mathcal{M}\) to podprzestrzeń rzeczywistej przestrzeni Hilberta[2] \(L^2 (\Omega, \mathcal{F}, P), \mathcal{M} \in L^2 \ \) [produktem tej podprzestrzeni jest wartość średnia \(\langle X,Y \rangle = E(X,Y) \ \)] oraz jeżeli \(X \in L^2 \) to oczekiwana wartość warunkowa zmiennej X

\( E_\mathcal{M}X = P_\mathcal{M}(X) = P_\mathcal{M}(X) \ \)

zdefiniowana jest poprzez rzutowanie zmiennej X na podprzestrzeń \(\mathcal{M}\).

Uwaga 
\( E_\mathcal{M}X \ \) ma sporo ciekawych własności, ale skoro nie są nam one potrzebne do toku zajęć to tymczasem je pominiemy...

Na podstawie twierdzenia o rzutowaniu oczekiwana wartość warunkowa \( E_{\mathcal{M}(Z_1,\dots, Z_n)}(X) \ \) to najlepszy średniokwadratowy predyktor zmiennej X w zbiorze \( \mathcal{M}(Z_1,\dots, Z_n) \ \) czyli najlepsza funkcja zmiennych \( \{ Z_1,\dots, Z_n \} \ \) przewidująca (prognozująca) X. Rzeczywiste oszacowanie tej funkcji jest z reguły bardzo trudne, dlatego bardzo często korzysta się z dopasowanej funkcji liniowej do zbioru danych \( \{ Z_1,\dots, Z_n \} \ \)

\( P_{\bar{sp}(1, Z_1,\dots, Z_n)}(X) = \sum_{i=0}^n \alpha_i Z_i, Z_0 = 1 \ \) gdzie \( \{ \alpha_0, \dots, \alpha_n \} \ \) spełniają warunek
\( \langle \sum_{i=0}^n \alpha_i Z_i, Z_j \rangle = \langle X, Z_j \rangle, ~~~j=0,1,\dots,n \ \)

lub też

\( \sum_{i=0}^n \alpha_i E(Z_iZ_j) = E(XZ_j), ~~~j=0,1,\dots,n \ \).

Jeżeli rozwiązanie równania na \( P_{\bar{sp}}X \) istnieje to istnieje też taki zbiór parametrów \( \{ \alpha_0, \dots, \alpha_n \} \ \) które minimalizują błąd szacunku X.

Uwaga 
Jeżeli stacjonarny proces ma średnią zerową to
\(P_{\bar{sp}(1,X_2,\dots,X_k)}(\cdot) = P_{\bar{sp}(X_2,\dots,X_k)}(\cdot) \ \)

Metoda znajdowania najlepszego liniowego predyktora

Definicja 
(Najlepszy liniowy predyktor X na podstawie zbioru Z)

Jeżeli \(X \in L^2 \ \) oraz \( \{ Z_\lambda, \lambda \in \Lambda \} \in L^2 \ \) to wtedy najlepszy liniowy predyktor X w oparciu o \(Z_\lambda \) jest zdefiniowany jako element \( \bar{sp} \{ Z_\lambda, \lambda \in \Lambda \} \ \) o najmniejszym średnio-kwadratowym odchyleniu od X.

Uwaga 
Twierdzenie rzutowe mówi nam, że takim predyktorem jest po prostu \( P_{\bar{sp} \{ Z_\lambda, \lambda \in \Lambda \}}X. \ \)

Alternatywna definicja funkcji autokorelacji częściowej

Niech \( \{ X_t \} \ \) będzie procesem stacjonarnym o średniej zerowej posiadającym funkcję autokowariancji \( \gamma (\cdot) \ \) taką, że \( \gamma (h) \to 0 \ \) gdy \( h \to \infty. \ \) Załóżmy, że \(\phi_{kj}, j=1,\dots,k,~k=1,2,\dots\) to współczynniki rozwinięcia

\( P_{\bar{sp}(X_1,\dots,X_k)}X_{k+1} = \sum_{j=1}^k \phi_{kj} X_{k+1-j}. \)

Wykorzystując równanie

\( \langle X_{k+1} - P_{\bar{sp}(X_1,\dots,X_k)}X_{k+1}, X_j \rangle = 0,~~~j=k,\dots,1 \ \)

możemy zdefiniować macierz

\( \begin{bmatrix} \rho(0) & \rho(1) & \rho(2) & \cdots & \rho(k-1) \\ \rho(1) & \rho(0) & \rho(1) & \cdots & \rho(k-2) \\ \rho(2) & \rho(1) & \rho(0) & \cdots & \rho(k-3) \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ \rho(k-1) & \rho(k-2) & \rho(k-3) & \cdots & \rho(0) \end{bmatrix} \begin{bmatrix} \phi_{k1} \\ \phi_{k2} \\ \phi_{k3} \\ \vdots \\ \phi_{kk}\end{bmatrix} = \begin{bmatrix} \rho(1) \\ \rho(2) \\ \rho(3) \\ \vdots \\ \rho(k)\end{bmatrix}. \)
Definicja 
Funkcja autokorelacji częściowej \( \alpha (k) \ \) rzędu k dla procesu \( \{ X \} \ \) dana jest wzorem
\( \alpha (k) = \phi_{kk},~~~k \ge 1. \ \)
Definicja 
Funkcja autokorelacji częściowej próby \( \hat{\alpha} (k) \ \) rzędu k dla realizacji procesu \( \{ X_1, X_2, \dots, X_n \} \ \) (zakładamy, że \( X_i \ne X_j\) dla i,j) dana jest wzorem
\( \hat{\alpha} (k) = \hat{\phi}_{kk},~~~1 \le k < n. \ \)

Poniżej podajemy skrypt Matlab / GNU Octave oparty właśnie na definicji alternatywnej.


Kod MS3 (Matlab / Octave)
function [z] = asc_fakc(x, stopien = 1)
 
  if ( isempty (x) || (!isvector(x)) ) 
    error ("autokor_czesc: pierwszy agrument musi byc niepustym wektorem");
    return;
  endif
 
  %dlugosc wektora
  lx = length(x);
 
  %funkcja autorelacji jest parzysta
  h = abs(stopien);
 
  if ( !isnumeric(stopien) || h == 0 || h > lx)
    error ("autokor_czesc: drugi argument musi byc liczba naturalna 0 < |h| <= dim(wektor)\n");
    return;
  endif
 
  G = [];
  for i=0:h
    G(i+1) = asc_kow(x,i);
  endfor
  G = G';
 
  M = [];
  for i=0:h-1
    for j=0:h-1
      M (i+1,j+1) = G(abs(i-j)+1);
    endfor
  endfor
  %M
 
  %M * P = G => P = M^-1 * G
  P = inv(M) * G(2:length(G));
 
  % \alpha = P(h)
  z = P(length(P));
 
endfunction

Przypisy 
  1. Słowo nad zbiorem skończonym (alfabetem) S to ciąg elementów zbioru S
  2. rzetelnego czytelnika odsyłamy do rozdziału 5 pozycji Matematyka w fizyce klasycznej i kwantowej, F.W. Byron, R.W. Fuller, (PWN 1974), a nieco mniej cierpliwego do artykułu na Wikipedii