Pada kesempatan ini akan dibahas terkait tutorial pemodelan ARIMA menggunakan R. Adapun tahapan yang dilakukan pada kesempatan ini meliputi:

  1. Spesifikasi Model

  2. Estimasi Model

  3. Diagnotstik Model

  4. Forecasting

Spesifikasi Model

Spesifikasi model dimaksudkan untuk mengetahui gambaran umum data. Gambaran umum data meliputi indektifikasi pola data dan uji stasioneritas data. Selain itu, pada tahapan ini akan dilakukan spesifikasi kandidat model yang sesuai dengan data menggunakan plot ACF, PACF, dan EACF.

library(TSA)
## 
## Attaching package: 'TSA'
## The following objects are masked from 'package:stats':
## 
##     acf, arima
## The following object is masked from 'package:utils':
## 
##     tar
library(readxl)
data.ts<- read_excel("homework_s12.xlsx", sheet = 1)
plot.ts(data.ts[,3])

#Uji FOrmal Stasioneritas Data
tseries::adf.test(ts(data.ts[,3]))
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts(data.ts[, 3])
## Dickey-Fuller = -4.8408, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary
#Spesifikasi Model 
acf(ts(data.ts[,3]))

pacf(ts(data.ts[,3]))

eacf(ts(data.ts[,3]))
## 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 o o  o  o  o 
## 1 x o o o o o o o o o o  o  o  o 
## 2 x o o o o o o o o o o  o  o  o 
## 3 x x x o o o o o o o o  o  o  o 
## 4 x x x o o o o o o o o  o  o  o 
## 5 x x o x o o o o o o o  o  o  o 
## 6 x o x x o o o o o o o  o  x  o 
## 7 x x x o o o o o o o o  o  o  o

Stasioneritas data

Berdasarkan plot data dan ACF (tidak turun secara melandai) menunjukkan bahwa data stasioner. Selanjutnya, berdasarkan hasil uji formal, diperoleh nilai p-value < 0.05 yang berarti TOLAK H0. Penolakan H0 menunjukkan bahwa data stasioner.

Spesifikasi model

Berdasarkan hasil di atas, didatkan hasil pada masing-masing plot sebagai berikut:

ACF(tails off); PACF(cut off after lag 1); EACF(p=0, q=1 dan p=1, q=1)

Kandidat model untuk data 3 yaitu AR(1), ARMA(1,1) dan MA(1)

Model Estimation

model.ar<-arima(ts(data.ts[,3]), order = c(1,0,0), method = "ML")
model.arma<-arima(ts(data.ts[,3]), order = c(1,0,1), method = "ML")
model.ma<-arima(ts(data.ts[,3]), order = c(0,0,1), method = "ML")
#Pengujian parameter
lmtest::coeftest(model.ar)
## 
## z test of coefficients:
## 
##             Estimate Std. Error   z value  Pr(>|z|)    
## ar1        -0.333613   0.095612   -3.4892 0.0004844 ***
## intercept 100.155529   0.074569 1343.1344 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lmtest::coeftest(model.ma)
## 
## z test of coefficients:
## 
##             Estimate Std. Error   z value  Pr(>|z|)    
## ma1        -0.272846   0.086518   -3.1536  0.001613 ** 
## intercept 100.154904   0.073277 1366.7974 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lmtest::coeftest(model.arma)
## 
## z test of coefficients:
## 
##             Estimate Std. Error   z value Pr(>|z|)    
## ar1        -0.504031   0.214907   -2.3453  0.01901 *  
## ma1         0.188410   0.234172    0.8046  0.42106    
## intercept 100.155526   0.078294 1279.2155  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Hasil Pengujian paramater pada masing-maasing model menunjukkan hasil yang signifikan pada model AR(1) dan MA(1). Sementara itu, pada model ARMA(1,1) menunjukkan hasil bahwa parameter signifikan pada ar1. Hal tersebut menunjukkan bahwa kandidat model yang tepat untuk data 3 yaitu AR(1) dan MA(1). Selanjutnya, untuk mengetahui model terbaik ditentukan berdasarkan nilai AIC dan hasil pengujian asumsi.

Diagnostik Model

Ukuran kebaikan model yang digunakan pada tutorial ini yaitu AIC.

(aic.model<- data.frame(Model = c("AR(1)","MA(1)","ARMA(1,1)"),
                       AIC = c(model.ar$aic, model.ma$aic, model.arma$aic)))
##       Model      AIC
## 1     AR(1) 286.2491
## 2     MA(1) 288.5843
## 3 ARMA(1,1) 287.6902

Hasil perhitungan AIC menunjukkan bahwa nilai AIC model AR(1) lebih kecil dibandingkan nilai AIC model lainnya. Hal tersebut menunjukkan bahwa AR(1) adalah model terbaik. Selanjutnya akan dilakukan diagnostik model pada AR(1). Diagnostik model dimaksudkan untuk mengetahui apakah sisaan dari model yang dibangun telah memenuhi asumsi pemodelan. Adapun asumsi sisaan yang harus terpenuhi adalah sebagai berikut:

  1. Sisaan mengikuti sebaran normal n(0,\(\sigma^2\))

  2. Sisaan saling bebas (tidak ada autokorelasi)

Pengujian sisaan pada tutorial ini akan dilakukan secara eksploratif dan uji formal. Berikut hasil pengujian secara esdploratif.

sisaan <- model.ar$residuals
# Eksplorasi
par(mfrow=c(2,2))
qqnorm(sisaan)
qqline(sisaan, col = "blue", lwd = 2)
plot(c(1:length(sisaan)),sisaan)
acf(sisaan)
pacf(sisaan)

Berdasarkan hasil diatas, secara eksploratif terlihat bahwa sisaan mengikuti sebaran normal dengan nilai tengah sama dengan 0. Selanjutnya, ditinjau dari plot ACF dan PACF terlihat bahwa tidak ada lag yang signifikan. Hal tersebut menunjukkan bahwa tidak ada gejala autokorelasi pada sisaan. Selanjutnya, untuk memastikan kembali akan dilakukan uji asumsi secara formal.

# Uji formal normalitas data
ks.test(sisaan,"pnorm")
## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  sisaan
## D = 0.047946, p-value = 0.9756
## alternative hypothesis: two-sided
# Uji nilai tengah sisaan
t.test(sisaan, mu = 0, alternative = "two.sided")
## 
##  One Sample t-test
## 
## data:  sisaan
## t = -0.043639, df = 99, p-value = 0.9653
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.2021241  0.1934248
## sample estimates:
##    mean of x 
## -0.004349634
# Uji autokorelasi
Box.test(sisaan, lag = 23 ,type = "Ljung")
## 
##  Box-Ljung test
## 
## data:  sisaan
## X-squared = 19.149, df = 23, p-value = 0.6926
# Overfitting model
(model.ar2<- arima(ts(data.ts[,3]), order = c(2,0,0), method = "ML"))
## 
## Call:
## arima(x = ts(data.ts[, 3]), order = c(2, 0, 0), method = "ML")
## 
## Coefficients:
##           ar1     ar2  intercept
##       -0.3067  0.0867   100.1553
## s.e.   0.1003  0.1021     0.0811
## 
## sigma^2 estimated as 0.9764:  log likelihood = -140.77,  aic = 287.53
lmtest::coeftest(model.ar2)
## 
## z test of coefficients:
## 
##             Estimate Std. Error   z value  Pr(>|z|)    
## ar1        -0.306679   0.100322   -3.0569  0.002236 ** 
## ar2         0.086655   0.102144    0.8484  0.396233    
## intercept 100.155281   0.081099 1234.9755 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Hasil pengujian:

Normalitas Data

\(H_0\) : sisaan menyebar normal

\(H_1\) : sisaan tidak mengikuti sebaran normal

Hasil : \(p-value = 0.9756 > \alpha = 0.05\) yang berarti TAK TOLAK \(H_0\). Sisaan menyebar normal

Nilai Tengah Sisaan

\(H_0 : \mu = 0\)

\(H_1 : \mu \neq 0\)

Hasil : \(p-value = 0.9653 > \alpha = 0.05\) yang berarti TAK TOLAK \(H_0\). Nilai tengah sisaan sama dengan 0

Autokorelasi

\(H_0\) : tidak ada autokorelasi

\(H_1\) : terdapat autokorelasi

Hasil : \(p-value = 0.693 > \alpha = 0.05\) yang berarti TAK TOLAK \(H_0\). Tidak terdapat gejala autokorelasi

Kesimpulan : Asumsi terpenuhi

Forecasting

#Peramalan
(ramalan<- forecast::forecast(ts(data.ts[,3]),model=model.ar, h=10))
## Registered S3 methods overwritten by 'forecast':
##   method       from
##   fitted.Arima TSA 
##   plot.Arima   TSA
##     Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 101       100.5811 99.31008 101.8520 98.63727 102.5249
## 102       100.0136 98.67372 101.3534 97.96445 102.0627
## 103       100.2029 98.85560 101.5502 98.14239 102.2634
## 104       100.1397 98.79161 101.4878 98.07796 102.2015
## 105       100.1608 98.81259 101.5090 98.09889 102.2227
## 106       100.1538 98.80555 101.5020 98.09185 102.2157
## 107       100.1561 98.80790 101.5043 98.09419 102.2180
## 108       100.1553 98.80711 101.5036 98.09341 102.2173
## 109       100.1556 98.80737 101.5038 98.09367 102.2175
## 110       100.1555 98.80729 101.5037 98.09358 102.2174
plot(ramalan)