2022-11-16

Peramalan Harga Penutupan Emas (Juli 2012-2022) dengan Metode Pemulusan dan Pemodelan

Kelompok 8 Paralel 1

Nama Anggota :

1. Afris Setiya Intan Amanda     (G1401201018)

2. Alwi Sidiq             (G1401201020)

3. Mufti Habibie Alayubi       (G1401201027)

4. Cahya Ireno Anugrah      (G1401201034)

5. Siti Aisyah            (G1401201039)

Latar Belakang

Ketika pandemi Covid-19 melanda Indonesia, harga emas mengalami fluktuatif karena kondisi perekonomian global yang sedang tidak baik. Indonesia mengalokasikan dana untuk penanganan pandemi dan para investor berusaha mencari perlindungan dalam berinvestasi, sehingga permintaan terhadap emas cenderung meningkat signifikan. Selain itu, stimulus fiskal dan besarnya moneter bank sentral global dan pemerintah menarik perhatian publik dalam berinvestasi emas yang menyebabkan keragaman harga emas mengalami kenaikan yang cukup fluktuatif. Oleh karena itu, metode peramalan permintaan harga emas yang tepat dan akurat diperlukan untuk masa depan (Dewi et al. 2022). Topik ini sangat menarik untuk dianalisis karena emas dianggap memiliki tren harga yang selalu meningkat, walaupun kondisi perekonomian tidak stabil.

Tinjauan Pustaka

Metode Pemulusan

SMA dan DMA

Rata-rata bergerak tunggal atau single moving average (SMA) merupakan suatu metode peramalan yang dilakukan dengan mengambil sekelompok nilai pengamatan, mencari nilai rata-rata tersebut sebagai ramalan untuk periode yang akan datang (Gempur 2010). Double moving average merupakan suatu metode peramalan yang dilakukan dengan menghitung nilai rata-rata bergerak sebanyak dua kali, kemudian dilanjutkan dengan meramal menggunakan suatu persamaan tertentu. Metode ini digunakan untuk mengatasi galat sistematis yang terjadi apabila rata-rata bergerak data dipakai pada data yang memiliki tren. Metode ini merupakan pengembangan dari metode simple moving average (SMA) (Muzaki 2022).

SES dan DES

Metode single exponential smoothing merupakan metode yang digunakan pada peramalan jangka pendek, biasanya hanya satu bulan ke depan, yang mengasumsikan bahwa data berfluktuasi di sekitar nilai mean tetap tanpa pola trend atau pertumbuhan konsisten /stasioner (Ginantra dan Anandita 2019). Metode double exponential smoothing merupakan metode yang digunakan untuk meramalkan data yang mengalami trend naik. Jika data yang digunakan semakin banyak dalam perhitungan peramalan, maka persentase error peramalannya akan semakin kecil dan sebaliknya (Hanief dan Purwanto 2017).

Ukuran Kesalahan Forecasting

1. SSE

Sum Square Error (SSE) menjumlahkan kuadrat dari kesalahan.

2. MSE dan RMSE

Mean Squared Error (MSE) merata-ratakan kuadrat dari kesalahan menggunakan penyebut n tanpa memperhatikan derajat bebas model. MSE banyak dipakai dalam optimalisasi pembobotan karena memberikan ketelitian yang lebih baik. Namun, Root Mean Square Error (RMSE) diperoleh dengan mengakarkan hasil MSE.

3. MAD

Simpangan absolut rata-rata atau mean absolute deviation (MAD) merata-ratakan nilai absolut kesalahan peramalan dalam unit ukuran yang sama seperti data aktual.

4. MAPE

Mean absolute percentage error (MAPE) merata-ratakan persentase absolut kesalahan dan memberikan gambaran seberapa kesalahan peramalan dibanding dengan nilai sebenarnya. MAPE banyak dipakai dalam perbandingan data-data yang mempunyai skala interval waktu yang berbeda.

ARIMA

Kestasioneran Data Deret Waktu

Data deret waktu dikatakan stasioner ketika rata-rata, ragam, dan peragam pada setiap lag sama dan konstan pada setiap waktu ke-t. Jika data deret waktu tidak memenuhi kriteria tersebut, data dapat dikatakan tidak stasioner (Aktivani 2020).

Uji Stasioner

Stasioneritas merupakan suatu proses pembangkitan yang mendasari suatu deret berkala yang didasarkan pada nilai tengah dan dan nilai ragam konstan. Untuk mengatasi ketidakstasioneran data deret waktu dalam nilai tengah atau rata-rata, metode yang digunakan adalah pembedaan (differencing), sedangkan untuk data yang tidak stasioner dalam ragam atau varians maka data perlu ditransformasi (Mulyana 2004).
Uji stasioneritas dapat dilakukan dengan dua metode, sebagai berikut:

1. Corellogram Test

Fungsi autokorelasi atau autocorrelation function (ACF) dibentuk dengan himpunan autokorelasi antara lag-k atau korelasi antara Zt dan Zt+k. Namun, fungsi autokorelasi parsial atau partial autocorrelation function (PACF) yang merupakan korelasi antara Yt dan Yt+k, setelah menghilangkan efek antara Y yang terletak di antara Yt dan Yt+k, sehingga Yt disumsikan sebagai konstanta. plot ACF dan PACF inilah yang dinamakan corellogram yang digunakan untuk melihat kestasioneran data. Jika nilai koefisien ACF sama pada setiap lag atau mendekati nol, maka data stasioner. Sebaliknya, data tidak stasioner ketika nilai koefisien ACF relatif tinggi. Jika plot ACF terlihat seperti histogram yang menurun (pola eksponensial), maka autokorelasi signifikan, dan jika diikuti oleh plot PACF seperti histogram yang langsung terpotong pada lag-2, maka data tidak stasioner (Aktivani 2020).

2. Augmented Dickey Fuller Test

Uji akar-akar unit atau augmented dickey fuller test (ADF test) mengasumsikan bahwa residuals et yang bersifat independen dengan rata-rata nol, ragam konstan, dan saling bebas (non autokorelasi). Namun, dalam banyak kasus, residuals et sering kali terdapat gejala autokorelasi. Oleh karena itu, ADF test perlu dikembangkan terhadap data yang mengandung autokorelasi pada residual et. ADF Test telah mempertimbangkan kemungkinan adanya autokorelasi pada error term jika series yang digunakan non stasioner. Jika nilai mutlak tau statistik dari uji ADF > nilai kritis ADF tabel, maka tolak H0 dan series dikatakan stasioner (Aktivani 2020).

Apabila data tidak stasioner, harus dilakukan differencing sampai data tersebut stasioner. Differencing merupakan pendekatan untuk menghitung perubahan atau selisih nilai observasi yang membuat data padat stasioner. Jika data asli atau aktual suatu series saling berintegrasi atau data sudah stasioner, maka data tersebut berintegrasi pada lag dengan order 0 (Aktivani 2020).

Model ARIMA

Menurut Salwa et al. (2018), ARIMA merupakan model yang secara penuh mengabaikan peubah bebas dalam peramalanforecasting* untuk memodelkan data yang tidak stasioner, sehingga cocok untuk menganalisis data harga penutupan emas periode Juli 2012-2022. Pada model ARIMA terdapat model-model lain yang dapat dibentuk, antara lain:

  1. Model Autoregressive (AR)

    Model AR merupakan model yang diregresikan terhadap nilai-nilai sebelumnya dari variabel itu sendiri. Model AR dengan ordo p disingkat menjadi AR(p) atau ARIMA(p,0,0).

    \[Y_t=\phi_1 Y_{t-1}+\phi_2 Y_{t-2}+...+\phi_p Y_{t-p}+e_t\]

  2. Model Moving Average (MA)

    Model MA merupakan model dengan orde q yang ditulis sebagai MA(q) atau ARIMA(0,0,q).

    \[Y_t=e_t-\theta_1 Y_{t-1}-\theta_2 Y_{t-2}-...-\theta_q Y_{t-q}\]

  3. Model Autoregressive Moving Average (ARMA)

    Model ARMA merupakan model gabungan dari Autoregresive (AR) dan Moving Average (MA). Model ARMA dapat dituliskan sebagai ARMA(p,q) atau ARIMA(p,0,q).

    \[Y_t=\phi_1 Y_{t-1}+\phi_2 Y_{t-2}+...+\phi_p Y_{t-p}+e_t-\theta_1 Y_{t-1}-\theta_2 Y_{t-2}-...-\theta_q Y_{t-q}\]

  4. Model Autoregressive Integrated Moving Average (ARIMA)

    Data yang mengalami ketidakstasioneran perlu dilakukan proses pembedaan atau differencing agar data menjadi stasioner. Adanya pembedaan ini tidak dapat dijelaskan oleh ketiga model sebelumnya sehingga digunakan model campuran yang disebut ARIMA(p,d,q).

    \[\phi_p (B)(1-B)^d Y_t=\theta_q (B)e_t\]

    dengan: \(B^k {Y_t}= Y_{t-k}\) dan \(B^k {e_t}= e_{t-k}\)
    Adapun notasi model ARIMA dirumuskan, sebagai berikut:

\[ ARIMA(p,d,q) \] Dimana,

p = parameter model autoregressive (AR)
d = parameter differencing
q = parameter model moving average (MA)

Menurut Ariwibowo (2020), metode baku yang digunakan untuk pemilihan model ARIMA menggunakan korelogram, yaitu autocorrelation function (ACF) dan partial autocorrelation function (PACF). Korelogram tersebut efektif untuk mengidentifikasi data dengan model AR dan MA murni, tetapi sulit mengidentifikasi model ARMA. Oleh karena itu, digunakan alat grafis lain untuk mengidentfikasi model ARMA berupa extended ACF (EACF). Pengidentifikasian digunakan dengan melihat pola berbentuk segitiga nol yang dengan ujung “0” segitiga lancip sebagai orde dari model ARMA.

Overfitting

Overfitting digunakan untuk memperluas model-model yang selanjutnya dipilih model yang paling sederhana yang menghasilkan sebuah model dengan ragam terbaik yang dapat mengakomodasi model. (Habibi dan Wardhani 2018).

Uji Diagnostik

Pemeriksaan diagnostik model merupakan langkah lanjutan setelah menentukan model sementara dari suatu analisis dengan melakukan uji signifikansi parameter model. Beberapa uji diagnostik yang digunakan, sebagai berikut:

1. Sisaan Menyebar Normal

Normalitas merupakan kondisi dimana sisaan berdistribusi normal. Dampak dari ketidaknormalan tidak terlalu serius, namun apabila ketidaknormalan tersebut disertai dengan ragam yang heterogen, masalahnya bisa menjadi serius (Hidayat 2013). Secara eksploratif, kenormalan galat dapat dilihat dari histogram, boxplot, dan plot peluang normal (plot kuantil-kuantil atau Q-Q Plot). Pada Q-Q Plot, pola pencaran titik-titik dalam plot yang membentuk garis lurus mengindikasikan sebaran data didekati pola sebaran normal. Uji normalitas digunakan untuk memeriksa kondisi kenormalan sebaran sisaan data yang dilakukan dengan Jarque Bera test dan Shapiro Wilk test. Pengujian ini dilakukan akibat Kalman filter bergantung pada asumsi sisaan data berdistribusi normal (Bernard et al. 2010).

2. Sisaan Saling Bebas

Autokorelasi merupakan keadaan dimana pada model ada korelasi antara sisaan pada periode t dengan sisaan pada periode sebelumnya (t-1). Uji ini digunakan untuk mendeteksi independensi sisaan yang terdapat pada lag tertentu. Kedua lag dikatakan tidak berkorelasi, jika antara lag tidak ada korelasi yang cukup berarti. Secara eksploratif, sisaan dapat dieksplorasi melalui plot yang mana tebaran amatan acak atau tidak berpola maka sisaan diindikasikan saling bebas (Mardiatmoko 2020).

3. Nilai Harapan Sisaan Sama Dengan Nol

Jika nilai harapan sisaan sama dengan nol, akan menghasilkan penduga parameter yang tidak bias. Secara eksploratif, dapat dieksplorasi melalui plot sisaan, amatan berada di sekitar nol maka rataan sisaan sama dengan nol. Selain itu, uji formal yang dapat dilakukan dengan t-test.

Analisis Data

Packages dan Library

Packages yang digunakan dalam harga penutupan emas (Juli 2012-2022) dengan metode pemulusan dan pemodelan ARIMA, sebagai berikut:

library(readxl)
library(forecast)
library(TTR)
library(imputeTS)
library(tseries)
library(ggplot2)
library(dplyr)
library(graphics)
library(TSA)
library(tidyverse)
library(lubridate)
library(gridExtra)
library(ggfortify)
library(cowplot)
library(graphics)
library(lmtest)
library(stats)
library(MASS)
library(fpp2)
library(FinTS)
library(car)
library(nortest)
library(fGarch)
library(rugarch)
library(aTSA)
library(dynlm)
library(psych)
library(moments)

Deskripsi Data

Data analisis yang digunakan adalahdata sekunder bulanan yang bersumber dari wesite World Gold Council mengenai harga penutupan emas periode Juli 2012-2022. Data ini terdiri atas satu peubah numerik. Data dapat diakses melalui tautan sebagai berikut: data emas. Dataset yang digunakan pada analisis dapat diakses melalui tautan sebagai berikut: dataset. Data ini dipilih untuk dianalisis lebih lanjut mengenai forecasting dengan model ARIMA.

Deklarasi Data

mpdw <- read_excel("D:/Kuliah Semester 5/MPDW/Gold Price Lag.xlsx", sheet=2)
ts.mpdw <- ts(mpdw)
kableExtra::kable(head(mpdw) ,caption = 'Subset Data Gold Price 2012-2022')
Subset Data Gold Price 2012-2022
Date Rupiah(gram) Inflasi
2012-07-31 483458.0 4.56
2012-08-31 496198.6 4.58
2012-09-28 535828.7 4.31
2012-10-31 538825.7 4.61
2012-11-30 532202.6 4.32
2012-12-31 522901.9 4.30

Statistik Deskriptif

summary(mpdw)
##       Date                         Rupiah(gram)       Inflasi     
##  Min.   :2012-07-31 00:00:00.00   Min.   :416772   Min.   :1.320  
##  1st Qu.:2015-01-30 00:00:00.00   1st Qu.:498327   1st Qu.:2.830  
##  Median :2017-07-31 00:00:00.00   Median :556474   Median :3.490  
##  Mean   :2017-07-30 13:53:03.47   Mean   :610643   Mean   :4.075  
##  3rd Qu.:2020-01-31 00:00:00.00   3rd Qu.:688126   3rd Qu.:4.890  
##  Max.   :2022-07-29 00:00:00.00   Max.   :930061   Max.   :8.790
describe(mpdw)
## Warning in FUN(newX[, i], ...): no non-missing arguments to min; returning Inf
## Warning in FUN(newX[, i], ...): no non-missing arguments to max; returning -Inf
##              vars   n      mean        sd    median  trimmed      mad       min
## Date            1 121       NaN        NA        NA      NaN       NA       Inf
## Rupiah(gram)    2 121 610643.43 143621.78 556474.10 596115.7 97401.48 416772.23
## Inflasi         3 121      4.07      1.99      3.49      3.9     1.50      1.32
##                    max     range skew kurtosis       se
## Date              -Inf      -Inf   NA       NA       NA
## Rupiah(gram) 930061.39 513289.16 0.85    -0.78 13056.53
## Inflasi           8.79      7.47 0.75    -0.35     0.18
z <- mpdw$`Rupiah(gram)`
boxplot(z, xlab="Harga Penutupan Emas")

Time Series Plot

DATES = seq(as.POSIXlt("2012-07-30 00:00:00", tz="UTC"), 
            as.POSIXlt("2022-07-29 00:00:00", tz="UTC"), length.out=40000)
dataset = data.frame(
  Outcome = c(rpois(20,1000),rpois(20,2000)),
  Week_End_Date = rep(DATES,2),
  Group = rep(c(FALSE,TRUE),each=20)
)
dataset$Condition = ifelse(
  dataset$Week_End_Date < as.POSIXlt("2020-03-06 00:00:00", tz="UTC"),
  "Before COVID-19",
  "After COVID-19"
)
FRAME = dataset %>% group_by(Condition) %>% 
  summarize(xmin=min(Week_End_Date),xmax=max(Week_End_Date))

p <- ggplot(mpdw, aes(x=Date, y=`Rupiah(gram)`)) +
  geom_line(lwd=1.2,col="#E26A2C")
p +labs( x="Year",y = "Gold Price (Rupiah/gram)",
         title="Time Series Plot Indonesia Gold Price ",
        subtitle = "Periode Juli 2012 - Juli 2022")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )+
  geom_rect(data=FRAME,inherit.aes=FALSE,
            aes(xmin=xmin,xmax=xmax,ymax=+Inf,ymin=-Inf,fill=Condition),
            alpha=0.1)+
  scale_fill_manual(values=c("red","blue"))

Berdasarkan time series plot, kenaikan harga penutupan emas terjadi secara signifikan selama pandemi Covid-19. Hal ini dikarenakan oleh harga penutupan emas yang akan meningkat pada saat terjadi resesi ekonomi, dilema geopolitik, atau ketidakpastian makroekonomi (Sui et al. 2021). Plot ini mengalami penurunan ketika kondisi perekonomian mulai membaik. Pola data harga penutupan emas yang terbentuk dari data tidak stasioner yang cenderung siklis dan tren positif.

Splitting Data

Data aktual dibagi menjadi data training dan data testing sebelum dilakukan pemodelan ARIMA. Dalam hal ini yang menjadi data training adalah data ke-1 pada Juli 2012 sampai data ke-100 pada November 2020, sedangkan yang menjadi data testing adalah data ke-101 pada Desember 2020 sampai data ke-120 pada Juli 2022. Kemudian, data diubah menjadi format data time series.

Harga.emas <- mpdw$`Rupiah(gram)`
training<-mpdw[1:100,2]
testing<-mpdw[101:120,2]

training.ts<-ts(training)
testing.ts<-ts(testing,start=101)

Time Series Plot Data Training dan Testing

DATES = seq(as.POSIXlt("2012-07-30 00:00:00", tz="UTC"), 
            as.POSIXlt("2022-07-29 00:00:00", tz="UTC"), length.out=40000)
dataset = data.frame(
  Outcome = c(rpois(20,1000),rpois(20,2000)),
  Week_End_Date = rep(DATES,2),
  Group = rep(c(FALSE,TRUE),each=20)
)
dataset$Spliting = ifelse(
  dataset$Week_End_Date < as.POSIXlt("2020-12-31 00:00:00", tz="UTC"),
  "Training",
  "Testing"
)
FRAME = dataset %>% group_by(Spliting) %>% 
  summarize(xmin=min(Week_End_Date),xmax=max(Week_End_Date))

p <- ggplot(mpdw, aes(x=Date, y=`Rupiah(gram)`)) +
  geom_line(lwd=1.2,col="#E26A2C")
p +labs( x="Year",y = "Gold Price (Rupiah/gram)",
         title="Spliting Plot Training & Testing ",
        subtitle = "Gold Price")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )+
  geom_rect(data=FRAME,inherit.aes=FALSE,
            aes(xmin=xmin,xmax=xmax,ymax=+Inf,ymin=-Inf,fill=Spliting),
            alpha=0.1)+
  scale_fill_manual(values=c("green","orange"))

Berdasarkan time series plot data training dan testing mengenai harga penutupan emas, terlihat bahwa data tidak konstan karena plot data tidak menyebar di sekitar rata-rata serta ragam yang konstan.Pada plot ini menyebar membentuk suatu kecenderungan menaik. Jadi, data training dan testing mengenai harga penutupan emas tidak stasioner dalam rata-rata dan ragam.

Ubah Format Data Date

golddates<-data.frame(x=seq(as.Date("2012-07-31"), as.Date("2022-08-29"), by="months"))
mpdw <- mpdw[,-1]
mpdw <- cbind(golddates,mpdw)
colnames(mpdw) <- c("Date","GoldPrice")
str(mpdw)
## 'data.frame':    121 obs. of  3 variables:
##  $ Date     : Date, format: "2012-07-31" "2012-08-31" ...
##  $ GoldPrice: num  483458 496199 535829 538826 532203 ...
##  $ NA       : num  4.56 4.58 4.31 4.61 4.32 4.3 4.57 5.31 5.9 5.57 ...
mpdw <- mpdw[,-3]

Metode Pemulusan

Parameter Optimum

Untuk mencari parameter optimum untuk metode single moving average (SMA) dan double moving average (DMA) dibuat sebuah fungsi optimasi dengan memanfaatkan iterasi dan evaluasi metrik tiap-tiap parameter.

#Fungsi Mencari Parameter Optimum
pemulusan <- function(x,min,max,method=c("SMA","DMA"),
                      metrik=c("SSE","MSE","RMSE","MAD","MAPE")){
  df.master <- data.frame()
  if(method=="SMA"){
    for(i in min:max) {
      sma <- TTR::SMA(x,i)
      ramal <- c(NA,sma)
      df <- cbind(aktual=c(x,NA),mulus=c(sma,NA),ramal)
      error <- df[,1]-df[,3]
      sse <- sum(error^2,na.rm=T)
      mse <- mean(error^2,na.rm=T)
      rmse <- sqrt(mse)
      mad <- mean(abs(error),na.rm=T)
      rerror <- error/df[,1]*100
      mape <- mean(abs(rerror),na.rm=T)
      ak <- data.frame(n=i,SSE=sse,MSE=mse,RMSE=rmse,MAD=mad,MAPE=mape)
      df.master <- rbind(df.master,ak)
    }
  }
  else if(method=="DMA") {
    for(i in min:max) {
      df_ts_sma <- TTR::SMA(x,  n=i)
      df_ts_dma <- TTR::SMA(df_ts_sma, n=i)
      At <- 2*df_ts_sma-df_ts_dma
      Bt <- df_ts_sma-df_ts_dma
      pemulusan_dma <- At+Bt
      ramal_dma <- c(NA, pemulusan_dma)
      df_dma <- cbind(df_aktual=c(x,NA), pemulusan_dma=c(pemulusan_dma,NA), ramal_dma)
      error.dma <- df_dma[, 1] - df_dma[, 3]
      SSE.dma <- sum(error.dma^2, na.rm = T)
      MSE.dma <- mean(error.dma^2, na.rm = T)
      RMSE.dma <- sqrt(mean(error.dma^2, na.rm = T))
      MAD.dma <- mean(abs(error.dma), na.rm = T)
      r.error.dma <- (error.dma/df_dma[, 1])*100 
      MAPE.dma <- mean(abs(r.error.dma), na.rm = T)
      ak <- data.frame(n=i,SSE=SSE.dma,MSE=MSE.dma,RMSE=RMSE.dma,
                       MAD=MAD.dma,MAPE=MAPE.dma)
      df.master <- rbind(df.master,ak)
    }
  }
  opt <- df.master$n[which.min(df.master[,metrik])]
  return(list(Akurasi=df.master,n.optimum=paste("n optimum adalah",opt)))
}

  1. Single Moving Average (SMA)

pemulusan(training.ts,min=2,max=10,method = "SMA",metrik = "RMSE")
## $Akurasi
##    n          SSE        MSE     RMSE      MAD     MAPE
## 1  2  71292639721  727475916 26971.76 18685.96 3.196290
## 2  3  86064038922  887258133 29786.88 20895.64 3.576735
## 3  4  99709281991 1038638354 32227.91 22523.35 3.841194
## 4  5 116053193716 1221612565 34951.57 24230.63 4.097525
## 5  6 135354016771 1439936349 37946.49 26215.03 4.394960
## 6  7 156253446748 1680144589 40989.57 28567.13 4.754873
## 7  8 177727646973 1931822250 43952.50 30357.02 5.002425
## 8  9 198805732325 2184678377 46740.54 31605.91 5.148045
## 9 10 216638269636 2407091885 49062.12 32485.32 5.216134
## 
## $n.optimum
## [1] "n optimum adalah 2"

Parameter optimum pada data training yang didapat untuk single moving average (SMA) berdasarkan nilai RMSE terkecil adalah n = 2.

  1. Double Moving Average (DMA)

pemulusan(training.ts,min=2,max=10,method = "DMA",metrik = "RMSE")
## $Akurasi
##    n          SSE        MSE     RMSE      MAD     MAPE
## 1  2  67312008515  693938232 26342.71 18830.46 3.156857
## 2  3  94458997939  994305241 31532.61 23160.58 3.980181
## 3  4  96336335105 1035874571 32185.01 23572.41 4.109238
## 4  5 104026903502 1143152786 33810.54 24364.05 4.252934
## 5  6 112293594797 1261725784 35520.78 25858.43 4.500373
## 6  7 113492628686 1304512973 36118.04 27344.49 4.824881
## 7  8  99133249999 1166273529 34150.75 26540.86 4.688529
## 8  9  91160398954 1098318060 33140.88 26194.97 4.607729
## 9 10  84378213381 1041706338 32275.48 25242.17 4.402280
## 
## $n.optimum
## [1] "n optimum adalah 2"

Parameter optimum pada data training yang didapat untuk double moving average (DMA) berdasarkan nilai RMSE terkecil adalah n = 2.

Berdasarkan perbandingan kedua metode, parameter optimum didapatkan nilai RMSE dari DMA lebih kecil daripada nilai RMSE SMA, sehingga metode pemulusan moving average akan dilanjutkan menggunakan metode DMA.

Untuk mencari parameter optimum single exponential smoothing (SES) dan double exponential smoothing (DES) dibuat sebuah fungsi optimasi dengan memanfaatkan iterasi dan evaluasi metrik tiap-tiap parameter (alpha dan beta).

#Fungsi Mencari Parameter Optimum
pemulusan2 <- function(x,alfa=NULL,beta=NULL,method=c("SES","DES"),
                       metrik=c("SSE","MSE","RMSE")){
  df.master <- data.frame()
  if(method=="SES"){
    for(i in alfa){
      df_ses <- HoltWinters(x, alpha = i, beta=F, gamma=F)
      sse <- df_ses$SSE 
      mse <- sse/length(x)
      rmse <-sqrt(mse)
      ak <- data.frame(Alpha=i,SSE=sse,MSE=mse,RMSE=rmse)
      df.master <- rbind(df.master,ak)
      datases <- data.frame(x, c(NA, df_ses$fitted[,1]))
      colnames(datases) = c("y","yhat")
    }
  }
  else if(method=="DES"){
    for (i in alfa) {
      for (j in beta) {
        df_des <- HoltWinters(x, alpha = i, beta=j, gamma=F)
        sse <- df_des$SSE 
        mse <- sse/length(x)
        rmse <-sqrt(mse)
        ak <- data.frame(Alpha=i,Beta=j,SSE=sse,MSE=mse,RMSE=rmse)
        df.master <- rbind(df.master,ak)
        datades <- data.frame(x, c(NA, NA, df_des$fitted[,1]))
        colnames(datades) = c("y","yhat")
      }
    }
  }
  if(method=="SES"){
    opt <- df.master$Alpha[which.min(df.master[,metrik])]
    return(list(Akurasi=df.master,n.optimum=paste("Alpha optimum adalah",opt)))
  }
  else if(method=="DES"){
    opt <- df.master$Alpha[which.min(df.master[,metrik])]
    opt2 <- df.master$Beta[which.min(df.master[,metrik])]
    return(list(Akurasi=df.master,n.optimum=paste("Alpha optimum adalah",opt,
                                                  ",Beta optimum adalah",opt2)))
  }
}

  1. Single Exponential Smoothing (SES)

pemulusan2(training.ts,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           method="SES",metrik="RMSE")
## $Akurasi
##   Alpha          SSE        MSE     RMSE
## 1   0.1 319932001361 3199320014 56562.53
## 2   0.2 167606453102 1676064531 40939.77
## 3   0.3 113747726334 1137477263 33726.51
## 4   0.4  87910102218  879101022 29649.64
## 5   0.5  73319442782  733194428 27077.56
## 6   0.6  64135818613  641358186 25325.05
## 7   0.7  57854216619  578542166 24052.90
## 8   0.8  53250703266  532507033 23076.11
## 9   0.9  49681612303  496816123 22289.37
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9"

Parameter optimum pada data training yang didapat dari metode SES berdasarkan nilai RMSE terkecil adalah α = 0.9.

  1. Double Exponential Smoothing (DES)

pemulusan2(training.ts,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           beta=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),method="DES",metrik="RMSE")
## $Akurasi
##    Alpha Beta          SSE        MSE     RMSE
## 1    0.1  0.1 275115382664 2751153827 52451.44
## 2    0.1  0.2 210063165784 2100631658 45832.65
## 3    0.1  0.3 181178385285 1811783853 42565.05
## 4    0.1  0.4 171370258093 1713702581 41396.89
## 5    0.1  0.5 171537452893 1715374529 41417.08
## 6    0.1  0.6 180801620939 1808016209 42520.77
## 7    0.1  0.7 189617892086 1896178921 43545.14
## 8    0.1  0.8 192271597859 1922715979 43848.79
## 9    0.1  0.9 190947668144 1909476681 43697.56
## 10   0.2  0.1 132667303194 1326673032 36423.52
## 11   0.2  0.2 113968292951 1139682930 33759.19
## 12   0.2  0.3 111773476093 1117734761 33432.54
## 13   0.2  0.4 114215352767 1142153528 33795.76
## 14   0.2  0.5 118136826055 1181368261 34371.04
## 15   0.2  0.6 123287008634 1232870086 35112.25
## 16   0.2  0.7 128161083515 1281610835 35799.59
## 17   0.2  0.8 130890190867 1308901909 36178.75
## 18   0.2  0.9 131120770276 1311207703 36210.60
## 19   0.3  0.1  93072419012  930724190 30507.77
## 20   0.3  0.2  86896897413  868968974 29478.28
## 21   0.3  0.3  88128641095  881286411 29686.47
## 22   0.3  0.4  90893491975  908934920 30148.55
## 23   0.3  0.5  93433332216  934333322 30566.87
## 24   0.3  0.6  95059976061  950599761 30831.80
## 25   0.3  0.7  96018056842  960180568 30986.78
## 26   0.3  0.8  96897695252  968976953 31128.39
## 27   0.3  0.9  98105773440  981057734 31321.84
## 28   0.4  0.1  75464704681  754647047 27470.84
## 29   0.4  0.2  73528823790  735288238 27116.20
## 30   0.4  0.3  75354642928  753546429 27450.80
## 31   0.4  0.4  77565731915  775657319 27850.63
## 32   0.4  0.5  79526062183  795260622 28200.37
## 33   0.4  0.6  81385847679  813858477 28528.20
## 34   0.4  0.7  83426388280  834263883 28883.63
## 35   0.4  0.8  85802562007  858025620 29292.07
## 36   0.4  0.9  88552014008  885520140 29757.69
## 37   0.5  0.1  65682876256  656828763 25628.67
## 38   0.5  0.2  65617692423  656176924 25615.95
## 39   0.5  0.3  67721549740  677215497 26023.36
## 40   0.5  0.4  70048834501  700488345 26466.74
## 41   0.5  0.5  72425943347  724259433 26912.07
## 42   0.5  0.6  74994376116  749943761 27385.10
## 43   0.5  0.7  77860402881  778604029 27903.48
## 44   0.5  0.8  81048945294  810489453 28469.10
## 45   0.5  0.9  84531066410  845310664 29074.23
## 46   0.6  0.1  59508881494  595088815 24394.44
## 47   0.6  0.2  60474184459  604741845 24591.50
## 48   0.6  0.3  62857320059  628573201 25071.36
## 49   0.6  0.4  65465015070  654650151 25586.13
## 50   0.6  0.5  68237287903  682372879 26122.27
## 51   0.6  0.6  71238632432  712386324 26690.57
## 52   0.6  0.7  74479456341  744794563 27290.92
## 53   0.6  0.8  77909412309  779094123 27912.26
## 54   0.6  0.9  81437471121  814374711 28537.25
## 55   0.7  0.1  55209791359  552097914 23496.76
## 56   0.7  0.2  56799180875  567991809 23832.58
## 57   0.7  0.3  59386009034  593860090 24369.24
## 58   0.7  0.4  62173188155  621731882 24934.55
## 59   0.7  0.5  65110687786  651106878 25516.80
## 60   0.7  0.6  68192077351  681920774 26113.61
## 61   0.7  0.7  71359135809  713591358 26713.13
## 62   0.7  0.8  74513164829  745131648 27297.10
## 63   0.7  0.9  77541961859  775419619 27846.36
## 64   0.8  0.1  51942129688  519421297 22790.82
## 65   0.8  0.2  53894705082  538947051 23215.23
## 66   0.8  0.3  56560280196  565602802 23782.41
## 67   0.8  0.4  59358690218  593586902 24363.64
## 68   0.8  0.5  62222172007  622221720 24944.37
## 69   0.8  0.6  65096799519  650967995 25514.07
## 70   0.8  0.7  67897375077  678973751 26057.12
## 71   0.8  0.8  70531519716  705315197 26557.77
## 72   0.8  0.9  72927123269  729271233 27005.02
## 73   0.9  0.1  49273228094  492732281 22197.57
## 74   0.9  0.2  51403233221  514032332 22672.28
## 75   0.9  0.3  54030140555  540301406 23244.38
## 76   0.9  0.4  56703800664  567038007 23812.56
## 77   0.9  0.5  59346238746  593462387 24361.08
## 78   0.9  0.6  61893515218  618935152 24878.41
## 79   0.9  0.7  64279268937  642792689 25353.36
## 80   0.9  0.8  66458160436  664581604 25779.48
## 81   0.9  0.9  68421069639  684210696 26157.42
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9 ,Beta optimum adalah 0.1"

Parameter optimum pada data training yang didapat dari metode DES berdasarkan nilia RMSE terkecil adalah α = 0.9 dan β = 0.1.

Berdasarkan perbandingan metrik nilai RMSE kedua metode dengan parameter optimum, metode DES dinyatakan lebih baik dibandingkan metode SES.

Forecasting

Double Moving Average

#Single Moving  Average dengan n=2
df_ts_sma <- TTR::SMA(training.ts,  n=2)
#Double Moving  Average dengan n=2
df_ts_dma <- TTR::SMA(df_ts_sma, n=2)
At <- 2*df_ts_sma-df_ts_dma
Bt <- df_ts_sma-df_ts_dma
pemulusan_dma <- At+Bt
ramal_dma <- c(NA, pemulusan_dma)
df_dma <- cbind(df_aktual=c(training.ts,NA), pemulusan_dma=c(pemulusan_dma,NA), ramal_dma)
head(df_dma)
##      df_aktual pemulusan_dma ramal_dma
## [1,]  483458.0            NA        NA
## [2,]  496198.6            NA        NA
## [3,]  535828.7      542199.0        NA
## [4,]  538825.7      558640.8  542199.0
## [5,]  532202.6      533701.1  558640.8
## [6,]  522901.9      519590.4  533701.1

Eksploratif: DMA Plot

ts.plot(mpdw[,2], xlab="periode waktu", ylab="Yt", col="blue", lty=1,
        xlim=c(60,120),ylim=c(600000,1000000))
points(mpdw[,2],col="black")
lines(pemulusan_dma,col="red",lwd=2,lty=1)
lines(ramal_dma,col="green",lwd= 2,lty=1)
title("Rataan Bergerak  Berganda N=2", cex.main=1, font.main=4 ,col.main="black")
legend("topright", c("Data Aktual","Pemulusan","Ramalan"),lty=1,col=c ("blue","red","green"))

Berdasarkan times series plot metode DMA di atas, dapat diketahui n yang optimum adalah 2. Hal ini dikarenakan oleh data training DMA dengan n = 2 lebih mendekati sebaran pola dari data aktual dibandingkan hasil peralamannya.

Double Exponential Smoothing

#DES alpha=0.9 dan beta=0.1
df_des <- HoltWinters(training.ts, alpha = 0.9, beta=0.1, gamma=F)
datades <- data.frame(training.ts, c(NA, NA, df_des$fitted[,1]))
colnames(datades) = c("y","yhat")
head(datades)
##          y     yhat
## 1 483458.0       NA
## 2 496198.6       NA
## 3 535828.7 508939.1
## 4 538825.7 548300.4
## 5 532202.6 554081.1
## 6 522901.9 546729.3
ramal_des <- forecast::forecast(df_des,h=3)
(df_ramal_des <- data.frame(ramal_des))
##     Point.Forecast    Lo.80    Hi.80    Lo.95     Hi.95
## 101       916747.0 887865.0 945629.0 872575.8  960918.2
## 102       931748.9 891107.3 972390.5 869592.9  993904.9
## 103       946750.8 895518.8 997982.8 868398.2 1025103.3

Eksploratif: DES Plot

#Plot 
ts.plot(mpdw[,2], xlab="periode  waktu", ylab="Yt", col="blue", lty=1,
        xlim=c(70,120),ylim=c(250000,1000000))
points(mpdw[,2],col="white")
lines(datades[,2], col="red",lwd=2) #nilai dugaan
lines(ts(df_ramal_des$Point.Forecast,start = 101), col="green",lwd=2) #nilai dugaan
title("Double Exponential Smoothing Alpha=0.9 Beta = 0.1",cex.main=1, font.main=4 ,col.main="black")
legend("bottomright", c("Data aktual", "Fitted DES","Ramalan"), lty=1,col=c ("blue","red","green"))

Berdasarkan eksplorasi dengan time series plot di atas, plot data training hasil fitted menggunakan metode pemulusan DES dengan α = 0.9 dan β = 0.1 menampilkan plot yang berpola tidak terlalu jauh dari plot data aktual dibandingkan dengan data peramalannya.

Perbandingan Nilai Akurasi

pemulusan(testing.ts,min=2,max=10,method = "SMA",metrik = "RMSE")
## $Akurasi
##    n         SSE       MSE     RMSE      MAD     MAPE
## 1  2 10139865228 563325846 23734.49 19809.27 2.343293
## 2  3 12010271147 706486538 26579.81 22716.22 2.685046
## 3  4 11856156942 741009809 27221.50 22271.09 2.625281
## 4  5 10432520867 695501391 26372.36 20113.23 2.344098
## 5  6 10562837839 754488417 27467.95 20387.10 2.363185
## 6  7 10175325671 782717359 27977.09 20423.60 2.358827
## 7  8 10105103595 842091966 29018.82 21645.20 2.496351
## 8  9 10003329995 909393636 30156.15 23271.85 2.683494
## 9 10  9879990416 987999042 31432.45 25131.91 2.896108
## 
## $n.optimum
## [1] "n optimum adalah 2"
pemulusan(testing.ts,min=2,max=10,method = "DMA",metrik = "RMSE")
## $Akurasi
##    n         SSE        MSE      RMSE       MAD      MAPE
## 1  2 10452970383  614880611 24796.786 20698.171 2.4488844
## 2  3 15975351016 1065023401 32634.696 26793.426 3.1502609
## 3  4 13915890752 1070453135 32717.780 28652.887 3.3833594
## 4  5 12963135616 1178466874 34328.805 30623.288 3.6007957
## 5  6 11264328490 1251592054 35377.847 29508.189 3.4094721
## 6  7  9155224532 1307889219 36164.751 31012.679 3.5421124
## 7  8  7367720006 1473544001 38386.769 32685.156 3.6948425
## 8  9  1837115191  612371730 24746.146 17501.813 1.9696508
## 9 10    14320934   14320934  3784.301  3784.301 0.4360753
## 
## $n.optimum
## [1] "n optimum adalah 10"
pemulusan2(testing.ts,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           method="SES",metrik="RMSE")
## $Akurasi
##   Alpha         SSE       MSE     RMSE
## 1   0.1 14056054409 702802720 26510.43
## 2   0.2 12407908901 620395445 24907.74
## 3   0.3 11012806633 550640332 23465.73
## 4   0.4 10023869365 501193468 22387.35
## 5   0.5  9289086603 464454330 21551.20
## 6   0.6  8696917416 434845871 20852.96
## 7   0.7  8195312549 409765627 20242.67
## 8   0.8  7764395852 388219793 19703.29
## 9   0.9  7398137316 369906866 19232.96
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9"
pemulusan2(testing.ts,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           beta=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),method="DES",metrik="RMSE")
## $Akurasi
##    Alpha Beta         SSE        MSE     RMSE
## 1    0.1  0.1 56252227558 2812611378 53034.06
## 2    0.1  0.2 40270514965 2013525748 44872.33
## 3    0.1  0.3 29894089848 1494704492 38661.41
## 4    0.1  0.4 23229340722 1161467036 34080.30
## 5    0.1  0.5 18990026159  949501308 30813.98
## 6    0.1  0.6 16320736145  816036807 28566.36
## 7    0.1  0.7 14665745126  733287256 27079.28
## 8    0.1  0.8 13673071768  683653588 26146.77
## 9    0.1  0.9 13125546550  656277328 25617.91
## 10   0.2  0.1 25965374709 1298268735 36031.50
## 11   0.2  0.2 18044060239  902203012 30036.69
## 12   0.2  0.3 14558131095  727906555 26979.74
## 13   0.2  0.4 13177333286  658866664 25668.40
## 14   0.2  0.5 12855637167  642781858 25353.14
## 15   0.2  0.6 13133980019  656699001 25626.14
## 16   0.2  0.7 13802648320  690132416 26270.37
## 17   0.2  0.8 14752039262  737601963 27158.83
## 18   0.2  0.9 15913701319  795685066 28207.89
## 19   0.3  0.1 16936696515  846834826 29100.43
## 20   0.3  0.2 13354572252  667728613 25840.45
## 21   0.3  0.3 12483262111  624163106 24983.26
## 22   0.3  0.4 12731436014  636571801 25230.37
## 23   0.3  0.5 13535227840  676761392 26014.64
## 24   0.3  0.6 14661153863  733057693 27075.04
## 25   0.3  0.7 15980328439  799016422 28266.88
## 26   0.3  0.8 17396646013  869832301 29492.92
## 27   0.3  0.9 18821118203  941055910 30676.63
## 28   0.4  0.1 13400305871  670015294 25884.65
## 29   0.4  0.2 11866507388  593325369 24358.27
## 30   0.4  0.3 12018196705  600909835 24513.46
## 31   0.4  0.4 12806225691  640311285 25304.37
## 32   0.4  0.5 13864121680  693206084 26328.81
## 33   0.4  0.6 15001561545  750078077 27387.55
## 34   0.4  0.7 16078903194  803945160 28353.93
## 35   0.4  0.8 16984920473  849246024 29141.83
## 36   0.4  0.9 17644785790  882239289 29702.51
## 37   0.5  0.1 11596603216  579830161 24079.66
## 38   0.5  0.2 11045083010  552254150 23500.09
## 39   0.5  0.3 11557158928  577857946 24038.68
## 40   0.5  0.4 12390989796  619549490 24890.75
## 41   0.5  0.5 13266926171  663346309 25755.51
## 42   0.5  0.6 14031520659  701576033 26487.28
## 43   0.5  0.7 14600650320  730032516 27019.11
## 44   0.5  0.8 14954639651  747731983 27344.69
## 45   0.5  0.9 15129225248  756461262 27503.84
## 46   0.6  0.1 10447554946  522377747 22855.58
## 47   0.6  0.2 10358276076  517913804 22757.72
## 48   0.6  0.3 10933646514  546682326 23381.24
## 49   0.6  0.4 11628784870  581439243 24113.05
## 50   0.6  0.5 12247627902  612381395 24746.34
## 51   0.6  0.6 12714582385  635729119 25213.67
## 52   0.6  0.7 13029182087  651459104 25523.70
## 53   0.6  0.8 13237731610  661886580 25727.16
## 54   0.6  0.9 13402039009  670101950 25886.33
## 55   0.7  0.1  9602951628  480147581 21912.27
## 56   0.7  0.2  9722285744  486114287 22048.00
## 57   0.7  0.3 10256141967  512807098 22645.24
## 58   0.7  0.4 10808323689  540416184 23246.85
## 59   0.7  0.5 11263407715  563170386 23731.21
## 60   0.7  0.6 11607579227  580378961 24091.06
## 61   0.7  0.7 11874809516  593740476 24366.79
## 62   0.7  0.8 12110869422  605543471 24607.79
## 63   0.7  0.9 12350525227  617526261 24850.08
## 64   0.8  0.1  8938105751  446905288 21140.13
## 65   0.8  0.2  9148637723  457431886 21387.66
## 66   0.8  0.3  9625345882  481267294 21937.81
## 67   0.8  0.4 10084150928  504207546 22454.57
## 68   0.8  0.5 10465419725  523270986 22875.12
## 69   0.8  0.6 10781241265  539062063 23217.71
## 70   0.8  0.7 11064153260  553207663 23520.37
## 71   0.8  0.8 11341512861  567075643 23813.35
## 72   0.8  0.9 11627137036  581356852 24111.34
## 73   0.9  0.1  8403941855  420197093 20498.71
## 74   0.9  0.2  8655800158  432790008 20803.61
## 75   0.9  0.3  9089925860  454496293 21318.92
## 76   0.9  0.4  9501276553  475063828 21795.96
## 77   0.9  0.5  9860834409  493041720 22204.54
## 78   0.9  0.6 10185835328  509291766 22567.49
## 79   0.9  0.7 10499397091  524969855 22912.22
## 80   0.9  0.8 10816469906  540823495 23255.61
## 81   0.9  0.9 11142916712  557145836 23603.94
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9 ,Beta optimum adalah 0.1"

Perbandingan Nilai Akurasi Model

Berdasarkan data testing yang digunakan, model single exponential smoothing (SES) merupakan model terbaik yang dilihat dari nilai RMSE terkecil.

Pemodelan ARIMA

Tahapan Pemodelan

Metode Box-Jenskins

Cek Tails Off / Cuts Off

ggAcf(training.ts,col="black",lwd=2)+labs( x="Lag",y = "ACF",
         title="Plot ACF Gold Price",
        subtitle = "(Juli 2012-2022)")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5))

ggPacf(training.ts,col="black",lwd=2)+labs( x="Lag",y = "PACF",
         title="Plot PACF Gold Price  ",
        subtitle = "(Juli 2012-2022)",col="white")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )

Berdasarkan plot ACF pada data harga penutupan emas, nilai korelasi setiap lag sebagaimana terlihat pada plot di atas menurun menuju nilai nol secara perlahan (tails off slowly). Artinya, data harga penutupan emas tidak stasioner.

Augmented Dickey-Fuller Test

Hipotesis yang diuji:
\(H_0\): \(\rho = 0\) (data tidak stasioner)
\(H_1\): \(\rho \neq 0\) (data stasioner)
dengan menggunakan α=5%.

adf.test(training.ts)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag  ADF p.value
## [1,]   0 2.18   0.990
## [2,]   1 1.46   0.962
## [3,]   2 2.35   0.990
## [4,]   3 2.08   0.990
## [5,]   4 1.89   0.984
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 1.472   0.990
## [2,]   1 0.565   0.987
## [3,]   2 2.781   0.990
## [4,]   3 2.521   0.990
## [5,]   4 2.428   0.990
## Type 3: with drift and trend 
##      lag    ADF p.value
## [1,]   0 -0.360   0.986
## [2,]   1 -1.295   0.868
## [3,]   2  0.432   0.990
## [4,]   3  0.352   0.990
## [5,]   4  0.239   0.990
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

Berdasarkan hasil augmented dickey-Fuller test (ADF test) didapatkan p-value > α, maka tak tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa data stasioner pada taraf nyata 5%, sehingga perlu dilakukan differencing.

Diferrencing Pertama

train.diff<-diff(training.ts,differences = 1) 

traindif <- as.data.frame(train.diff)
xdif <- as.Date(mpdw[2:100,1])
xdif <- as.data.frame(xdif)
data1 <- cbind(xdif,traindif)

ggplot(data1, aes(x=xdif, y=`Rupiah(gram)`)) +
  geom_line(lwd=1.2,col="red3")+
  ggtitle("Plot Data Difference Train Gold Price")+
  xlab("Time Period")+ylab("Data Difference Train Gold Price")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )

Setelah dilakukan differencing satu kali d=1, pola data harga penutupan emas sudah stasioner dilihat dari times series plot di atas.

Cek Stasioneritas Setelah Differencing

Stasioner dalam Rataan

Hipotesis yang diuji:
\(H_0\): \(\rho = 0\) (data tidak stasioner dalam rataan)
\(H_1\): \(\rho \neq 0\) (data stasioner dalam rataan)
dengan menggunakan α=5%.

adf.test(train.diff)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -7.37    0.01
## [2,]   1 -8.44    0.01
## [3,]   2 -5.38    0.01
## [4,]   3 -4.01    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -7.54    0.01
## [2,]   1 -8.77    0.01
## [3,]   2 -5.71    0.01
## [4,]   3 -4.33    0.01
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -7.81    0.01
## [2,]   1 -9.74    0.01
## [3,]   2 -6.65    0.01
## [4,]   3 -5.40    0.01
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

Diperoleh p-value = 0.01 < α = 0.05, maka tolak H0. Artinya, cukup bukti untuk mengatakan bahwa data stasioner dalam rataan setelah dilakukan differencing sebanyak satu kali pada taraf nyata 5%.

ACF & PACF Plot

ggAcf(train.diff,col="black",lwd=2)+labs( x="Lag",y = "ACF",
         title="Plot ACF Gold Price  ",
        subtitle = "(Juli 2012-2022)")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )

ggPacf(train.diff,col="black",lwd=2)+labs( x="Lag",y = "PACF",
         title="Plot PACF Gold Price  ",
        subtitle = "(Juli 2012-2022)",col="white")+
  theme_bw()+
  theme(
    plot.title = element_text(size = 14L,
                              face = "bold",
                              hjust = 0.5),
    plot.subtitle = element_text(size = 11L,
                              face = "plain",
                              hjust = 0.5)
  )

Berdasarkan plot ACF dan PACF di atas,terlihat bahwa nilai korelasi antara data dengan lag seperti gambar di atas tidak turun secara perlahan, dimana diperoleh cuts off pada lag ke-2. Berdasarkan hasil eksplorasi di atas, kandidat model yang dapat dibentuk secara berurutan adalah ARIMA(0,1,2) dan ARIMA(2,1,0).

Extended ACF Matriks

eacf(train.diff)
## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x x o o o o o o o o o  o  o  o 
## 1 x x o o o o o o o o o  o  o  o 
## 2 x x o o o o o o o o o  o  o  o 
## 3 x o o o o o o o o o o  o  o  o 
## 4 x o o o o o o o o o o  o  o  o 
## 5 x o o o x o o o o o o  o  o  o 
## 6 o x x x o o o o o o o  o  o  o 
## 7 o x o x o o o o o o o  o  o  o

Kandidat Model EACF

Dengan menggunakan matriks EACF, dapat diperoleh kandidat model lainnya.

Identifikasi Model

Berdasarkan plot ACF, PACF, dan matrisk EACF, diperoleh tiga kandidat model beserta orde parameternya, sebagai berikut:

ARIMA(2,1,0), dengan orde p = 2, d = 1, dan q = 0

ARIMA(0,1,2), dengan orde p = 0, d = 1, dan q = 2

ARIMA(1,1,2), dengan orde p = 1, d = 1, dan q = 2

Signifikansi Model

model1 <- Arima(training.ts,order = c(2,1,0),method ="ML")
model2 <- Arima(training.ts,order = c(0,1,2),method ="ML")
model3 <- Arima(training.ts,order = c(1,1,2),method ="ML")
#ARIMA(2,1,0)
coeftest(model1)    #signifikan semua
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.367742   0.094573  3.8884 0.0001009 ***
## ar2 -0.334850   0.095580 -3.5033 0.0004595 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#ARIMA(0,1,2)
coeftest(model2)    #ma2 tidak signifikan
## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ma1  0.38379    0.10289  3.7301 0.0001914 ***
## ma2 -0.17113    0.10521 -1.6265 0.1038404    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#ARIMA(1,1,2)
coeftest(model3)    #tidak signifikan semua
## 
## z test of coefficients:
## 
##     Estimate Std. Error z value Pr(>|z|)  
## ar1  0.23134    0.34151  0.6774  0.49816  
## ma1  0.16412    0.32532  0.5045  0.61393  
## ma2 -0.26783    0.15247 -1.7567  0.07898 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Perbandingan Kebaikan Kandidat Model

Model <- c("ARIMA (2,1,0)","ARIMA (0,1,2)","ARIMA (1,1,2)")
AIC <- c(model1$aic,model2$aic,model3$aic)
BIC <- c(model1$bic,model2$bic,model3$bic)
Akurasi <- data.frame(Model,AIC,BIC)
kableExtra::kable(Akurasi)
Model AIC BIC
ARIMA (2,1,0) 2244.889 2252.674
ARIMA (0,1,2) 2245.924 2253.709
ARIMA (1,1,2) 2247.509 2257.889
paste("Model yang terbaik adalah model",Akurasi$Model[which.min(Akurasi[,"AIC"])])
## [1] "Model yang terbaik adalah model ARIMA (2,1,0)"

Berdasarkan nilai AIC terkecil, dugaan model terbaik sementara adalah ARIMA(2,1,0). Namun, prioritas pertama dalam penentuan model terbaik dapat dilihat dari signifikansi pendugaan parameter model terlebih dahulu sebelum melihat nilai AIC terkecil, hanya ada satu model tentatif dengan semua parameter yang berpengaruh signifikan. Oleh karena itu, dugaan model terbaik berdasarkan pemenuhan kedua kriteria tersebut adalah model ARIMA (2,1,0).

Overfitting

Overfitting dilakukan dengan menambah 1 ordo p atau q, sehingga model perluasan yang digunakan adalah model ARIMA(3,1,0) dan ARIMA(2,1,1). Kedua kandidat dibandingkan dengan model ARIMA (2,1,0), lalu dilakukan pemeriksaan signifikansi model.

model4 <- Arima(training.ts, order=c(3,1,0),method="ML") #ar3 tidak signifikan
lmtest::coeftest(model4)
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.412895   0.099532  4.1483 3.349e-05 ***
## ar2 -0.389999   0.103153 -3.7808 0.0001563 ***
## ar3  0.143884   0.107048  1.3441 0.1789153    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
model5 <- Arima(training.ts, order=c(2,1,1),method="ML") #ar1 ma1 tidak signifikan
lmtest::coeftest(model5)
## 
## z test of coefficients:
## 
##     Estimate Std. Error z value Pr(>|z|)  
## ar1  0.13257    0.24082  0.5505  0.58196  
## ar2 -0.27628    0.12544 -2.2024  0.02764 *
## ma1  0.26780    0.24445  1.0955  0.27328  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Berdasarkan hasil uji signifikansi di atas, kedua kandidat model perluasan overfitting ARIMA (2,1,0) , yaitu ARIMA (3,1,0) dan ARIMA (2,1,1) tidak signifikan, sehingga model yang dipilih untuk tahapan analisis selanjutnya tetap model ARIMA (2,1,0).

Uji Diagnostik

Diagnostik Model: Eksploratif

Analisis sisaan dapat dilihat secara eksploratif menggunakan Q-Q plot, residuals plot, ACF dan PACF plot.

sisaan <- model1$residuals
par(mfrow=c(2,2))
qqnorm(sisaan)
qqline(sisaan, col = "red", lwd =1)
plot(c(1:length(sisaan)),sisaan)
acf(sisaan)
pacf(sisaan)

Berdasarkan hasil eksplorasi menggunakan Q-Q plot, terlihat bahwa sisaan berdistribusi mengikuti garis normal, sehingga dapat dikatakan bahwa sisaan menyebar normal. Kemudian, plot sisaan yang diperoleh menunjukkan bahwa sisaan memiliki pola acak dan tersebar di sekitar nilai nol. Sedangkan pada plot ACF dan PACF, nilai awal amatan tidak melewati garis signifikan, atau dapat dikatakan bahwa sisaan saling bebas.

Diagnostik Model: Uji Formal

1. Sisaan Menyebar Normal

Uji formal ini dilakukan dengan Jarque Bera test dan Shapiro-Wilk test.

Hipotesis yang diuji:
\(H_0\): Sisaan Menyebar Normal
\(H_1\): Sisaan Tidak Menyebar Normal

jarque.bera.test(sisaan)
## 
##  Jarque Bera Test
## 
## data:  sisaan
## X-squared = 45.156, df = 2, p-value = 1.565e-10
shapiro.test(sisaan)
## 
##  Shapiro-Wilk normality test
## 
## data:  sisaan
## W = 0.93848, p-value = 0.000156

Berdasarkan Jarque-Bera test, diperoleh p-value (1.565e-10) < α (0.05). Selain itu, hasil Shapiro-Wilk test, diperoleh p-value (0.000156) < α (0.05) maka tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa sisaan menyebar normal pada taraf nyata 5%. Namun, asumsi normalitas tidak wajib terpenuhi.

2. Sisaan Saling Bebas

Uji formal ini dilakukan dengan LJung-Box test.

Hipotesis yang diuji:
\(H_0\): Sisaan antara lag saling bebas
\(H_1\): Sisaan antara lag tidak saling bebas

Box.test(sisaan, type = "Ljung")
## 
##  Box-Ljung test
## 
## data:  sisaan
## X-squared = 0.0013169, df = 1, p-value = 0.9711

Berdasarkan LJung-Box test, diperoleh p-value (0.9711) > α (0.05), maka tak tolak H0. Artinya, cukup bukti untuk menyatakan bahwa sisaan antara lag saling bebas atau dapat dikatakan tidak ada autokorelasi antara sisaan lag pada taraf nyata 5%.

3. Nilai Tengah Sisaan Sama dengan Nol

Uji formal ini dilakukan dengan t-test.

Hipotesis yang diuji:
\(H_0\): Nilai tengah sisaan sama dengan nol
\(H_1\): Nilai tengah sisaan tidak sama dengan nol

t.test(sisaan, mu = 0, conf.level = 0.95) 
## 
##  One Sample t-test
## 
## data:  sisaan
## t = 2.0969, df = 99, p-value = 0.03856
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##   216.8209 7855.9831
## sample estimates:
## mean of x 
##  4036.402

Berdasarkan t-test, diperoleh p-value (0.03856) < α (0.05), maka tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa nilai tengah sisaan sama dengan nol pada taraf nyata 5%. Namun, asumsi nilai harapan sama dengan nol tidak wajib terpenuhi.

4. Sisaan Homogen

Uji formal ini dilakukan dengan ARCH test.

Hipotesis yang diuji:
\(H_0\): Sisaan homogen
\(H_1\): Sisaan tidak homogen

ArchTest(sisaan)
## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  sisaan
## Chi-squared = 36.928, df = 12, p-value = 0.0002295

Berdasarkan ARCH test, diperoleh p-value (0.0002295) < α (0.05), maka tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa sisaan homogen pada taraf nyata 5%.

Auto ARIMA

model.auto <- auto.arima(training.ts,trace = T,d = 1)
## 
##  ARIMA(2,1,2) with drift         : 2246.161
##  ARIMA(0,1,0) with drift         : 2258.767
##  ARIMA(1,1,0) with drift         : 2254.504
##  ARIMA(0,1,1) with drift         : 2246.609
##  ARIMA(0,1,0)                    : 2260.422
##  ARIMA(1,1,2) with drift         : 2246.65
##  ARIMA(2,1,1) with drift         : 2244.425
##  ARIMA(1,1,1) with drift         : 2246.93
##  ARIMA(2,1,0) with drift         : 2242.92
##  ARIMA(3,1,0) with drift         : 2244.138
##  ARIMA(3,1,1) with drift         : 2245.166
##  ARIMA(2,1,0)                    : 2245.141
## 
##  Best model: ARIMA(2,1,0) with drift
summary(model.auto)
## Series: training.ts 
## ARIMA(2,1,0) with drift 
## 
## Coefficients:
##          ar1      ar2     drift
##       0.3373  -0.3641  4130.944
## s.e.  0.0936   0.0944  1890.488
## 
## sigma^2 = 381662419:  log likelihood = -1117.25
## AIC=2242.49   AICc=2242.92   BIC=2252.88
## 
## Training set error measures:
##                    ME     RMSE      MAE        MPE     MAPE      MASE
## Training set 112.8012 19141.47 13862.54 -0.1805612 2.441015 0.9053651
##                    ACF1
## Training set 0.03682325

Setelah melakukan analisis pemodelan ARIMA dengan metode Box-Jenskins, dapat diperiksa dengan menggunakan fungsi auto.arima(). Kedua metode ini menghasilkan output yang sama dimana model terbaik adalah ARIMA (2,1,0) with drift (konstanta). Model ARIMA (2,1,0) memiliki konstanta sebesar 4130.944 dan konstanta galat bakunya sebesar 1890.488.

Forecasting ARIMA

ramalan <- forecast::forecast(Arima(training.ts, order=c(2,1,0),method="ML",include.drift = TRUE), h=20) 

data.ramalan <- ramalan$mean
plot(ramalan,lwd=2)

ts.plot(mpdw$GoldPrice,col="black",lwd=2,main="Perbandingan Peramalan dan Aktual",ylim=c(400000,1000000))
lines(data.ramalan, col = "red",lwd=2)

Berdasarkan plot di atas, terlihat bahwa hasil forecasting mengalami tren naik.

perbandingan.temp<-matrix(data=c(testing.ts[1:5], data.ramalan[1:5]), nrow = 5, ncol = 2)
colnames(perbandingan.temp)<-c("Aktual","Hasil Forecast")
head(perbandingan.temp)
##        Aktual Hasil Forecast
## [1,] 851021.7       901821.8
## [2,] 842547.3       913629.9
## [3,] 841661.6       920419.4
## [4,] 815544.9       922652.5
## [5,] 794819.7       925175.7

Berdasarkan perbandingan nilai data aktual dan forecasting, terjadi perbedaan nilai data yang cukup signifikan. Hal ini kurang sejalan dengan data aktual yang cenderung menurun, tetapi data forecasting cenderung meningkat, sehingga perbandingan kedua data terpaut cukup jauh.

accuracy(testing.ts,data.ramalan)
##                ME     RMSE      MAE      MPE     MAPE      ACF1 Theil's U
## Test set 106527.7 109463.8 106527.7 11.20824 11.20824 0.5810819  22.32879

Berdasarkan analisis yang telah dilakukan, pendugaan model terbaik pada data harga penutupan emas adalah ARIMA (2,1,0) dilihat dari signifikasi pada semua pendugaan parameter model dan nilai AIC terkecil. Karena berdasarkan uji asumsi, asumsi homogentitas sisaan tidak terpenuhi maka direkomendasikan menggunakan model ARIMA ARCH-GARCH.

Pemodelan ARCH-GARCH

Heteroskedastik juga mempengaruhi akurasi selang kepercayaan peramalan. Efisiensi dan akurasi pendugaan dapat dilakukan dengan menyusun model yang mempertimbangkan adanya heteroskedastik. Keberadaan masalah heteroskedastik dapat diuji menggunakan the Portmanteau Q test atau Engle Lagrange Multiplier tests (LM Test) untuk menguji apakah ragam sisaan adalah fungsi dari kuadrat residual pada waktu-waktu sebelumnya (lag time). Penentuan pemodelan sisaan menggunakan ARCH atau GARCH melihat banyaknya lag/ordo yang nyata. Jika ordo yang signifikan kurang dari lima (singkat), maka menggunakan model ARCH. Sebaliknya, ordo yang nyata berjumlah lebih dari lima (lama), menggunakan model GARCH.

Mean Model

model.mean <- stats::arima(mpdw$GoldPrice, order = c(2,1,0), include.mean = FALSE, method = "ML")
checkresiduals(model.mean)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(2,1,0)
## Q* = 2.6243, df = 8, p-value = 0.9557
## 
## Model df: 2.   Total lags used: 10
kurtosis(model.mean$residuals)
## [1] 5.006638

Berdasarkan LJung-Box test, diperoleh p-value (0.9557) > α (0.05), maka tak tolak H0. Artinya, cukup bukti untuk menyatakan bahwa sisaan saling bebas pada model mean ARIMA (2,1,0) pada taraf nyata 5%, sama halnya pada visualisasi ACF plot. Selain itu, berdasarkan time series plot, sisaan menyebar heterogen. Berdasarkan grafik histogram, sisaan cenderung leptokurtic (nilai kurtosis (5.006638) > 3) yang mana data banyak mengumpul di titik nol.

Pemeriksaan efek ARCH

ehatsq <- ts(resid(model1)^2)
goldprice.ARCH <- dynlm(ehatsq~L(ehatsq))
summary(goldprice.ARCH)
## 
## Time series regression with "ts" data:
## Start = 2, End = 100
## 
## Call:
## dynlm(formula = ehatsq ~ L(ehatsq))
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -627652986 -316037792 -258187477  -33624055 5366806050 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)   
## (Intercept) 3.181e+08  9.551e+07   3.331  0.00122 **
## L(ehatsq)   1.785e-01  9.985e-02   1.788  0.07689 . 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 869600000 on 97 degrees of freedom
## Multiple R-squared:  0.03191,    Adjusted R-squared:  0.02193 
## F-statistic: 3.197 on 1 and 97 DF,  p-value: 0.07689

Diperoleh p-value (0.07689) > α (0.05), maka tak tolak H0. Artinya, output di atas menunjukkan tidak terdapat efek ARCH pada sisaan mean model pada ARIMA (2,1,0) sehingga tidak cocok dilanjutkan dengan melakukan pemodelan sisaan melalui variance model. Pemodelan ARIMAX akan dicobakan pada analisis selanjutnya.

Proses Short or Long Memory

aTSA::arch.test(model.mean)
## ARCH heteroscedasticity test for residuals 
## alternative: heteroscedastic 
## 
## Portmanteau-Q test: 
##      order   PQ  p.value
## [1,]     4 28.8 8.40e-06
## [2,]     8 34.3 3.60e-05
## [3,]    12 37.8 1.65e-04
## [4,]    16 39.1 1.06e-03
## [5,]    20 41.4 3.31e-03
## [6,]    24 44.1 7.45e-03
## Lagrange-Multiplier test: 
##      order    LM  p.value
## [1,]     4 77.14 1.11e-16
## [2,]     8 12.20 9.42e-02
## [3,]    12  7.00 7.99e-01
## [4,]    16  4.78 9.94e-01
## [5,]    20  3.09 1.00e+00
## [6,]    24  2.16 1.00e+00

Berdasarkan the Portmanteau Q test diperoleh p-value < α hingga lag ke-24, maka tolak H0. Artinya, ada indikasi proses long memory dimana terdapat ordo yang nyata sampai lag ke-24 sehingga penggunaan model ARCH menjadi kurang tepat. Jadi, variance modelakan dilanjutkan dengan menggunakan pemodelan GARCH.

for (i in 1:24) {
  ArchTest <- ArchTest(model.mean$residuals, lags=i, demean=TRUE)
  cat("P Value LM Test lag ke", i,"adalah" , ArchTest$p.value, "\n") }
## P Value LM Test lag ke 1 adalah 0.1586672 
## P Value LM Test lag ke 2 adalah 0.1453683 
## P Value LM Test lag ke 3 adalah 0.1953917 
## P Value LM Test lag ke 4 adalah 2.674232e-05 
## P Value LM Test lag ke 5 adalah 1.662288e-05 
## P Value LM Test lag ke 6 adalah 2.658517e-05 
## P Value LM Test lag ke 7 adalah 9.573099e-06 
## P Value LM Test lag ke 8 adalah 1.432823e-05 
## P Value LM Test lag ke 9 adalah 2.786684e-05 
## P Value LM Test lag ke 10 adalah 6.069207e-05 
## P Value LM Test lag ke 11 adalah 0.0001309282 
## P Value LM Test lag ke 12 adalah 0.0002635894 
## P Value LM Test lag ke 13 adalah 0.0005519993 
## P Value LM Test lag ke 14 adalah 0.0009063845 
## P Value LM Test lag ke 15 adalah 0.001500025 
## P Value LM Test lag ke 16 adalah 0.001676533 
## P Value LM Test lag ke 17 adalah 0.002921707 
## P Value LM Test lag ke 18 adalah 0.004557954 
## P Value LM Test lag ke 19 adalah 0.007711867 
## P Value LM Test lag ke 20 adalah 0.01097896 
## P Value LM Test lag ke 21 adalah 0.01526735 
## P Value LM Test lag ke 22 adalah 0.02411995 
## P Value LM Test lag ke 23 adalah 0.03331705 
## P Value LM Test lag ke 24 adalah 0.05032072

Simpulan

….

Daftar Pustaka

Aktivani S. 2020. Uji stasioneritas data inflasi kota padang periode 2014-2019. Jurnal Statistika. 20(2):83-90.

Anderson DR, Burnham KP, White GC. 1998. Comparison of akaike information criterion and consistent akaike information criterion for model selection and statistical inference from capture-recapture studies. Journal of Applied Statistics. 25(2): 263-282.

Ariwibowo A. 2020 Mei 14. Analisis value at risk menggunakan pendekatan threshold generalized autoregressive conditional heteroscedasticity dan generalized pareto distribution (Studi kasus: harga penutupan saham PT Hotel X Tbk. Periode: 2 Januari 2017 - 21 November 2019) (skripsi). Jakarta: UIN Syarif Hidayatullah Jakarta. https://repository.uinjkt.ac.id/dspace/handle/123456789/51563.

Dewi DM, Nafi MZ, Nasrudin. 2022. Analisis peramalan harga emas di Indonesia pada masa pandemi Covid-19 untuk investasi. Jurnal Litbang Sukowati. 5(10): 38—50.

Gempur. 2010. Rata-rata bergerak tunggal(single moving average) [internet]. [diacu 2022 Oktober 20]. Tersedia dari: https://exponensial.wordpress.com/tag/rata-rata-bergerak-tunggal/#:~:text=Rata-rata%20bergerak%20tunggal%20%28Moving%20average%29%20untuk%20periode%20t,ini%20digunakan%20untuk%20memprediksi%20nilai%20pada%20periode%20berikutnya.

Ginantra NLWSR, Anandita IBG. 2019. Penerapan metode single exponential smoothing dalam peramalan penjualan barang. Jurnal Sains Komputer & Informatika (J-SAKTI). 3(2): 433—441.

Habibi MA, Wardhani LP. 2018. Estimasi parameter pada model negatif binomial generalized autoregressive moving average (GARMA) dengan algoritma IRLS (Studi kasus peramalan jumlah kecelakaan di jalan tol Gempol-Surabaya). Jurnal Sains dan Seni ITS. 7(2): 46-51.

Hamidah SN, Salam N, Susanti DS. 2013. Teknik peramalan menggunakan metode pemulusan eksponensial holt-winters. Jurnal Matematika Murni dan Terapan “Epsilon”. 7(2): 26—33.

Hidayat A. 2013. Penjelasan tentang uji normalitas dan metode perhitungan [internet]. [diacu 2022 Oktober 20]. Tersedia dari: https://www.statistikian.com/2013/01/uji-normalitas.html.

Mardiatmoko G. 2020. Pentingnya uji asumsi klasik pada analisis regresi linier berganda (Studi kasus penyusunan persamaan allometrik kenari muda (Canarium indicum L.)). Jurnal Ilmu Matematika Terapan. 14(3):333-342.

Mulyana. 2004. Buku Ajar Analisis Deret Waktu. Bandung: Universitas Padjajaran.

Muzaki L. 2022. Contoh studi kasus metode double moving average, kelebihan kekurangan, dan cara menghitung [internet]. [diacu pada 2022 Oktober 20]. Tersedia dari: https://www.pengadaanbarang.co.id/2022/05/metode-double-moving-average.html.

Qudratullah, Farhan M. 2007. Bayesian information criterion (BIC) dalam pemilihan model terbaik feed forward neural network (FFNN) : Studi kasus data posisi dana simpanan tabungan bank Umum dan BPR di provinsi D.I. Yogyakarta [skripsi]: Universitas Gajah Mada.

Salwa N, Tatsara N, Amalia R, Zohra AF. 2018. Peramalan harga bitcoin menggunakan metode arima (autoregressive integrated moving average). Journal of Data Analysis. 1(1):21-31. doi:10.24815/jda.v1i1.11874.

Setyowati OAD. 2020. Peramalan harga cabai rawit di Provinsi Jawa Timur menggunakan metode arimax (skripsi). Surabaya: UIN Sunan Ampel. http://digilib.uinsby.ac.id/43026/.

Sui M, Rengifo EW, Court E. 2021. Gold, inflation and exchange rate in dollarized economies–A comparative study of Turkey, Peru and the United States. International Review of Economics & Finance. 71: 82—99. doi:10.1016/j.iref.2020.08.014.