1 Impor Library

Library yang digunakan pada analisis ini meliputi beberapa package untuk pemodelan deret waktu.

library(readxl)
library(tseries)
library(forecast)
library(tsoutliers)
library(lmtest)
library(TSA)

2 Impor Data

Data yang digunakan pada penerapan model ARIMA Outlier dan ARIMA intervensi adalah data Indeks Harga Saham Gabungan (IHSG) dari 12 Desember 2023 hingga 9 April 2025. Indeks Harga Saham Gabungan (IHSG) merupakan salah satu indeks pasar saham yang digunakan oleh Bursa Efek Indonesia. Nilai IHSG mengalami fluktuasi setiap hari sehingga berpotensi memiliki nilai pencilan (outlier). Selain itu, terdapat kebijakan pemerintah berupa pengesahan RUU TNI pada tanggal 20 Maret 2025 yang diduga memengaruhi pergerakan IHSG.

data_ihsg <- read_excel("D:/Semester 6/Metode Peramalan/Ghina Azziqra - ihsg 2.xlsx")

#pendefinisan variabel
ihsg <- as.numeric(
  data_ihsg$`Harga Saham Penutup`
)

tanggal <- data_ihsg$Tanggal

data.ihsg <- data.frame(
  tanggal,
  ihsg
)

3 Eksplorasi Data Deret Waktu

Pada tahap awal dilakukan eksplorasi data untuk melihat pola umum pergerakan IHSG dari waktu ke waktu. Eksplorasi ini bertujuan untuk mengidentifikasi adanya tren, fluktuasi, serta indikasi awal kestasioneran data.

plot(tanggal,
     ihsg,
     type = "l",
     col = "skyblue",
     lwd = 2,
     main = "Plot Time Series IHSG",
     xlab = "Tanggal",
     ylab = "Harga Penutup")

Berdasarkan plot time series tersebut, terlihat bahwa data IHSG mengalami fluktuasi yang cukup dinamis sepanjang periode pengamatan. Pada awal periode tahun 2024, nilai IHSG cenderung bergerak relatif stabil dengan kecenderungan meningkat secara bertahap. Pola data menunjukkan adanya tren menurun pada bagian akhir periode pengamatan, yang ditunjukkan dengan penurunan nilai IHSG hingga mendekati nilai terendah pada akhir data. Dari plot tersebut juga terlihat bahwa variabilitas data cenderung berubah sepanjang waktu, sehingga secara visual data diduga belum stasioner terhadap rata-rata.

4 Pemodelan ARIMA

Pemodelan ARIMA dilakukan untuk memperoleh model awal yang dapat menggambarkan pola data deret waktu IHSG sebelum mempertimbangkan pengaruh outlier dan intervensi.

4.1 Uji Stasioneritas

Stasioneritas merupakan asumsi penting dalam data deret waktu. Uji stasioneritas terhadap rata-rata dilakukan menggunakan Augmented Dickey-Fuller (ADF) Test. Hipotesis yang digunakan adalah sebagai berikut:

H0 : Data IHSG mengandung akar unit sehingga tidak stasioner terhadap rata-rata.

H1 : Data IHSG tidak mengandung akar unit sehingga stasioner terhadap rata-rata.

Kriteria pengujian: Tolak H0 apabila nilai p-value < 0,05.

adf.test(ihsg)

    Augmented Dickey-Fuller Test

data:  ihsg
Dickey-Fuller = -0.59149, Lag order = 6, p-value = 0.9773
alternative hypothesis: stationary

Nilai p-value tersebut > 0.05 sehingga terima H0 yang menyatakan bahwa data mengandung akar unit. Dengan demikian, data IHSG pada kondisi awal belum stasioner terhadap rata-rata. untuk mengatasi ketidakstasioneran tersebut, dilakukan proses differencing orde pertama, sebagai berikut.

\[ \Delta Y_t = Y_t - Y_{t-1} \] di mana:

\(Y_t\) : nilai deret waktu pada waktu ke-\(t\)
\(Y_{t-1}\) : nilai deret waktu pada waktu ke-\((t-1)\)
\(\Delta Y_t\) : hasil differencing orde pertama

diff_ihsg <- diff(ihsg)
adf.test(diff_ihsg)

    Augmented Dickey-Fuller Test

data:  diff_ihsg
Dickey-Fuller = -7.703, Lag order = 6, p-value = 0.01
alternative hypothesis: stationary

Berdasarkan ADF test differencing orde 1, terlihat bahwa p-value < 0.05. Hal ini menunjukkan bahwa H0 ditolak, sehingga data setelah differencing pertama telah stasioner terhadap rata-rata.

plot(diff_ihsg,
     type = "l",
     col = "maroon",
     main = "Plot Differencing IHSG")

4.2 Identifikasi Model ARIMA

Acf(diff_ihsg)

Pacf(diff_ihsg)

Berdasarkan plot ACF dari data hasil differencing, terlihat bahwa beberapa lag melewati batas signifikansi, khususnya pada lag 2, 7, 10, dan 19, yang menunjukkan adanya autokorelasi yang signifikan pada lag tersebut.

Sementara itu, plot PACF menunjukkan adanya beberapa lag yang signifikan, terutama pada lag 2, 10, dan 19, yang mengindikasikan adanya komponen autoregressive pada lag tersebut.

Diperoleh beberapa kandidat model ARIMA yang memungkinkan untuk digunakan dalam pemodelan data IHSG, yaitu:

• ARIMA([2,10,19],1,0) dengan parameter autoregressive tertentu yang dibatasi (restricted model), dan

• ARIMA(0,1,[2,7,10,19]) sebagai model pembanding dengan komponen moving average.

4.3 Estimasi Model ARIMA

model.restricted1 <- Arima(
  ihsg,
  order = c(19,1,0),
  fixed = c(
    0, NA,     # lag 2
    rep(0,7),
    NA,        # lag 10
    rep(0,8),
    NA         # lag 19
  )
)

coeftest(model.restricted1)

z test of coefficients:

      Estimate Std. Error z value Pr(>|z|)   
ar2  -0.123405   0.059480 -2.0747 0.038011 * 
ar10  0.170637   0.062387  2.7351 0.006236 **
ar19 -0.191442   0.063179 -3.0301 0.002444 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Berdasarkan hasil diatas, parameter autoregressive pada lag 2, 10, dan 19 memiliki nilai p-value < 0,05, sehingga dapat disimpulkan bahwa parameter tersebut signifikan dan berpengaruh dalam model.

model.restricted2 <- Arima(
  ihsg,
  order = c(0,1,19),
  fixed = c(
    0, NA,
    rep(0,4),
    NA,
    rep(0,2),
    NA,
    rep(0,8),
    NA
  )
)

coeftest(model.restricted2)

z test of coefficients:

      Estimate Std. Error z value Pr(>|z|)   
ma2  -0.096872   0.054646 -1.7727 0.076275 . 
ma7  -0.127676   0.055611 -2.2959 0.021683 * 
ma10  0.120953   0.058196  2.0784 0.037676 * 
ma19 -0.177836   0.060364 -2.9461 0.003219 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Berdasarkan estimasi model pembanding dengan komponen moving average , parameter pada lag 7, 10, dan 19 memiliki nilai p-value < 0,05 sehingga signifikan, sedangkan parameter pada lag 2 tidak signifikan karena p-value > 0.05.

Berdasarkan kedua model, model yang dapat digunakan adalah ARIMA([2,10,19],1,0) karena memiliki parameter yang signifikan.

4.4 Diagnostik Residual Model

Uji independensi residual dilakukan menggunakan uji Ljung-Box. Hipotesis yang digunakan dalam pengujian ini adalah sebagai berikut:

\[ H_0 : \rho_1 = \rho_2 = \rho_3 = \cdots = \rho_k = 0 \] \[ H_1 : \text{Minimal terdapat satu } \rho_k \neq 0 \]

dengan:

\(\rho_k\) : autokorelasi residual pada lag ke-\(k\)

Kriteria pengujian: Tolak \(H_0\) apabila nilai p-value < \(\alpha\) (0,05).

resid.arima1 <- residuals(model.restricted1)

Box.test(
  resid.arima1,
  lag = 10,
  type = "Ljung-Box"
)

    Box-Ljung test

data:  resid.arima1
X-squared = 11.319, df = 10, p-value = 0.3332

Uji normalitas residual dilakukan menggunakan uji Jarque-Bera. Hipotesis yang digunakan adalah sebagai berikut:

\[ H_0 : \text{Residual berdistribusi normal} \]

\[ H_1 : \text{Residual tidak berdistribusi normal} \]

jarque.bera.test(resid.arima1)

    Jarque Bera Test

data:  resid.arima1
X-squared = 223.39, df = 2, p-value < 2.2e-16

Berdasarkan hasil diagnostik residual tersebut, dapat disimpulkan bahwa residual model tidak memiliki autokorelasi karena p-value > 0.05, sehingga residual dapat dikatakan bersifat white noise. Namun, residual tidak berdistribusi normal.

Meskipun residual belum memenuhi asumsi normalitas, model ARIMA([2,10,19],1,0) tetap digunakan sebagai model dasar karena telah memenuhi asumsi tidak adanya autokorelasi residual dan memiliki parameter yang signifikan secara statistik.

5 Pemodelan ARIMA dengan Outlier

Pemodelan ARIMA dengan mempertimbangkan outlier dilakukan untuk meningkatkan akurasi model dengan memperhitungkan adanya pengamatan ekstrem yang dapat mempengaruhi estimasi parameter model.

5.1 Deteksi Outlier

Deteksi outlier dilakukan dengan mempertimbangkan beberapa jenis outlier, yaitu Additive Outlier (AO), Innovational Outlier (IO), Level Shift (IS), dan Temporary Change (TC).

prediksi.arima <- ihsg - resid.arima1

outlier <- tso(
  prediksi.arima,
  types = c('AO','IO','IS','TC')
)

print(outlier)
Series: prediksi.arima 
Regression with ARIMA(0,1,0) errors 

Coefficients:
          TC151     AO292      AO311
      -251.9815  225.5150  -400.6111
s.e.    61.6684   47.2975    66.8888

sigma^2 = 4518:  log likelihood = -1742.81
AIC=3493.62   AICc=3493.75   BIC=3508.57

Outliers:
  type ind time coefhat  tstat
1   TC 151  151  -252.0 -4.086
2   AO 292  292   225.5  4.768
3   AO 311  311  -400.6 -5.989

Berdasarkan hasil deteksi tersebut, terdapat tiga outlier yang terdiri dari satu Temporary Change (TC) pada waktu ke-151 serta dua Additive Outlier (AO) pada waktu ke-292 dan ke-311. Nilai t yang relatif besar menunjukkan bahwa outlier tersebut signifikan dan berpengaruh terhadap pola data IHSG.

5.2 Estimasi Model ARIMA dengan Outlier

Setelah outlier teridentifikasi, selanjutnya dibentuk model ARIMA dengan memasukkan variabel outlier sebagai variabel eksogen (xreg). Model yang digunakan adalah ARIMA([2,10,19],1,0) dengan variabel outlier berupa TC151, AO292, dan AO311.

n <- length(ihsg)

# AO (pulse)
ao1 <- rep(0, n); ao1[292] <- 1
ao2 <- rep(0, n); ao2[311] <- 1

# TC (decay)
tc1 <- rep(0, n)

for(i in 151:n){
  tc1[i] <- 0.7^(i-151)
}

xreg <- cbind(tc1, ao1, ao2)
dim(xreg)
[1] 311   3
#fixed vector
fixed_vec <- rep(0, 22)  
# 19 AR + 3 outlier

# Lag AR signifikan
fixed_vec[2] <- NA
fixed_vec[10] <- NA
fixed_vec[19] <- NA

# Outlier (3 parameter)
fixed_vec[20:22] <- NA

length(fixed_vec)
[1] 22
sum(is.na(fixed_vec))
[1] 6
#model
arima.tc <- arima(
  ihsg,
  order = c(19,1,0),
  xreg = xreg,
  method = "ML",
  transform.pars = FALSE,
  include.mean = FALSE,
  fixed = fixed_vec
)

arima.tc

Call:
arima(x = ihsg, order = c(19, 1, 0), xreg = xreg, include.mean = FALSE, transform.pars = FALSE, 
    fixed = fixed_vec, method = "ML")

Coefficients:
      ar1      ar2  ar3  ar4  ar5  ar6  ar7  ar8  ar9    ar10  ar11  ar12  ar13
        0  -0.1689    0    0    0    0    0    0    0  0.1421     0     0     0
s.e.    0   0.0614    0    0    0    0    0    0    0  0.0638     0     0     0
      ar14  ar15  ar16  ar17  ar18     ar19      tc1        ao1       ao2
         0     0     0     0     0  -0.1866  17.0703  -151.9409  -13.1014
s.e.     0     0     0     0     0   0.0631  64.1072    51.2453   71.8610

sigma^2 estimated as 5060:  log likelihood = -1762.38,  aic = 3536.77
coeftest(arima.tc)

z test of coefficients:

        Estimate  Std. Error z value Pr(>|z|)   
ar2    -0.168882    0.061405 -2.7503 0.005954 **
ar10    0.142069    0.063750  2.2285 0.025845 * 
ar19   -0.186591    0.063060 -2.9589 0.003087 **
tc1    17.070311   64.107161  0.2663 0.790025   
ao1  -151.940943   51.245316 -2.9650 0.003027 **
ao2   -13.101404   71.860975 -0.1823 0.855335   
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Berdasarkan hasil tersebut, diperoleh bahwa parameter autoregressive pada lag 2, 10, dan 19 memiliki nilai p-value < 0,05, sehingga parameter tersebut signifikan secara statistik. Selain itu, variabel outlier AO(292) juga menunjukkan signifikansi statistik dengan nilai p-value sebesar 0,003027.

Namun demikian, variabel outlier TC(151) dan AO(311) memiliki nilai p-value > 0,05, sehingga tidak signifikan secara statistik. Meskipun demikian, variabel tersebut tetap dipertahankan dalam model karena secara keseluruhan model menunjukkan performa yang lebih baik dibandingkan model tanpa mempertimbangkan outlier.

5.3 Diagnostik Residual Model ARIMA dengan Outlier

Dilakukan pengujian diagnostik residual untuk memastikan bahwa residual model memenuhi asumsi model ARIMA. Pengujian dilakukan menggunakan uji Ljung-Box dan Jarque-Bera.

resid.tc <- residuals(arima.tc)

Box.test(
  resid.tc,
  lag = 10,
  type = "Ljung-Box"
)

    Box-Ljung test

data:  resid.tc
X-squared = 11.046, df = 10, p-value = 0.3539
jarque.bera.test(resid.tc)

    Jarque Bera Test

data:  resid.tc
X-squared = 253.9, df = 2, p-value < 2.2e-16

Berdasarkan hasil diagnostik residual tersebut, dapat disimpulkan bahwa residual model tidak memiliki autokorelasi, sehingga residual dapat dikatakan bersifat white noise. Namun, residual tidak berdistribusi normal.

Meskipun demikian, model tetap dapat digunakan karena telah memenuhi asumsi independensi residual dan mampu mengakomodasi pengaruh outlier pada data. Model ARIMA([2,10,19],1,0) dengan outlier tetap dapat digunakan.

6 Pemodelan ARIMA Intervensi

Pemodelan ARIMA intervensi dilakukan untuk menganalisis pengaruh suatu kejadian tertentu terhadap pola data IHSG.

6.1 Pemodelan Pra Intervensi

Model pra intervensi dibentuk menggunakan data sebelum terjadinya intervensi, yaitu sebanyak 150 pengamatan pertama.

#data pra intervensi
ihsg.pra <- ts(ihsg[1:150])
ihsg.pra
Time Series:
Start = 1 
End = 150 
Frequency = 1 
  [1] 7125.31 7075.34 7176.02 7190.99 7119.52 7187.85 7219.67 7209.62 7237.52
 [10] 7237.52 7245.92 7303.89 7272.80 7323.59 7279.09 7359.76 7350.62 7283.58
 [19] 7200.20 7227.30 7219.96 7241.14 7224.00 7242.79 7200.63 7252.97 7227.40
 [28] 7247.93 7256.23 7227.82 7178.04 7137.09 7157.17 7192.22 7207.94 7201.70
 [37] 7238.79 7198.62 7247.41 7235.15 7297.67 7209.74 7303.28 7335.54 7296.70
 [46] 7352.60 7349.02 7339.64 7295.10 7283.82 7285.32 7328.64 7316.11 7311.91
 [55] 7276.75 7247.46 7329.80 7373.96 7381.91 7421.21 7433.31 7328.05 7302.45
 [64] 7336.75 7331.13 7338.35 7350.15 7377.76 7365.66 7310.09 7288.81 7205.06
 [73] 7236.98 7166.84 7254.40 7286.88 7164.81 7130.84 7166.81 7087.32 7073.82
 [82] 7110.81 7174.53 7155.29 7036.08 7155.78 7234.20 7234.20 7117.42 7134.72
 [91] 7135.89 7123.61 7088.79 7099.26 7083.76 7179.83 7246.70 7317.24 7266.69
[100] 7186.04 7222.38 7176.42 7253.63 7140.23 7034.14 6970.74 7036.19 7099.31
[109] 6947.67 6974.90 6897.95 6921.55 6855.69 6850.10 6831.56 6734.83 6726.92
[118] 6819.32 6879.98 6889.17 6882.70 6905.64 6967.95 7063.58 7139.63 7125.14
[127] 7196.75 7220.89 7253.37 7250.98 7269.80 7287.04 7300.41 7327.58 7278.86
[136] 7224.29 7224.22 7321.07 7294.50 7321.98 7313.86 7262.76 7240.28 7288.17
[145] 7288.90 7241.86 7255.76 7325.98 7308.12 7059.65

6.1.1 Uji Stasioneritas Data Pra Intervensi

adf.test(ihsg.pra)

    Augmented Dickey-Fuller Test

data:  ihsg.pra
Dickey-Fuller = -2.6061, Lag order = 5, p-value = 0.3239
alternative hypothesis: stationary
diff.pra <- diff(ihsg.pra)
adf.test(diff.pra)

    Augmented Dickey-Fuller Test

data:  diff.pra
Dickey-Fuller = -4.3277, Lag order = 5, p-value = 0.01
alternative hypothesis: stationary

Hasil uji stasioneritas menggunakan Augmented Dickey-Fuller (ADF) menunjukkan bahwa data pra intervensi belum stasioner pada kondisi awal dengan nilai p-value sebesar 0,3239. Setelah dilakukan differencing orde pertama, diperoleh nilai p-value sebesar 0,01, sehingga data hasil differencing orde 1 telah stasioner terhadap rata-rata.

6.1.2 Identifikasi Model

Acf(diff.pra)

Pacf(diff.pra)

Berdasarkan plot tersebut, terlihat bahwa lag 18 signifikan pada plot PACF. Maka, diperoleh kandidat model ARIMA([18],1,0) sebagai model pra intervensi.

fixed_vec_pra <- rep(0,18)
fixed_vec_pra[18] <- NA

arima.pra <- arima(
  ihsg.pra,
  order = c(18,1,0),
  fixed = fixed_vec_pra,
  method = "ML",
  include.mean = FALSE,
  transform.pars = FALSE
)

arima.pra

Call:
arima(x = ihsg.pra, order = c(18, 1, 0), include.mean = FALSE, transform.pars = FALSE, 
    fixed = fixed_vec_pra, method = "ML")

Coefficients:
      ar1  ar2  ar3  ar4  ar5  ar6  ar7  ar8  ar9  ar10  ar11  ar12  ar13  ar14
        0    0    0    0    0    0    0    0    0     0     0     0     0     0
s.e.    0    0    0    0    0    0    0    0    0     0     0     0     0     0
      ar15  ar16  ar17     ar18
         0     0     0  -0.2162
s.e.     0     0     0   0.0877

sigma^2 estimated as 3097:  log likelihood = -810.7,  aic = 1623.4

Berdasarkan hasil diatas, parameter autoregressive pada lag 18 memiliki nilai estimasi sebesar −0,2162, yang menunjukkan adanya pengaruh autokorelasi pada lag tersebut terhadap data pra intervensi.

6.1.3 Diagnostik Residual Model Pra Intervensi

resid.pra <- residuals(arima.pra)

Box.test(
  resid.pra,
  lag = 10,
  type = "Ljung-Box"
)

    Box-Ljung test

data:  resid.pra
X-squared = 10.729, df = 10, p-value = 0.379
jarque.bera.test(resid.pra)

    Jarque Bera Test

data:  resid.pra
X-squared = 32.087, df = 2, p-value = 1.077e-07

Pada diagnostik menunjukkan bahwa residual model tidak memiliki autokorelasi dan dapat dikatakan bersifat white noise. Namun demikian, hasil uji Jarque-Bera menunjukkan bahwa residual tidak memenuhi asumsi normalitas.

Model pra intervensi yang telah terbentuk selanjutnya digunakan untuk menghasilkan nilai residual yang digunakan dalam identifikasi waktu terjadinya intervensi.

6.2 Estimasi Model ARIMA Intervensi

Estimasi model ARIMA intervensi dilakukan menggunakan model dasar ARIMA([18],1,0) dengan memasukkan variabel intervensi sebagai variabel transfer function.

intervensi <- c(
  rep(0,150),
  1,
  rep(0,length(ihsg)-151)
)

model.intervensi <- arimax(
  ihsg,
  order = c(18,1,0),
  fixed = c(
    rep(0,17),   
    NA,         
    NA, NA        
  ),
  
  xtransf = intervensi,
  transfer = list(c(1,0)),
  
  method = "ML"
)

model.intervensi

Call:
arimax(x = ihsg, order = c(18, 1, 0), fixed = c(rep(0, 17), NA, NA, NA), method = "ML", 
    xtransf = intervensi, transfer = list(c(1, 0)))

Coefficients:
      ar1  ar2  ar3  ar4  ar5  ar6  ar7  ar8  ar9  ar10  ar11  ar12  ar13  ar14
        0    0    0    0    0    0    0    0    0     0     0     0     0     0
s.e.    0    0    0    0    0    0    0    0    0     0     0     0     0     0
      ar15  ar16  ar17     ar18  T1-AR1    T1-MA0
         0     0     0  -0.0048  0.0000  -17.3954
s.e.     0     0     0   0.0653  0.0151   53.0092

sigma^2 estimated as 5618:  log likelihood = -1778.11,  aic = 3562.22
coeftest(model.intervensi)

z test of coefficients:

          Estimate  Std. Error z value Pr(>|z|)
ar18   -4.7705e-03  6.5321e-02 -0.0730   0.9418
T1-AR1  7.1066e-13  1.5121e-02  0.0000   1.0000
T1-MA0 -1.7395e+01  5.3009e+01 -0.3282   0.7428

Berdasarkan hasil diatas, parameter intervensi memiliki nilai p-value > 0,05, sehingga parameter tersebut tidak signifikan secara statistik. Hal ini menunjukkan bahwa pengaruh intervensi yang dimasukkan dalam model belum memberikan kontribusi yang signifikan terhadap perubahan pola data IHSG.

6.3 Diagnostik Residual Model ARIMA Intervensi

Pengujian diagnostik residual dilakukan untuk memastikan bahwa residual model memenuhi asumsi independensi setelah model ARIMA intervensi terbentuk.

resid.intervensi <- residuals(model.intervensi)

Box.test(
  resid.intervensi,
  lag = 4,
  type = "Ljung-Box"
)

    Box-Ljung test

data:  resid.intervensi
X-squared = 5.0366, df = 4, p-value = 0.2836
jarque.bera.test(resid.intervensi)

    Jarque Bera Test

data:  resid.intervensi
X-squared = 733.84, df = 2, p-value < 2.2e-16

Berdasarkan hasil diagnostik tersebut, residual model bersifat white noise dan tidak berdistribusi normal.

Dengan demikian, model ARIMA intervensi yang telah terbentuk selanjutnya akan dibandingkan dengan model lainnya pada tahap pemilihan model terbaik.

7 Pemilihan Model Terbaik

Pemilihan model terbaik dilakukan dengan membandingkan beberapa model yang telah dibentuk sebelumnya, yaitu model ARIMA dasar, ARIMA dengan outlier, dan ARIMA intervensi. Perbandingan dilakukan menggunakan beberapa kriteria evaluasi model, yaitu Akaike Information Criterion (AIC), Root Mean Square Error (RMSE), dan Mean Absolute Percentage Error (MAPE).

# Fungsi RMSE
rmse <- function(actual, predicted){
  sqrt(mean((actual - predicted)^2))
}

# Fungsi MAPE
mape <- function(actual, predicted){
  mean(abs((actual - predicted)/actual))*100
}

# Prediksi dari masing-masing model
pred_arima <- fitted(model.restricted1)

pred_outlier <- fitted(arima.tc)

pred_intervensi <- fitted(model.intervensi)


rmse_arima <- rmse(ihsg, pred_arima)
rmse_outlier <- rmse(ihsg, pred_outlier)
rmse_intervensi <- rmse(ihsg, pred_intervensi)

mape_arima <- mape(ihsg, pred_arima)
mape_outlier <- mape(ihsg, pred_outlier)
mape_intervensi <- mape(ihsg, pred_intervensi)

#ringkasan
ringkasan_model <- data.frame(
  Model = c(
    "ARIMA",
    "ARIMA + Outlier",
    "ARIMA + Intervensi"
  ),
  
  AIC = c(
    AIC(model.restricted1),
    AIC(arima.tc),
    AIC(model.intervensi)
  ),
  
  RMSE = c(
    rmse_arima,
    rmse_outlier,
    rmse_intervensi
  ),
  
  MAPE = c(
    mape_arima,
    mape_outlier,
    mape_intervensi
  )
)

ringkasan_model
               Model      AIC     RMSE      MAPE
1              ARIMA 3541.839 72.05631 0.7560614
2    ARIMA + Outlier 3538.765 71.02092 0.7405569
3 ARIMA + Intervensi 3564.216 74.83586 0.7673116

Beradasarkan nilai tersebut, nilai AIC terkecil diperoleh pada model ARIMA + Outlier yang menunjukkan bahwa model tersebut memiliki tingkat kesesuaian terbaik dengan data dibandingkan model lainnya. Selain itu, model ARIMA dengan outlier juga memiliki nilai RMSE dan MAPE yang menunjukkan bahwa model tersebut memiliki tingkat kesalahan prediksi yang lebih kecil dibandingkan model ARIMA dasar maupun model ARIMA intervensi.

Sebaliknya, model ARIMA intervensi memiliki nilai AIC, RMSE, dan MAPE yang paling besar, sehingga model tersebut kurang optimal dalam menggambarkan pola data dibandingkan model lainnya.

Maka dari itu, berdasarkan kriteria AIC, RMSE, dan MAPE, model ARIMA([2,10,19],1,0) dengan outlier dipilih sebagai model terbaik yang akan digunakan dalam tahap peramalan data IHSG pada periode selanjutnya.

8 Peramalan

Setelah diperoleh model terbaik berdasarkan kriteria pemilihan model, selanjutnya dilakukan proses peramalan menggunakan model ARIMA([2,10,19],1,0) dengan outlier.

Dalam penelitian ini, hasil peramalan dinyatakan dalam bentuk periode waktu, bukan tanggal kalender. Hal ini dikarenakan data IHSG merupakan data deret waktu yang hanya tersedia pada hari perdagangan saham, sehingga terdapat beberapa tanggal yang tidak memiliki pengamatan, seperti pada akhir pekan dan hari libur. Oleh karena itu, penggunaan periode waktu lebih tepat untuk menjaga konsistensi interval pengamatan dalam analisis deret waktu.

Peramalan dilakukan untuk 30 periode ke depan menggunakan model ARIMA dengan variabel outlier yang telah ditentukan sebelumnya. Hasil peramalan menghasilkan nilai ramalan (point forecast) beserta batas interval kepercayaan sebesar 80% dan 95%.

p <- 19
d <- 1
q <- 0

# fixed vector
fixed_arima <- c(
  0, NA, rep(0,7),   # ar2
  NA,                 # ar10
  rep(0,8),
  NA,                 # ar19
  NA, NA, NA          # tc1, ao1, ao2
)

# xreg training
xreg_train <- cbind(
  tc1 = tc1,
  ao1 = ao1,
  ao2 = ao2
)

# refit model
arima.tc2 <- Arima(
  ihsg,
  order = c(p,d,q),
  xreg = xreg_train,
  include.mean = FALSE,
  fixed = fixed_arima,
  transform.pars = FALSE,
  method = "ML"
)

summary(arima.tc2)
Series: ihsg 
Regression with ARIMA(19,1,0) errors 

Coefficients:
      ar1      ar2  ar3  ar4  ar5  ar6  ar7  ar8  ar9    ar10  ar11  ar12  ar13
        0  -0.1689    0    0    0    0    0    0    0  0.1421     0     0     0
s.e.    0   0.0614    0    0    0    0    0    0    0  0.0638     0     0     0
      ar14  ar15  ar16  ar17  ar18     ar19      tc1        ao1       ao2
         0     0     0     0     0  -0.1866  17.0703  -151.9409  -13.1014
s.e.     0     0     0     0     0   0.0631  64.1072    51.2453   71.8610

sigma^2 = 5160:  log likelihood = -1762.38
AIC=3538.77   AICc=3539.14   BIC=3564.92

Training set error measures:
                    ME     RMSE      MAE         MPE      MAPE      MASE
Training set -4.050069 71.02092 52.35529 -0.06783506 0.7405569 0.9664456
                   ACF1
Training set 0.03777222
h <- 30
last_index <- length(tc1)

# lanjutkan decay TC
tc1_future <- 0.7^((last_index+1):(last_index+h) - 151)

# AO masa depan = 0
ao1_future <- rep(0,h)
ao2_future <- rep(0,h)

xreg_future <- cbind(
  tc1 = tc1_future,
  ao1 = ao1_future,
  ao2 = ao2_future
)

colnames(xreg_future)
[1] "tc1" "ao1" "ao2"
colnames(xreg_future) <- colnames(xreg_train)

######## FORECAST ########
forecast.outlier <- forecast(
  arima.tc2,
  h = h,
  xreg = xreg_future
)

forecast.outlier
    Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
312       6032.841 5940.782 6124.900 5892.049 6173.633
313       6031.792 5901.601 6161.983 5832.682 6230.902
314       6029.612 5878.603 6180.621 5798.663 6260.560
315       6019.296 5850.010 6188.582 5760.395 6278.197
316       6015.659 5828.789 6202.530 5729.866 6301.453
317       6005.731 5802.795 6208.668 5695.367 6316.096
318       6043.267 5825.607 6260.927 5710.384 6376.150
319       6074.970 5843.521 6306.419 5720.999 6428.941
320       6003.689 5759.204 6248.175 5629.781 6377.598
321       6042.577 5785.715 6299.438 5649.741 6435.412
322       6045.496 5772.711 6318.281 5628.308 6462.684
323       6025.717 5737.888 6313.545 5585.521 6465.912
324       6047.956 5747.145 6348.767 5587.905 6508.007
325       6067.923 5754.666 6381.179 5588.838 6547.007
326       6049.768 5724.240 6375.296 5551.915 6547.620
327       6000.812 5663.458 6338.166 5484.874 6516.750
328       6002.071 5653.356 6350.786 5468.758 6535.385
329       6110.841 5751.123 6470.558 5560.700 6660.981
330       6103.309 5732.903 6473.715 5536.822 6669.796
331       6080.809 5703.629 6457.988 5503.962 6657.655
332       6082.691 5698.510 6466.872 5495.137 6670.245
333       6084.088 5691.906 6476.269 5484.298 6683.877
334       6088.854 5689.016 6488.692 5477.355 6700.353
335       6092.134 5685.062 6499.205 5469.571 6714.696
336       6090.602 5676.363 6504.840 5457.079 6724.125
337       6076.089 5654.746 6497.432 5431.700 6720.478
338       6070.611 5642.297 6498.926 5415.561 6725.662
339       6101.815 5666.653 6536.977 5436.293 6767.337
340       6094.414 5652.508 6536.321 5418.576 6770.252
341       6085.403 5637.662 6533.144 5400.642 6770.164
plot(
  forecast.outlier,
  main = "Peramalan IHSG 30 Hari ke Depan (ARIMA + Outlier)",
  xlab = "Waktu",
  ylab = "IHSG",
  col = "red",
  type = "l"
)

Berdasarkan grafik tersebut, nilai IHSG pada periode mendatang diperkirakan mengalami pergerakan yang relatif stabil dengan kecenderungan fluktuasi kecil di sekitar nilai 6000–6100. Area bayangan yang semakin melebar menunjukkan bahwa ketidakpastian ramalan meningkat pada periode yang lebih jauh.

Secara umum, model ARIMA([2,10,19],1,0) dengan outlier mampu memberikan gambaran pergerakan IHSG pada periode mendatang dengan tingkat kesalahan yang relatif kecil berdasarkan nilai RMSE sebesar 71,02092 dan MAPE sebesar 74,056%.