\[ y_t = a_t + \Theta_1 a_{t-12} + \Theta_2 a_{t-24} + \Theta_3 a_{t-36} \]
#Data
set.seed(123)
at = rnorm(1000, mean = 0, sd = 2)
theta1 = 0.6
theta2 = -0.4
theta3 = 0.3

n <- length(at)
y <- numeric(n)
# Inisialisasi awal (Agar indeks negatif tidak error)
y[1:36] <- 0
#model
for (t in 37:n) {
  y[t] <- at[t] +
    theta1 * at[t - 12] +
    theta2 * at[t - 24] +
    theta3 * at[t - 36]}
#plot
ts.plot(y)

#ACF dan PACF
par(mfrow = c(1, 2))
acf(y)
pacf(y)