Подсетимо се да постоје строго и слабо стационарни модели. Овде се бавимо слабо стационарним моделима, тј. моделима који имају константно очекивање и дисперзију у времену и \(cov(x_{t},x_{t+k})\) зависи само од кашњења \(|k|\). Стационарност је идеализација модела и ми увек настојимо да “очистимо” податке од трендова и сезонских компоненти да бисмо могли да им доделимо стационарни модел.
\(MA(q)\) модел
\(MA\) модел реда \(q\) је модел који се може записати у облику:
\[x_{t}= e_{t}-\theta_{1}e_{t-1}-\theta_{2}e_{t-2}-\dots-\theta_{q}e_{t-q},\] где је \(e_{t}\) бели шум са очекивањем 0 и дисперзијом \(\sigma^2_w\) (nајчешће је то Гаусов бели шум). У \(MA\) моделима тренутна вредност серије зависи од претходних грешака. Претходна једначина се може записати као: \[x_t = (1-\theta_1L - \theta_2L^2-\dots-\theta_qL^q)e_t = \phi_q(L)e_t,\] где је \(\phi_q(L)\) полином реда \(q\), а \(\theta_1,\theta_2,\dots,\theta_q\) су параметри покретних просека.
Слично као код ауторегресионих модела, моделима покретних просека можемо доделити карактеристичну једначину: \[g^q-\theta_1g^{q-1}-\dots-\theta_q=0.\]
Директно из једначине модела можемо видети следеће особине:
1. Средња вредност: \(E(x_{t})=0\)
2. Дисперзија: \(D(x_{t})=(1+\theta_{1}^2+\dots+\theta_{q}^2)\sigma^2\)
3. Аутокорелациона функција: \[\rho(k) =
\begin{cases}1, \; k=0\\
\frac{\sum_{i=0}^{q-k}\theta_i\theta_{i+k}}{\sum_{i=0}^q \theta_i^2}, \;
k=1, \dots,q\\
0, \; k>q\end{cases}\]
Приметимо да је \(MA(q)\) специјалан случај линеарног модела са коначно много коефицијената различитих од нуле, када узмемо да су коефицијенти у том моделу \(\theta_{j}=0\) за \(j>q\).
\(MA(q)\) модел је увек стационаран (за свако \(q\) и све коефицијенте \(\theta\)). Објашњење за то је некорелисаност \(x_{t}\) и свих пре \(x_{t-q}\), јер се не појављују исти шумови.
Линеарни модел је заправо MA(\(\infty\)): \[x_{t}=\displaystyle\sum_{i=0}^{\infty} \theta_{i}e_{t-i},\] где ред \(\displaystyle\sum_{i=0}^{\infty} \theta_{i}\) апсолутно конвергира.
Подсетимо се \(AR\) модела. Ауторегресиона серија реда \(q\) (ознакa \(AR(q)\)) је облика: \[x_{t}= e_{t}+\theta_{1}x_{t-1}+\theta_{2}x_{t-2}+\dots+\theta_{q}x_{t-q}\]
Видимо да њена тренутна вреност зависи од њених претходних вредности. Сада ћемо видети да су те зависности (од грешака и вредности модела) на неки начин исте. Следећа теорема нам даје веома важно својство стационарних процеса.
Волдова теорема о декомпозицији Сваки стационарни модел се може приказати у облику \(MA(\infty)\) модела.
Закључујемо да се сваки стационарни \(AR(q)\) модел може конвертовати у \(MA(\infty)\) и обрнуто, сваки \(MA(q)\) у \(AR(\infty)\) модел ако важи услов инвертибилности. \(MA(q)\) је инвертибилан и има јединствен инвертибилнии модел ако су нуле полинома са коефицијентима \(\theta\) ван јединичног корена (Више о овоме теоријски на предавањима).
- Сада се питамо када онда користимо \(MA\) моделе.
Не постоји правило, а одговор на ово питање налазимо у искуству. Свакако што се тиче стационарних процеса ако нам је један модел одговарајући, онда ће и други бити скоро сигурно.
Посматрајмо сада најједноставнији модел покретних просека првог реда \(MA(1)\): \[x_{t}= e_{t}-\theta_{1}e_{t-1},\] или другачије записано \[x_t = (1-\theta_1L)e_t.\]
Аутокорелациона функција овог модела дата са: \[\rho_k = \frac{\gamma_k}{\gamma_0} = \begin{cases}1, & k=0 \\ \frac{\theta_{1}}{1+\theta_{1}^2}, & k=1\\ 0, & k>1.\end{cases}\]
Можемо приметити да \(\rho_1= \frac{\theta_{1}}{1+\theta_{1}^2} \in \left[-\frac{1}{2},\frac{1}{2}\right]\), као и да аутокорелациони коефицијенти не одређују једнозначно \(MA(1)\) модел, јер се добија иста вредност \(\rho_1\) за моделе са параметрима \(\theta\) и \(\frac{1}{\theta}\) (види услов инвертибилности који нам каже који модел користимо).
Теоријски ACF ћемо израчунати у R-у:
rho <- function(lag, theta){#lag je kašnjenje u ACF, theta koeficijenti modela sa dodatom jedinicom
q <- length(theta) - 1
if(lag > q){
ACF <- 0
}
else{
s1 <- 0
s2 <- 0
for (i in 1:(q-lag+1)) {
s1 <- s1 + theta[i] * theta[i+lag]}
for (i in 1:(q+1)) {
s2 <- s2 + theta[i]^2}
ACF <- s1 / s2
}
return(ACF)
}Нацртајмо сада добијене теоријске вредности ACF за један \(MA(3)\) модел са параметрима \(\beta_1=0.7\), \(\beta_2=0.5\) и \(\beta_3=0.2\) и за један \(MA(3)\) модел са параметрима \(\beta_1=-0.7\), \(\beta_2=0.5\) и \(\beta_3=-0.2\):
## [1] -0.25+1.561249i -2.00+0.000000i -0.25-1.561249irho_k <- rep(0, 10) #ovde ćemo smestiti vrednosti funkcije
for (k in 1:10){
rho_k[k] <- rho(k, theta)}
rho_k## [1] 0.6460674 0.3595506 0.1123596 0.0000000 0.0000000 0.0000000 0.0000000
## [8] 0.0000000 0.0000000 0.0000000## Warning: package 'ggplot2' was built under R version 4.2.3data<-as.data.frame(cbind(lag=0:10, rho_k=c(1,rho_k)))
ggplot(data,aes(x=lag,y=rho_k))+
geom_segment(aes(y=0, x=lag,xend=lag,yend=rho_k),col="red")+
geom_point()
## [1] -0.25+1.561249i -2.00+0.000000i -0.25-1.561249irho_k <- rep(0, 10) #ovde ćemo smestiti vrednosti funkcije
for (k in 1:10){
rho_k[k] <- rho(k, theta)}
rho_k## [1] -0.6460674 0.3595506 -0.1123596 0.0000000 0.0000000 0.0000000
## [7] 0.0000000 0.0000000 0.0000000 0.0000000library(ggplot2)
data<-as.data.frame(cbind(lag=0:10, rho_k=c(1,rho_k)))
ggplot(data,aes(x=lag,y=rho_k))+
geom_segment(aes(y=0, x=lag,xend=lag,yend=rho_k),col="red")+
geom_point()
Примери
Симулирајмо један \(MA(3)\) модел.
set.seed(1)
b <- c(0.8, 0.6, 0.4)
x <- w <- rnorm(1000)
for (t in 4:1000){
for(j in 1:3) x[t] <- x[t] + b[j] * w[t - j]
}
plot(x, type = 'l')
Функција у R-у која симулира \(MA\)
модел је arima.sim. Њени обавезни елементи су листа са
елементима:
1. order = (0,0,q), где је q ред MA модела
2. maкоефицијенти \(\theta\) у моделу
3. n представља број серија који желимо да
симулирамо.
Генеришимо сада два \(MA(1)\) модела:
set.seed(5)
ma1<-arima.sim(list(order=c(0,0,1), ma=0.5), n=100)
ma2<-arima.sim(list(order=c(0,0,1), ma=-0.5), n=100)
par(mfrow = c(1,2))
plot.ts(ma1,main="MA(1) sa koef. 0.5")
plot.ts(ma2,main="MA(1) sa koef. -0.5")
Нацртајмо сада графике узорачких (оцењених) функција аутокорелације.
Плава испрекидана линија нам говори да су већ од \(lag=2\) вредности ван критичне области, па
за њих прихватамо хипотезу да је аутокорелацја нула што је очекивано,
јер су ово MA(1) модели. Видимо да има неких вредности за \(lag\geq2\) које су ван плавог испрекиданог
појаса, али ако оне само мало излазе то занемарујемо, јер тада би опет
прихватили хипотезу ако само мало повећамо праг значајности.
Упоредимо сада са теоријским вредностима за ACF:
theta1 <- c(1,0.5) #ponavljamo kod od gore
rho_k1 <- rep(0, 50)
for (k in 1:50){
rho_k1[k] <- rho(k, theta1)}
theta2 <- c(1,-0.5)
rho_k2 <- rep(0, 50)
for (k in 1:50){
rho_k2[k] <- rho(k, theta2)}
par(mfrow = c(1,2))
plot(0:50,c(1,rho_k1), type = "h", xlab = "lag", ylab = "ma1")
plot(0:50,c(1,rho_k2), type = "h",xlab = "lag", ylab = "ma2")
ma1<-arima.sim(list(order=c(0,0,1), ma=0.5), n=100)
dist(rbind(rho_k1, acf(ma1,lag.max = 50,plot = F)$acf[1:50]), method = "euclidian")## rho_k1
## 1.068124ma1<-arima.sim(list(order=c(0,0,1), ma=0.5), n=1000)
dist(rbind(rho_k1, acf(ma1,lag.max = 50,plot = F)$acf[1:50]), method = "euclidian")## rho_k1
## 0.7399258Можемо приметити да је са повећањем обима узорка узорачка ACF ближа
теоријској.
Фитовање МА модела
Функција arima се користи за фитовање \(MA(q)\) модела за неке податке тако што је
параметар order=(0,0,q) где је \(q\) ред модела. За разлику од функције
ar функција arima иницијално не одузима средњу
вредности, па ћемо имати и оцењени intercept. Једна од
повратних вредности су коефицијенти који се оцењују методом најмањих
квадрата. У R-у се нумеричким методама тражи вредност коефицијената које
минимизују:
\[ S(\hat\theta_{1},\dots,\hat\theta_{q})=\displaystyle\sum_{t=1}^{n} {\big(x_{t}-(\hat\theta_{1}\hat{e_{t-1}}+\dots+\hat\theta_{q}\hat{e_{t-q}})\big)^2},\]
где су оцењене грешке заправо резидуали који се добијали итеративном методом, а за њихову почетну вредност узете су нуле.
##
## Call:
## arima(x = x, order = c(0, 0, 3))
##
## Coefficients:
## ma1 ma2 ma3 intercept
## 0.7898 0.5665 0.3959 -0.0322
## s.e. 0.0307 0.0351 0.0320 0.0898
##
## sigma^2 estimated as 1.068: log likelihood = -1452.41, aic = 2914.83- Како бирамо ред \(q\)?
За \(MA(q)\) модел теоријске вредности ACF су једнаке 0 за задршке веће од \(q\). Зато би узорачке требало да буду блиске 0. Ако серију уклапамо у \(MA\) модел, има смисла за ред модела изабрати ред оне узорачке аутокорелације која је последња значајна. Некад и код серија тачно симулираних из неког \(MA(q)\) модела, имаћемо значајне корелације и реда већег од \(q\). Битно је и колико те вредности испадају из граница. Ако је само мало, онда са мало већим нивоом значајности вредности узорачке ACF упадају у област прихватања хипотезе.
Посматрајмо сада базу са прошлог часа - индекси потрошачких цена.
data<-read.table("CPIAUCSL.txt")
data.ts<-ts(data$V2,start = c(1947,1), end = c(2015,3), frequency = 12 )
plot(data.ts)
Када нацртамо корелограм резидуала видимо да модел није баш задовољавајући, јер не добијамо бели шум, али није ни тако лош што је последица тога да смо га ми на прошлом часу фитовали линеарно тј. за MA(\(\infty\)).
## Warning in arima(data.ts, order = c(0, 0, 10)): possible convergence problem:
## optim gave code = 1##
## Call:
## arima(x = data.ts, order = c(0, 0, 10))
##
## Coefficients:
## ma1 ma2 ma3 ma4 ma5 ma6 ma7 ma8
## 3.3654 6.8150 10.4034 12.9768 13.5686 11.9834 8.8217 5.1968
## s.e. 0.0320 0.0972 0.1847 0.2643 0.3071 0.2963 0.2377 0.1544
## ma9 ma10 intercept
## 2.2301 0.5169 100.6168
## s.e. 0.0767 0.0256 3.2114
##
## sigma^2 estimated as 1.446: log likelihood = -1322.62, aic = 2669.25Напомена: У R-у постоји функција ma
која изравњава временску серију методом покретних просека. Не треба
мешати са моделовањем и функцијом arima!
Задатак
Уклопити временску серију pounds_nz у одговарајући \(MA\) модел.
Нека је \(\{x_t: t=1,\dots,n\}\) стационарна временска серија за коју је \(E(x_t)=\mu\), \(D(x_t)=\sigma^2\) и \(Cor(x_t,x_{t+k})=\rho_k\).
Израчунати \(D(\bar{x})\) ако је \(\{x_t\}\) модел \(MA(1)\) дат са \(x_t=e_t+\frac{1}{2}e_{t-1}.\)
Израчунати \(D(\bar{x})\) ако је \(\{x_t\}\) модел \(MA(1)\) дат са \(x_t=e_t-\frac{1}{2}e_{t-1}.\)
Упоредити претходне резултате са дисперзијом узорачке средине добијене за бели шум у случају \(\rho_k=0\) \((k>0)\). Који од ова три модела би имао најтачнију оцену за \(\mu\) засновану на дисперзији узорачке средине?