Pendugaan Parameter, Diagnostik Model, dan Peramalan
Kheni Hikmah Lestari
Packages yang digunakakan
library(ggplot2)
library(tsibble)##
## Attaching package: 'tsibble'## The following objects are masked from 'package:base':
##
## intersect, setdiff, unionlibrary(tseries)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(MASS)
library(forecast)
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':
##
## tarlibrary(TTR)
library(aTSA)##
## Attaching package: 'aTSA'## The following object is masked from 'package:forecast':
##
## forecast## The following objects are masked from 'package:tseries':
##
## adf.test, kpss.test, pp.test## The following object is masked from 'package:graphics':
##
## identifylibrary(graphics)
library(rio)Penyiapan Data
#Input data
data = import("https://raw.githubusercontent.com/khenihikmah/praktikummpdw/main/Pertemuan1/Data%20Penumpang%20Pesawat.csv")
data## Month Passengers
## 1 1 13224
## 2 2 10789
## 3 3 6434
## 4 4 6433
## 5 5 4002
## 6 6 3540
## 7 7 7214
## 8 8 5294
## 9 9 4570
## 10 10 7315
## 11 11 7210
## 12 12 7801
## 13 13 10677
## 14 14 13367
## 15 15 8047
## 16 16 7286
## 17 17 5440
## 18 18 6423
## 19 19 10193
## 20 20 8830
## 21 21 8170
## 22 22 10606
## 23 23 9898
## 24 24 12807
## 25 25 18810
## 26 26 14584
## 27 27 8525
## 28 28 7987
## 29 29 7595
## 30 30 8146
## 31 31 14233
## 32 32 11762
## 33 33 10758
## 34 34 13648
## 35 35 11995
## 36 36 12600
## 37 37 16191
## 38 38 13642
## 39 39 14349
## 40 40 14261
## 41 41 12606
## 42 42 12565
## 43 43 15547
## 44 44 13013
## 45 45 12823
## 46 46 15167
## 47 47 17486
## 48 48 18651
## 49 49 21201
## 50 50 22568
## 51 51 19143
## 52 52 18401
## 53 53 15417
## 54 54 14021
## 55 55 17872
## 56 56 14866
## 57 57 16354
## 58 58 15542
## 59 59 15501
## 60 60 17671
## 61 61 28195
## 62 62 23112
## 63 63 15725
## 64 64 14801
## 65 65 13887
## 66 66 12924
## 67 67 17184
## 68 68 14702
## 69 69 13758
## 70 70 16287
## 71 71 14964
## 72 72 22199
## 73 73 27937
## 74 74 28596
## 75 75 19387
## 76 76 15145
## 77 77 14393
## 78 78 13101
## 79 79 16482
## 80 80 14852
## 81 81 16740
## 82 82 21313
## 83 83 20606
## 84 84 24804
## 85 85 35429
## 86 86 36264
## 87 87 23066
## 88 88 23386
## 89 89 17894
## 90 90 17799
## 91 91 22409
## 92 92 16729
## 93 93 17522
## 94 94 21301
## 95 95 24107
## 96 96 32365
## 97 97 45646
## 98 98 37151
## 99 99 29907
## 100 100 30467
## 101 101 24784
## 102 102 22659
## 103 103 33972
## 104 104 26180
## 105 105 23997
## 106 106 30318
## 107 107 26439
## 108 108 32959
## 109 109 43647
## 110 110 45129
## 111 111 33745
## 112 112 29896
## 113 113 22655
## 114 114 21210
## 115 115 32686
## 116 116 27857
## 117 117 23826
## 118 118 31235
## 119 119 30048
## 120 120 37872
## 121 121 50330
## 122 122 47541
## 123 123 31185
## 124 124 31263
## 125 125 32392
## 126 126 29729
## 127 127 36012
## 128 128 33170
## 129 129 30961
## 130 130 38057
## 131 131 37204
## 132 132 39130#Mengambil data Passengers
data <- data$Passengers
#Mengubah tipe data menjadi time series
data.ts<-ts(data)Eksplorasi Data
Plot Data Penuh
plot.ts(data.ts, xlab="waktu", ylab="Jumlah Passengers", main="Plot Data Passengers")
Berdasarkan plot data deret waktu, terlihat bahwa data memiliki pola
siklik. Berdasarkan pola data, pembagian data latih dan data uji
ditetapkan dengan proporsi 80%:20%.
Plot Data Latih
datatrain<-data[1:106]
train.ts<-ts(datatrain)
plot.ts(train.ts, lty=1, xlab="waktu", ylab="Jumlah Passengers", main="Plot Data Passengers")
Berdasarkan plot data deret waktu pada data latih, terlihat bahwa data cenderung memiliki trend yang naik. Hal ini mengindikasikan bahwa data tidak stasioner dalam rataan.
Plot Data Uji
datatest<-data[107:132]
test.ts<-ts(datatest)
plot.ts(test.ts, lty=1, xlab="waktu", ylab="Jumlah Passengers", main="Plot Data Passengers")
Uji Stasioneritas Data
Plot ACF
acf(train.ts)
Uji ADF
tseries::adf.test(train.ts)## Warning in tseries::adf.test(train.ts): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: train.ts
## Dickey-Fuller = -4.6475, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary\(H_0\) : Data tidak stasioner dalam rataan
\(H_1\) : Data stasioner dalam rataan
Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.
Plot Box-Cox
index <- seq(1:106)
bc = boxcox(train.ts~index, lambda = seq(-1,1,by=0.01))
#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda## [1] 0.24#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]## [1] 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.11 0.12 0.13 0.14 0.15 0.16
## [16] 0.17 0.18 0.19 0.20 0.21 0.22 0.23 0.24 0.25 0.26 0.27 0.28 0.29 0.30 0.31
## [31] 0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.40 0.41 0.42 0.43 0.44 0.45 0.46Plot Boxcox menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 0.24 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0.02 dan batas atas 0.46. Selang tersebut tidak memuat nilai satu sehingga dapat dikatakan bahwa data yang digunakan tidak stasioner dalam ragam.
Penanganan Ketidakstasioneran Data
train.diff<-diff(train.ts,differences = 1)
plot.ts(train.diff, lty=1, xlab="waktu", ylab="Jumlah Passengers", main="Plot Data Passengers")
Berdasarkan plot data deret waktu, terlihat bahwa data sudah stasioner dalam rataan ditandai dengan data bergerak pada nilai tengah tertentu (tidak terdapat trend ataupun musiman pada data)
Plot ACF
acf(train.diff)
Berdasarkan plot tersebut, terlihat bahwa plot ACF cuts off pada lag ke 2. Hal ini menandakan data sudah stasioner dalam rataan dan ketidakstasioneran data telah berhasil tertangani.
Uji ADF
tseries::adf.test(train.diff)## Warning in tseries::adf.test(train.diff): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: train.diff
## Dickey-Fuller = -6.7933, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary\(H_0\) : Data tidak stasioner dalam rataan
\(H_1\) : Data stasioner dalam rataan
Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) atau data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF, sehingga dalam hal ini ketidakstasioneran data sudah berhasil ditangani dan dapat dilanjutkan ke pemodelan
Identifikasi Model
Plot ACF
acf(train.diff)
Berdasarkan plot tersebut, terlihat bahwa plot ACF cenderung cuts off pada lag ke 2, sehingga jika plot PACF dianggap tails of, maka model tentatifnya adalah ARIMA(0,1,2).
Plot PACF
pacf(train.diff)
Berdasarkan plot tersebut, terlihat bahwa plot PACF cenderung cuts off pada lag ke 2, sehingga jika plot ACF dianggap tails of, maka model tentatifnya adalah ARIMA(2,1,0).
Jika baik plot ACF maupun plot PACF keduanya dianggap tails of, maka model yang terbentuk adalah ARIMA(1,1,1)
Plot EACF
eacf(train.diff)## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o x o x o x o x o o o x o o
## 1 o x o x o o o x o o o x o o
## 2 x x o x x o x x o o o x x x
## 3 x x x x o o o x o o o x o o
## 4 x x x x o o o x o o o x x o
## 5 x o x o x o o o o o o o x o
## 6 x o x x o o o o o o o o o o
## 7 x o x o o o x o o o o o o oIdentifikasi model menggunakan plot EACF dilakukan dengan melihat ujung segitiga pada pola segitiga nol. Dalam hal ini model tentatif yang terbentuk adalah ARIMA(0,1,2), ARIMA(1,1,0), dan ARIMA(0,1,2)
Pendugaan Parameter Model Tentatif
ARIMA(1,1,0)
model1.da=Arima(train.diff, order=c(1,1,0),method="ML")
summary(model1.da) #AIC=2093.26## Series: train.diff
## ARIMA(1,1,0)
##
## Coefficients:
## ar1
## -0.3573
## s.e. 0.0920
##
## sigma^2 = 31311723: log likelihood = -1044.63
## AIC=2093.26 AICc=2093.38 BIC=2098.55
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 85.432 5542.139 4090.014 -3457.36 3766.799 0.8999882 -0.1836536lmtest::coeftest(model1.da) #seluruh parameter signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.357317 0.092002 -3.8838 0.0001028 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ARIMA(0,1,2)
model2.da=Arima(train.diff, order=c(0,1,2),method="ML")
summary(model2.da) #AIC=2044.75 ## Series: train.diff
## ARIMA(0,1,2)
##
## Coefficients:
## ma1 ma2
## -0.9431 -0.0569
## s.e. 0.1493 0.1471
##
## sigma^2 = 18646470: log likelihood = -1019.37
## AIC=2044.75 AICc=2044.99 BIC=2052.68
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 301.2277 4256.021 3083.701 -2626.017 2809.041 0.6785539
## ACF1
## Training set -0.02564734lmtest::coeftest(model2.da) #ada parameter yang tidak signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.943075 0.149292 -6.3170 2.667e-10 ***
## ma2 -0.056923 0.147124 -0.3869 0.6988
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ARIMA(1,1,2)
model3.da=Arima(train.diff, order=c(1,1,2),method="ML")
summary(model3.da) #AIC=2040.37## Series: train.diff
## ARIMA(1,1,2)
##
## Coefficients:
## ar1 ma1 ma2
## -0.5997 -0.1969 -0.8029
## s.e. 0.1817 0.1370 0.1361
##
## sigma^2 = 17707925: log likelihood = -1016.19
## AIC=2040.37 AICc=2040.78 BIC=2050.95
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 266.7534 4127.146 3039.496 -2742.671 2960.346 0.6688267 -0.1143309lmtest::coeftest(model3.da) #ada parameter yang tidak signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.59971 0.18171 -3.3005 0.0009653 ***
## ma1 -0.19685 0.13695 -1.4374 0.1506142
## ma2 -0.80291 0.13611 -5.8991 3.655e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ARIMA(2,1,2)
model4.da=Arima(train.diff, order=c(2,1,2),method="ML")
summary(model4.da) #AIC=2036.73## Series: train.diff
## ARIMA(2,1,2)
##
## Coefficients:
## ar1 ar2 ma1 ma2
## -0.4959 -0.2670 -0.4146 -0.5854
## s.e. 0.1673 0.1012 0.1544 0.1530
##
## sigma^2 = 16845088: log likelihood = -1013.36
## AIC=2036.73 AICc=2037.34 BIC=2049.95
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 288.9185 4005.364 3024.152 -2479.662 2731.728 0.6654504
## ACF1
## Training set 0.003638669lmtest::coeftest(model4.da) #seluruh parameter signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.49586 0.16732 -2.9635 0.0030415 **
## ar2 -0.26703 0.10118 -2.6392 0.0083101 **
## ma1 -0.41462 0.15442 -2.6850 0.0072531 **
## ma2 -0.58538 0.15302 -3.8255 0.0001305 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Berdasarkan pendugaan parameter di atas, nilai AIC terkecil dimiliki oleh model ARIMA(2,1,2) dan seluruh parameternya signifikan. Model ARIMA(1,1,2) memiliki nilai AIC kedua terkecil namun ada parameter yang tidak signifikan, sedangkan model ARIMA(0,1,2) juga ada parameter yang tidak signifikan sehingga model yang dipilih adalah model ARIMA(2,1,2).
Analisis Sisaan
Model terbaik hasil identifikasi kemudian dicek asumsi sisaannya. Sisaan model ARIMA harus memenuhi asumsi normalitas, kebebasan sisaan, dan kehomogenan ragam. Diagnostik model dilakukan secara eksplorasi dan uji formal.
Eksplorasi Sisaan
#Eksplorasi
sisaan.da <- model4.da$residuals
par(mfrow=c(2,2))
qqnorm(sisaan.da)
qqline(sisaan.da, col = "blue", lwd = 2)
plot(c(1:length(sisaan.da)),sisaan.da)
acf(sisaan.da)
pacf(sisaan.da) 
par(mfrow = c(1,1))Berdasarkan plot kuantil-kuantil normal, secara eksplorasi ditunjukkan sisaan menyebar normal ditandai dengan titik titik yang cenderung mengikuti garis \(45^{\circ}\). Kemudian dapat dilihat juga lebar pita sisaan yang cenderung sama menandakan bahwa sisaan memiliki ragam yang homogen. Plot ACF dan PACF sisaan ARIMA(2,1,2) juga tidak signifikan pada 20 lag awal yang menandakan saling bebas. Kondisi ini akan diuji lebih lanjut dengan uji formal.
Uji Formal
#1) Sisaan Menyebar Normal
tseries::jarque.bera.test(sisaan.da) #tak tolak H0 < sisaan tidak menyebar normal##
## Jarque Bera Test
##
## data: sisaan.da
## X-squared = 6.6939, df = 2, p-value = 0.03519Selain dengan eksplorasi, asumsi tersebut dapat diuji menggunakan uji formal. Pada tahapan ini uji formal yang digunakan untuk normalitas adalah uji Jarque Bera . Hipotesis pada uji Jarque Bera adalah sebagai berikut.
\(H_0\) : Sisaan menyebar normal
\(H_1\) : Sisaan tidak menyebar normal
Berdasarkan uji KS tersebut, didapat p-value sebesar 0.03519 yang lebih dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa sisaan tidak menyebar normal. Hal ini tidak sesuai dengan hasil eksplorasi menggunakan plot kuantil-kuantil normal.
#2) Sisaan saling bebas/tidak ada autokorelasi
Box.test(sisaan.da, type = "Ljung") #tak tolak H0 > sisaan saling bebas##
## Box-Ljung test
##
## data: sisaan.da
## X-squared = 0.0014303, df = 1, p-value = 0.9698Selanjutnya akan dilakukan uji formal untuk kebebasan sisaan menggunakan uji Ljung-Box. Hipotesis yang digunakan adalah sebagai berikut.
\(H_0\) : Sisaan saling bebas
\(H_1\) : Sisaan tidak tidak saling bebas
Berdasarkan uji Ljung-Box tersebut, didapat p-value sebesar 0.9698 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa sisaan saling bebas. Hal ini sesuai dengan hasil eksplorasi.
#3) Sisaan homogen
Box.test((sisaan.da)^2, type = "Ljung") #tak tolak H0 < sisaan homogen##
## Box-Ljung test
##
## data: (sisaan.da)^2
## X-squared = 4.309, df = 1, p-value = 0.03791Hipotesis yang digunakan untuk uji kehomogenan ragam adalah sebagai berikut.
\(H_0\) : Ragam sisaan homogen
\(H_1\) : Ragam sisaan tidak homogen
Berdasarkan uji Ljung-Box terhadap sisaan kuadrat tersebut, didapat p-value sebesar 0.03791 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa ragam sisaan tidak homogen.
#4) Nilai tengah sisaan sama dengan nol
t.test(sisaan.da, mu = 0, conf.level = 0.95) #tak tolak h0 > nilai tengah sisaan sama dengan 0##
## One Sample t-test
##
## data: sisaan.da
## t = 0.73754, df = 104, p-value = 0.4625
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -487.907 1065.744
## sample estimates:
## mean of x
## 288.9185Terakhir, dengan uji-t, akan dicek apakah nilai tengah sisaan sama dengan nol. Hipotesis yang diujikan sebagai berikut.
\(H_0\) : nilai tengah sisaan sama dengan 0
\(H_1\) : nilai tengah sisaan tidak sama dengan 0
Berdasarkan uji-ttersebut, didapat p-value sebesar 0.4625 yang lebih besar dari taraf nyata 5% sehingga tak tolak \(H_0\) dan menandakan bahwa nilai tengah sisaan sama dengan nol. Hal ini sama dengan hasil eksplorasi.
Overfitting
Tahapan selanjutnya adalah overfitting dilakukan dengan menaikkan orde AR(p) dan MA(q) dari model ARIMA(2,1,2) untuk melihat apakah terdapat model lain yang lebih baik dari model saat ini. Kandidat model overfitting adalah ARIMA(3,1,2) dan ARIMA(2,1,3).
Model ARIMA(3,1,2)
model.overfit1=Arima(train.diff, order=c(3,1,2),method="ML")
summary(model.overfit1) #2038.03 ## Series: train.diff
## ARIMA(3,1,2)
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2
## -0.3789 -0.2394 0.1063 -0.5142 -0.4858
## s.e. 0.2002 0.1032 0.1194 0.1802 0.1789
##
## sigma^2 = 16922761: log likelihood = -1013.01
## AIC=2038.03 AICc=2038.9 BIC=2053.9
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 274.5238 3994.464 2975.148 -2483.256 2729.907 0.6546673
## ACF1
## Training set 0.008051894lmtest::coeftest(model.overfit1) #ar1 dan ar3 tidak signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.37893 0.20024 -1.8924 0.058442 .
## ar2 -0.23940 0.10317 -2.3203 0.020323 *
## ar3 0.10632 0.11938 0.8906 0.373127
## ma1 -0.51423 0.18023 -2.8532 0.004328 **
## ma2 -0.48577 0.17887 -2.7158 0.006612 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Model ARIMA(2,1,3)
model.overfit2=Arima(train.diff, order=c(2,1,3),method="ML")
summary(model.overfit2) #2019.41## Series: train.diff
## ARIMA(2,1,3)
##
## Coefficients:
## ar1 ar2 ma1 ma2 ma3
## -1.0049 -0.9804 0.0728 -0.1813 -0.8915
## s.e. 0.0329 0.0237 0.0905 0.0948 0.0557
##
## sigma^2 = 13806197: log likelihood = -1003.7
## AIC=2019.41 AICc=2020.27 BIC=2035.27
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 247.5603 3607.946 2639.201 -1893.645 2143.475 0.5807438 0.1594988lmtest::coeftest(model.overfit2) #ma1 dan ma2 tidak signifikan##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -1.004934 0.032881 -30.5630 < 2e-16 ***
## ar2 -0.980361 0.023729 -41.3154 < 2e-16 ***
## ma1 0.072771 0.090465 0.8044 0.42116
## ma2 -0.181302 0.094845 -1.9116 0.05593 .
## ma3 -0.891456 0.055734 -15.9948 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#model yang dipilih adalah model awal, yaitu ARIMA(2,1,2)Peramalan
Peramalan dilakukan menggunakan fungsi forecast() .
Contoh peramalan berikut ini dilakukan untuk 30 periode ke depan.
#---FORECAST---#
ramalan.da <- forecast::forecast(model4.da, h = 26)
ramalan.da## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 107 -1179.220086 -6464.152 4105.712 -9261.826 6903.386
## 108 -825.791158 -6136.321 4484.739 -8947.546 7295.964
## 109 1001.739974 -4541.813 6545.293 -7476.393 9479.873
## 110 1.163169 -5589.946 5592.272 -8549.700 8552.026
## 111 9.304452 -5583.446 5602.055 -8544.070 8562.679
## 112 272.451382 -5323.203 5868.105 -8285.363 8830.266
## 113 139.793215 -5457.278 5736.864 -8420.188 8699.775
## 114 135.305092 -5461.938 5732.548 -8424.939 8695.550
## 115 172.954265 -5424.120 5770.029 -8387.032 8732.941
## 116 155.483985 -5441.694 5752.662 -8404.661 8715.629
## 117 154.093356 -5443.110 5751.297 -8406.090 8714.277
## 118 159.448001 -5437.720 5756.616 -8400.682 8719.578
## 119 157.164183 -5440.016 5754.344 -8402.984 8717.313
## 120 156.866788 -5440.317 5754.051 -8403.288 8717.021
## 121 157.624102 -5439.555 5754.803 -8402.523 8717.771
## 122 157.327993 -5439.853 5754.509 -8402.821 8717.477
## 123 157.272597 -5439.909 5754.454 -8402.877 8717.423
## 124 157.379136 -5439.801 5754.560 -8402.770 8717.528
## 125 157.341100 -5439.840 5754.522 -8402.808 8717.490
## 126 157.331511 -5439.849 5754.512 -8402.818 8717.481
## 127 157.346422 -5439.834 5754.527 -8402.803 8717.496
## 128 157.341589 -5439.839 5754.522 -8402.808 8717.491
## 129 157.340004 -5439.841 5754.521 -8402.809 8717.489
## 130 157.342081 -5439.839 5754.523 -8402.807 8717.491
## 131 157.341474 -5439.839 5754.522 -8402.808 8717.491
## 132 157.341220 -5439.840 5754.522 -8402.808 8717.491data.ramalan.da <- ramalan.da$mean
plot(ramalan.da)
Berdasarkan hasil plot ramalan di atas, dapat dilihat bahwa ramalan ARIMA(0,1,2) cenderung stabil hingga akhir periode. Selanjutnya, dapat dicari nilai akurasi antara hasil ramalan dengan data uji sebagai berikut.
pt_1 <- train.ts[106] #nilai akhir data latih
hasil.forc.Diff <- data.ramalan.da
hasil <- diffinv(hasil.forc.Diff, differences = 1) + pt_1
#has.1 sama hasilnta dengan: cumsum(c(pt_1,hasil.forc.Diff))
ts.plot(train.ts,hasil)
perbandingan.da<-matrix(data=c(head(test.ts, n=26), hasil[-1]),
nrow = 26, ncol = 2)
colnames(perbandingan.da)<-c("Aktual","Hasil Forecast")
perbandingan.da## Aktual Hasil Forecast
## [1,] 26439 29138.78
## [2,] 32959 28312.99
## [3,] 43647 29314.73
## [4,] 45129 29315.89
## [5,] 33745 29325.20
## [6,] 29896 29597.65
## [7,] 22655 29737.44
## [8,] 21210 29872.75
## [9,] 32686 30045.70
## [10,] 27857 30201.18
## [11,] 23826 30355.28
## [12,] 31235 30514.73
## [13,] 30048 30671.89
## [14,] 37872 30828.76
## [15,] 50330 30986.38
## [16,] 47541 31143.71
## [17,] 31185 31300.98
## [18,] 31263 31458.36
## [19,] 32392 31615.70
## [20,] 29729 31773.03
## [21,] 36012 31930.38
## [22,] 33170 32087.72
## [23,] 30961 32245.06
## [24,] 38057 32402.40
## [25,] 37204 32559.74
## [26,] 39130 32717.09accuracy(ts(hasil[-1]), head(test.ts, n=26))## ME RMSE MAE MPE MAPE ACF1 Theil's U
## Test set 2950.942 7587.468 5380.307 4.864096 15.03574 0.5171162 1.003985Dari hasil diperoleh nilai MAPE 15.03574 baik untuk menggambarkan model deret waktu.