Time Series: Midterm

library(stats)
library(MARSS)
library(forecast)
library(datasets)

1 AR

Mari kita mulai dengan mensimulasikan beberapa model AR(1) dan membandingkan perilaku mereka. Pertama, mari kita pilih model dengan koefisien AR yang kontras. Ingat bahwa agar model AR(1) stasioner, \(\varphi < |1|\), jadi kita akan mencoba 0,1 dan 0,8. Kami akan kembali mengatur benih angka acak sehingga kami akan mendapatkan jawaban yang sama.

set.seed(343)
## list description for AR(1) model with small coef
AR.sm <- list(order = c(1, 0, 0), ar = 0.1, sd = 0.1)
## list description for AR(1) model with large coef
AR.lg <- list(order = c(1, 0, 0), ar = 0.8, sd = 0.1)
## simulate AR(1)
AR1.sm <- arima.sim(n = 50, model = AR.sm)
AR1.lg <- arima.sim(n = 50, model = AR.lg)

## setup plot region
par(mfrow = c(1, 2))
## get y-limits for common plots
ylm <- c(min(AR1.sm, AR1.lg), max(AR1.sm, AR1.lg))
## plot the ts
plot.ts(AR1.sm, ylim = ylm, ylab = expression(italic(x)[italic(t)]), 
    main = expression(paste(phi, " = 0.1")))
plot.ts(AR1.lg, ylim = ylm, ylab = expression(italic(x)[italic(t)]), 
    main = expression(paste(phi, " = 0.9")))

Sepertinya seri waktu dengan koefisien AR yang lebih kecil lebih “berombak” dan tampaknya tetap lebih dekat ke 0 sedangkan seri waktu dengan koefisien AR yang lebih besar tampaknya lebih berkeliaran. Ingatlah bahwa sebagai koefisien dalam model AR(1) ke 0, model mendekati urutan WN, yang stasioner dalam rata-rata dan varians. Namun, ketika koefisien pergi ke 1, model mendekati jalan acak, yang tidak stasioner baik dalam rata-rata atau varians. Selanjutnya, mari kita hasilkan dua model AR(1) yang memiliki koefisien magnitudo yang sama, tetapi tanda-tanda yang berlawanan, dan bandingkan hasilnya.

set.seed(343)
## list description for AR(1) model with small coef
AR.pos <- list(order = c(1, 0, 0), ar = 0.5, sd = 0.1)
## list description for AR(1) model with large coef
AR.neg <- list(order = c(1, 0, 0), ar = -0.5, sd = 0.1)
## simulate AR(1)
AR1.pos <- arima.sim(n = 50, model = AR.pos)
AR1.neg <- arima.sim(n = 50, model = AR.neg)

## setup plot region
par(mfrow = c(1, 2))
## get y-limits for common plots
ylm <- c(min(AR1.pos, AR1.neg), max(AR1.pos, AR1.neg))
## plot the ts
plot.ts(AR1.pos, ylim = ylm, ylab = expression(italic(x)[italic(t)]), 
    main = expression(paste(phi[1], " = 0.5")))
plot.ts(AR1.neg, ylab = expression(italic(x)[italic(t)]), main = expression(paste(phi[1], 
    " = -0.5")))

Sekarang tampak seperti kedua seri waktu bervariasi di sekitar rata-rata dengan jumlah yang sama, tetapi model dengan koefisien negatif menghasilkan seri waktu yang jauh lebih “sawtooth”.Ternyata setiap model AR(1) dengan $ -1 < < 0$ akan menunjukkan osilasi 2 poin yang Anda lihat di sini.

2 MA

Kita bisa mensimulasikan MA \((q)\) proses seperti yang kita lakukan untuk AR \((p)\) proses menggunakan \(arima.sim()\). Berikut adalah 3 yang berbeda dengan \(θ\)’s yang kontras:

set.seed(343)
## list description for MA(1) model with small coef
MA.sm <- list(order = c(0, 0, 1), ma = 0.2, sd = 0.1)
## list description for MA(1) model with large coef
MA.lg <- list(order = c(0, 0, 1), ma = 0.7, sd = 0.1)
## list description for MA(1) model with large coef
MA.neg <- list(order = c(0, 0, 1), ma = -0.4, sd = 0.1)
## simulate MA(1)
MA1.sm <- arima.sim(n = 50, model = MA.sm)
MA1.lg <- arima.sim(n = 50, model = MA.lg)
MA1.neg <- arima.sim(n = 50, model = MA.neg)

## setup plot region
par(mfrow = c(1, 3))
## plot the ts
plot.ts(MA1.sm, ylab = expression(italic(x)[italic(t)]), main = expression(paste(theta, 
    " = 0.2")))
plot.ts(MA1.lg, ylab = expression(italic(x)[italic(t)]), main = expression(paste(theta, 
    " = 0.7")))
plot.ts(MA1.neg, ylab = expression(italic(x)[italic(t)]), main = expression(paste(theta, 
    " = -0.4")))

Berbeda dengan proses AR(1), model MA(1) tidak menunjukkan perilaku yang berbeda secara radikal dengan perubahan \(θ\) .Ini seharusnya tidak terlalu mengejutkan mengingat bahwa mereka hanya kombinasi linear dari white noise.

3 ARMA

Mari kita lihat contoh cara kerja \(arima()\). Pertama kita akan mensimulasikan model ARMA (1,1) dan kemudian memperkirakan parameter untuk melihat seberapa baik kita dapat memulihkannya. Selain itu, kita akan menambahkan konstanta untuk membuat rata-rata bukan nol, yang menurut \(arima()\) sebagai penyadapan dalam outputnya.

set.seed(343)
## ARMA(1,1) description for arim.sim()
ARMA11 <- list(order = c(1, 0, 1), ar = c(0.4), ma = c(0.4))
## mean of process
mu <- 5
## simulated process (+ mean)
ARMA.sim <- arima.sim(n = 100, model = ARMA11) + mu
## estimate parameters
arima(x = ARMA.sim, order = c(1, 0, 1))

## 
## Call:
## arima(x = ARMA.sim, order = c(1, 0, 1))
## 
## Coefficients:
##          ar1     ma1  intercept
##       0.0443  0.5819     4.9044
## s.e.  0.1739  0.1466     0.1373
## 
## sigma^2 estimated as 0.6933:  log likelihood = -123.82,  aic = 255.63