Import Data

Data yang digunakan adalah data Nilai Impor (Juta USD) bulanan pada rentang tahun 2010 - 2019 dengan jumlah 120 amatan yang bersumber dari Tabel Dinamis Subjek Ekspor-Impor BPS https://www.bps.go.id

Nilai impor merupakan nilai suatu kelompok barang yang diimpor (melewati batas negara Indonesia) tanpa menggunakan dokumen non PEB/PIB, dalam satuan US Dollar (USD).

#Load Package
library(readxl)
library(forecast)

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

library(TTR)
library(TSA)

## Registered S3 methods overwritten by 'TSA':
##   method       from    
##   fitted.Arima forecast
##   plot.Arima   forecast

## 
## Attaching package: 'TSA'

## The following objects are masked from 'package:stats':
## 
##     acf, arima

## The following object is masked from 'package:utils':
## 
##     tar

#---Input data---#
data1 <- read_excel("C:/Users/User/Documents/DIAH/SEMESTER 5/STA1341 METODE PERAMALAN DERET WAKTU/Impor.xlsx", sheet=3)
impor <- data1$Impor
# Membentuk Objek Time Series
impor.ts<-ts(data1$Impor, frequency=12, start=2010)

Pertemuan 5 - Cek Kestasioneran Data

### Eksplorasi Data - Stasioneritas Data
#Plot time series
plot(impor.ts,
     col = "navyblue",
     lwd = 1,
     type = "o",
     xlab = "Periode",
     ylab = "Nilai Impor (Juta USD)",
     main = "Time Series Plot Nilai Impor 2010 - 2019")

Plot time series di atas menunjukkan pergerakan nilai impor bulanan nasional dalam jutaan USD sejak 2010 hingga 2019 yang selalu berubah baik diamati berdasarkan bulan maupun tahun. Time series plot tersebut juga menunjukkan bahwa nilai Impor bulanan 2010-2019 terlihat tidak stasioner baik dalam nilai tengah maupun ragam dan cenderung terlihat memiliki pola siklik. Nilai Impor ini dapat dipengaruhi oleh fluktuasi ekonomi jangka panjang seperti yang berhubungan dengan siklus bisnis.

Untuk memperkuat argumen mengenai ketidakstasioneran data nilai Impor 2010-2019, plot ACF dan uji formal Augmented Dicky-Fuller Test dapat digunakan sebagai pertimbangan tambahan.

#Melihat Kestasioenran Data menggunakan plot ACF
library(ggplot2)
acf(data1$Impor, lag.max = 25, plot = F) %>% 
  autoplot() + 
  ggtitle("Plot ACF Nilai Impor 2010-2019") + 
  theme_classic()

Berdasarkan plot ACF di atas terlihat bahwa autokorelasi antar lag menurun secara perlahan. Hal ini mengindikasikan bahwa data tidak stasioner. Selanjutnya akan dilakukan uji formal menggunakan Augmented Dicky-Fuller Test untuk memastikan hal tersebut.

Hipotesis yang digunakan dalam uji ini yaitu:

H0: Data tidak stasioner

H1: Data stasioner

Penolakan terhadap H0 dapat diputuskan bilai p-value yang dihasilkan uji Augmented Dicky-Fuller lebih kecil dibandingkan alpha (5%).

Hasil perhitungan ADF dapat dilihat sebagai berikut.

#Uji formal - Augmented Dicky-Fuller Test
library(tseries)
adf.test(impor.ts)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  impor.ts
## Dickey-Fuller = -2.1404, Lag order = 4, p-value = 0.5182
## alternative hypothesis: stationary

Karena nilai p-value > 0.05, maka belum cukup bukti untuk menolak H0. Oleh karena itu dapat dikatakan bahwa data memang tidak stasioner.

Untuk membuat kesimpulan statistik tentang struktur proses stokastik berdasarkan catatan yang diamati dari proses itu, beberapa asumsi penyederhanaan (dan mungkin masuk akal) tentang struktur terkadang diperlukan. Salah satu asumsi paling penting adalah asumsi stasioneritas data (Cyrer JD dan Chan KS. 2008). Ketidakstasionaran data impor 2010-2019 di atas harus ditangani. Salah satu metode untuk menangani masalah ini adalah metode differencing. Hasil differencing pertama data nilai impor dapat dilihat di bawah ini.

#differencing pertama
impor.dif1<-diff(impor, differences = 1)
impor.dif1 <- ts(impor.dif1, frequency = 12, start=2010)

#identifikasi stasioneritas
plot.ts(impor.dif1,lty=1,xlab="Periode", ylab="Nilai Impor Differencing 1")

acf(impor.dif1, lag.max = 25, plot = F) %>% 
  autoplot() + 
  ggtitle("Plot ACF Diff 1") + 
  theme_classic()

adf.test(impor.dif1)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  impor.dif1
## Dickey-Fuller = -6.6573, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary

Pada plot time-series terlihat bahwa data yang sudah dilakukan differencing sudah menyebar di sekitar suatu konstanta sumbu y, sehingga dapat dikatakan bahwa data sudah terlihat stasioner. Akan tetapi, dari plot tersebut masih diragukan apabila ragam pada data tersebut sudah homogen.

Berdasarkan plot ACF terlihat bahwa autokorelasi antar lag sudak tidak menurun secara perlahan namun berupa tails off. Begiru pula dengan hasil pengujian ADF memperlihatkan bahwa data sudah stasioner (p-value < alpha). Dengan tiga pertimbangan ini, maka dapat disimpulkan bahwa data hasil differencing sudah stasioner.

Pertemuan 6 - Identifikasi p, d, q

Penentuan parameter model ARIMA dapat dilihat dari beberapa sudut pandang. Untuk parameter p dan q dapat dilihat dari plot ACF, plot PACF, dan plot EACF. Apabila dilakukan differencing ordo 1, maka nilai \(d=1\), dan seterusnya. Sebelum memulai pemodelan tersebut, data dibagi menjadi data training dan data testing. Partisi yang akan digunakan yaitu 80% (96 amatan) untuk data training dan 20% (24 amatan) untuk data testing.

training<-data1[1:96,3]
testing<-data1[97:120,3]

# Membentuk Objek Time Series
training.ts<-ts(training, frequency=12, start=2010)
testing.ts<-ts(testing,frequency=12, start=2018)

## Plot Data Training dan testing
ts.plot(impor.ts, xlab = "Periode", ylab ="Nilai Impor (Juta USD)", 
        main = "Plot Data Training dan Data Testing")
lines(training.ts, col = "blue")
lines(testing.ts, col="Red")
legend(2015,18000,c("Data Training","Data Testing"), 
       lty=8, col=c("blue","red"), cex=0.8)

Pemodelan akan didasarkan pada data training, kemudian validasi model akan didasarkan pada testing. Karena data dalam keadaan tidak stasioner, maka data training dan testing akan dikenakan proses differencing. Kemudian, plot ACF, PACF, dan EACF dari data training dapat dilihat sebagai berikut.

train_ts_diff <- diff(training.ts, differences=1)
test_ts_diff <- diff(testing.ts, differences=1)

acf.diff <- 
  acf(train_ts_diff, lag.max = 20, plot = F) %>% 
  autoplot() +
  ggtitle(NULL) +
  theme_classic()

pacf.diff <- 
  pacf(train_ts_diff, lag.max = 20, plot = F) %>% 
  autoplot() +
  ggtitle(NULL) +
  theme_classic()

vip::grid.arrange(acf.diff, pacf.diff, ncol = 2)

eacf(train_ts_diff)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x o o o o o o o o x o  x  x  o 
## 1 x o o o o o o o o x o  x  o  o 
## 2 x x o o o o o o o o o  x  o  o 
## 3 x o o o o o o o o o o  x  o  o 
## 4 o o o o o o o o o o o  x  o  o 
## 5 o x o o o o o o o o o  x  o  o 
## 6 x x o o o o o o o o o  x  o  o 
## 7 x x o o o o o o o o o  x  o  o

Berdasarkan plot ACF dan PACF, terlihat bahwa kedua plot seperti tails off sehingga model tentatif yang dapat terbentuk adalah ARIMA(p,1,q). Dari plot EACF terlihat bahwa model tentatif yang dapat terbentuk adalah model ARIMA(0,1,1), ARIMA(0,1,2), ARIMA(1,1,1), ARIMA(1,1,2), dan ARIMA(2,1,2). Untuk melihat model tentatif yang paling efektif untuk dilakukan peramalan, dapat dilihat terlebih dahulu model yang memiliki seluruh koefisien yang signifikan.

mod1.impor <- Arima(training.ts, order = c(0,1,1), method = "ML")
mod2.impor <- Arima(training.ts, order = c(0,1,2), method = "ML")
mod3.impor <- Arima(training.ts, order = c(1,1,1), method = "ML")
mod4.impor <- Arima(training.ts, order = c(1,1,2), method = "ML")
mod5.impor <- Arima(training.ts, order = c(2,1,2), method = "ML")
lmtest::coeftest(mod1.impor)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ma1 -0.58539    0.06974 -8.3939 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lmtest::coeftest(mod2.impor)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ma1 -0.66810    0.11408 -5.8567 4.723e-09 ***
## ma2  0.10134    0.10122  1.0012    0.3167    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lmtest::coeftest(mod3.impor)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ar1 -0.12591    0.14411 -0.8737    0.3823    
## ma1 -0.51952    0.11310 -4.5935 4.359e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lmtest::coeftest(mod4.impor)

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value  Pr(>|z|)    
## ar1  0.936690   0.058464  16.0216 < 2.2e-16 ***
## ma1 -1.689110   0.086425 -19.5443 < 2.2e-16 ***
## ma2  0.758193   0.078201   9.6954 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lmtest::coeftest(mod5.impor)

## 
## z test of coefficients:
## 
##     Estimate Std. Error  z value  Pr(>|z|)    
## ar1  0.69683    0.15101   4.6144 3.942e-06 ***
## ar2 -0.34765    0.16838  -2.0647   0.03895 *  
## ma1 -1.36495    0.10030 -13.6082 < 2.2e-16 ***
## ma2  0.77715    0.13571   5.7267 1.024e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Hanya model ARIMA(0,1,1), ARIMA(1,1,2), dan ARIMA(2,1,2) yang memiliki seluruh parameter signifikan dalam model yang ditunjukkan oleh nilai pr(>|z|) < alpha (5%). Sehingga model-model tersebut yang akan dipertimbangkan untuk digunakan dalam peramalan.