Analiza Szeregów Czasowych/Procesy stochastyczne

Z Skrypty dla studentów Ekonofizyki UPGOW

(Różnice między wersjami)
(Alternatywna definicja funkcji autokorelacji częściowej)
(Alternatywna definicja funkcji autokorelacji częściowej)
Linia 380: Linia 380:
</source>
</source>
|}
|}
-
[[Grafika:pacf_dane.png]]
+
 
-
[[Grafika:pacf_kowariancja.png]]
+
-
[[Grafika:pacf_korelacja.png]]
+
-
[[Grafika:pacf_pacf.png]]
+
<center><gallery widths=240 heights=180 perrow=2 caption="Przykład szeregu czasowego z obliczonymi numerycznie funkcjami: autokowariancji, autokorelacji i autokorelacji częściowej">
<center><gallery widths=240 heights=180 perrow=2 caption="Przykład szeregu czasowego z obliczonymi numerycznie funkcjami: autokowariancji, autokorelacji i autokorelacji częściowej">
Plik:pacf_dane.png|Szereg czasowy.
Plik:pacf_dane.png|Szereg czasowy.

Wersja z 23:41, 22 gru 2010

Analiza Szeregów Czasowych

Spis treści

Elementy teorii prawdopodobieństwa

Procesy stochastyczne

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[1] \(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) \ \)
  1. 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

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