Simulasi SARIMA (4,1,4)(4,1,4)12

Inisialisasi parameter

set.seed(123)
# Panjang data 
n <- 1200

# Inisialisasi vektor
# y adalah data sebelum differencing
y <- rep(0, n + 60) 
# wt = xt - xt-1 - xt-12 + xt-13
w <- rep(0, n + 60)
# White noise
e <- rnorm(n + 60)  
# Inisialisasi parameter
phi <- c(0.5, -0.3, 0.1, -0.05)
Phi <- c(0.6, -0.2, 0.1, -0.05)
theta <- c(0.4, 0.2, -0.1, -0.05)
Theta <- c(0.5, 0.3, -0.1, 0.05)

Simulasi ARIMA(4,1,4)(4,1,4)12

for (t in 61:(n + 60)) {
  
  # Komponen AR
  AR_sum <- 0
  for (i in 1:4) {
    AR_sum <- AR_sum + phi[i] * w[t - i]
  }
  for (j in 1:4) {
    AR_sum <- AR_sum + Phi[j] * w[t - 12 * j]
  }
  for (i in 1:4) {
    for (j in 1:4) {
      AR_sum <- AR_sum - (phi[i] * Phi[j] * w[t - i - 12 * j])
    }
  }
  
  # Komponen MA
  MA_sum <- 0
  for (i in 1:4) {
    MA_sum <- MA_sum - theta[i] * e[t - i]
  }
  for (j in 1:4) {
    MA_sum <- MA_sum - Theta[j] * e[t - 12 * j]
  }
  for (i in 1:4) {
    for (j in 1:4) {
      MA_sum <- MA_sum + theta[i] * Theta[j] * e[t - i - 12 * j]
    }
  }
  
  # wt
  w[t] <- AR_sum + e[t] + MA_sum
  
  # yt = yt-1 + yt-12 - yt-13 + wt
  y[t] <- y[t-1] + y[t-12] - y[t-13] + w[t]
}

Buang buffer awal (data ke 1-60) dan buat jadi seasonal 12

y <- ts(y[61:(n + 60)], frequency = 12)
w <- ts(w[61:(n + 60)],frequency = 12)

Plot dari data sebelum dan sesudah differencing

ts.plot(y)

ts.plot(w)

PACF dan ACF dari data sebelum dan sesudah differencing

par(mfrow = c(1,2))
acf(y)
pacf(y)

acf(w)
pacf(w)