Autoregresné procesy

Číslovanie príkladov sa odvoláva na slajdy, v ktorých sú príklady aj podrobnejšie popísané. Za každým príkladom nasledujú v slajdoch cvičenia.

Príklad 1: Analýza zadaného autoregresného procesu v R-ku

(a) Overenie stacionarity AR procesu

library(fArma) 

armaRoots(c(1.4,-0.85))

##       re      im   dist
## 1 0.8235  0.7059 1.0847
## 2 0.8235 -0.7059 1.0847

(b) Využitie R-ka na počítanie s komplexnými číslami

x <- complex(real=0.8235, imaginary=0.7059)
# vseobecne:
# korene <- armaRoots(c(1.4,-0.85))
# x <- complex(real= korene$re[1], imaginary= korene$im[1])
Arg(x)

## [1] 0.7086563

Opakovanie z prednášky: Ako sa z tohto vypočíta perióda?

(c) Výpočet ACF a PACF

acf1 <- ARMAacf(ar=c(1.4,-0.85), lag.max=10)  # ACF
pacf1 <- ARMAacf(ar=c(1.4,-0.85),lag.max=10, pacf="true") # PACF

acf1

##           0           1           2           3           4           5 
##  1.00000000  0.75675676  0.20945946 -0.35000000 -0.66804054 -0.63775676 
##           6           7           8           9          10 
## -0.32502500  0.08705824  0.39815279  0.48341440  0.33835029

pacf1

##  [1]  7.567568e-01 -8.500000e-01 -4.353501e-15  3.134088e-16  5.821933e-16
##  [6] -1.561284e-16 -6.281330e-16  3.716454e-16  2.414129e-16 -5.009978e-16

Kvôli numerickému výpočtu sú nenulové čísla aj tam, kde by mali byť nuly (ale sú veľmi blízke nule), uzitočná môže byť funkcia zapsmall.

zapsmall(pacf1)

##  [1]  0.7567568 -0.8500000  0.0000000  0.0000000  0.0000000  0.0000000
##  [7]  0.0000000  0.0000000  0.0000000  0.0000000

Graficky:

# alebo mozeme zrusit nulovy lag (vzdy sa rovna 1, nedava nam ziadnu informaciu)
# barplot(acf1[1:length(acf1)])
barplot(acf1)

barplot(pacf1)

(d) Vygenerovanie realizácie procesu

set.seed(123)
x <- arima.sim(n=50, list(ar = c(1.4,-0.85)), sd=0.1)
plot(x)

Príklad 2: Odhadovanie parametrov AR procesu

Využijeme dáta uložené v premennej x, ktoré boli vygenerované z procesu so známymi parametrami.

Budeme používať knižnicu astsa

library(astsa)

(a) Výberova ACF a PACF

Máme dve možnosti:

• štandardnú funkciu acf, resp. s parametrom type=partial

acf(x) # aj s nulovym lagom

acf(x, type="partial")

• funkciu acf2 z knižnice astsa

acf2(x) # bez nuloveho lagu, ACF aj PACF

##         ACF  PACF
##  [1,]  0.74  0.74
##  [2,]  0.25 -0.69
##  [3,] -0.18  0.10
##  [4,] -0.37 -0.01
##  [5,] -0.31  0.06
##  [6,] -0.11  0.01
##  [7,]  0.11  0.06
##  [8,]  0.21 -0.06
##  [9,]  0.21  0.14
## [10,]  0.15  0.05
## [11,]  0.10  0.09
## [12,]  0.05 -0.11
## [13,] -0.03 -0.04
## [14,] -0.16 -0.23
## [15,] -0.26  0.10
## [16,] -0.22  0.10
## [17,] -0.02  0.20
## [18,]  0.24  0.02

(b) Odhadovanie a testovanie AR modelu

Funkcia sarima z knižnice astsa.

sarima(x,2,0,0, details=FALSE) # ar(2) model

## $fit
## 
## Call:
## stats::arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, 
##     Q), period = S), xreg = xmean, include.mean = FALSE, optim.control = list(trace = trc, 
##     REPORT = 1, reltol = tol))
## 
## Coefficients:
##          ar1      ar2   xmean
##       1.3571  -0.7737  0.0049
## s.e.  0.0877   0.0883  0.0282
## 
## sigma^2 estimated as 0.006735:  log likelihood = 52.71,  aic = -97.42
## 
## $AIC
## [1] -3.880417
## 
## $AICc
## [1] -3.822639
## 
## $BIC
## [1] -4.765696

(c) Konštrukcia predikcií

Funkcia sarima.for z knižnice astsa

sarima.for(x,30,2,0,0)

Príklad 3: Práca s reálnymi dátami

• Načítanie dát

x <- read.table("http://www.iam.fmph.uniba.sk/institute/stehlikova/cr14/data/kava.txt")

Z vektora spravíme časový rad pomocou funkcie ts, ako parameter zadávame frekvenciu dát (1 = ročné, 4 = kvartálne, 12 = mesačné) a začiatočný/koncový čas v tvare vektora.

x <- ts(x, frequency=1, start=c(1910))
plot(x)

• Pracovať musíme so stacionárnymi dátami.

Zdá sa, že v dátach je trend, preto budeme pracovať s diferenciami, aby sme modelovali stacionárny časový rad. Typickým znakom radu, ktorý treba zdiferencovať je velmi pomaly klesajúca výberová ACF. Na prednáške (a potom aj na druhom cvičení pri PC) bude aj konkrétny štatistický test, ktorý slúži na zisťovanie, či dáta diferencovať.

• Výberová ACF a PACF

acf2(diff(x))

• Odhadnutie modelu a jeho testovanie

Odhadnite vhodný AR(p) model

sarima(x,p,1,0) # AR(p) model pre prve diferencie procesu 'x'

• Predikcie

Spravte predikcie na nasledujucich 10 rokov.