ARIMA
Raid Naufal
2024-11-05
Tahap Persiapan
library("forecast")## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary("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':
##
## tarlibrary("graphics")
library("tseries")
library("readxl")dataawal<-read_xlsx("D:/TUGAS/R/datalatihan.xlsx")
data1<-dataawal[-111:-120,]
data2<-dataawal[-1:-110,]persiapan
data.ts<-ts(data1)
head(data.ts)## Time Series:
## Start = 1
## End = 6
## Frequency = 1
## Sales (in thousands)
## [1,] 10618.1
## [2,] 10537.9
## [3,] 10209.3
## [4,] 10553.0
## [5,] 9934.9
## [6,] 10534.5ts.plot(data.ts, xlab="Time Period ", ylab="Sales", main= "Time Series Plot Data Sales")
points(data.ts)
persiapan sudah selesai, langkah selanjutnya adalah mengidentifikasi
model yang sesuai
Model ARIMA
adf.test(data.ts)## Warning in adf.test(data.ts): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: data.ts
## Dickey-Fuller = -5.9291, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationaryacf(data.ts)
pacf(data.ts, lag.max = 50)
eacf(data.ts)## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o x o o o o o o o o o o o
## 1 o o o o o o o o o o o o o o
## 2 o 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 o o x o o o o o o o o o o o
## 6 x o o o o o o o o o o o o o
## 7 x x o o x o o o o o o o o oPenentuan Model Terbaik
arima(data.ts, order=c(0,0,3), method="ML")##
## Call:
## arima(x = data.ts, order = c(0, 0, 3), method = "ML")
##
## Coefficients:
## ma1 ma2 ma3 intercept
## 0.0141 -0.0518 -0.2326 10374.4514
## s.e. 0.1002 0.1171 0.0967 14.9337
##
## sigma^2 estimated as 45079: log likelihood = -745.56, aic = 1499.12arima(data.ts, order=c(3,0,0), method="ML") ##
## Call:
## arima(x = data.ts, order = c(3, 0, 0), method = "ML")
##
## Coefficients:
## ar1 ar2 ar3 intercept
## 0.0904 0.0086 -0.2312 10374.9540
## s.e. 0.0927 0.0929 0.0926 17.9023
##
## sigma^2 estimated as 44745: log likelihood = -745.15, aic = 1498.3arima(data.ts, order=c(1,0,1), method="ML")##
## Call:
## arima(x = data.ts, order = c(1, 0, 1), method = "ML")
##
## Coefficients:
## ar1 ma1 intercept
## 0.0672 0.0246 10375.1957
## s.e. 0.4398 0.4314 22.7734
##
## sigma^2 estimated as 47350: log likelihood = -748.18, aic = 1502.36arima(data.ts, order=c(3,0,1), method="ML") ##
## Call:
## arima(x = data.ts, order = c(3, 0, 1), method = "ML")
##
## Coefficients:
## ar1 ar2 ar3 ma1 intercept
## 0.6498 -0.0437 -0.2197 -0.6094 10374.0433
## s.e. 0.2868 0.1149 0.1105 0.2950 12.8663
##
## sigma^2 estimated as 43484: log likelihood = -743.64, aic = 1497.28ARIMA_data <- Arima(data.ts, order=c(3,0,1), method="ML") Peramalan 10 Waktu Kedepan
hasil_forecastt <- forecast(ARIMA_data)
hasil_forecast <- as.data.frame(hasil_forecastt)
hasil_forecast## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 111 10371.02 10097.49 10644.55 9952.692 10789.35
## 112 10361.97 10088.22 10635.72 9943.305 10780.64
## 113 10362.54 10088.75 10636.33 9943.808 10781.27
## 114 10367.76 10086.66 10648.86 9937.849 10797.67
## 115 10373.12 10088.65 10657.58 9938.063 10808.17
## 116 10376.24 10090.72 10661.76 9939.579 10812.91
## 117 10376.89 10091.37 10662.41 9940.230 10813.56
## 118 10376.00 10090.29 10661.72 9939.037 10812.97
## 119 10374.71 10088.73 10660.69 9937.338 10812.08
## 120 10373.76 10087.68 10659.84 9936.242 10811.29plot(hasil_forecastt)
garis biru mewakili prediksi titik, area berbayang biru tua menunjukkan
interval prediksi 80% dan area berbayang biru muda menunjukkan interval
prediksi 95% untuk prediksi titik
Penilaian Akurasi dengan MAPE
data_forecast1 <- hasil_forecast$`Point Forecast`
data_forecast1 <- as.data.frame(data_forecast1)
head(data_forecast1)## data_forecast1
## 1 10371.02
## 2 10361.97
## 3 10362.54
## 4 10367.76
## 5 10373.12
## 6 10376.24dataaktual <- as.data.frame(data2)
dataaktual## Sales (in thousands)
## 1 10210.1
## 2 10352.5
## 3 10423.8
## 4 10519.3
## 5 10596.7
## 6 10650.0
## 7 10741.6
## 8 10246.0
## 9 10354.4
## 10 10155.4akurasi.arima <- accuracy(data_forecast1$data_forecast1,dataaktual$`Sales (in thousands)`)
akurasi.arima## ME RMSE MAE MPE MAPE
## Test set 53.57804 193.8149 161.3917 0.4821229 1.538676Simple Exponential Smoothing
plot(data.ts,main="Plot semua data")
points(data.ts)
data1.ses<-HoltWinters(data.ts,alpha=0.1, beta=FALSE, gamma=FALSE)
data1.ses## Holt-Winters exponential smoothing without trend and without seasonal component.
##
## Call:
## HoltWinters(x = data.ts, alpha = 0.1, beta = FALSE, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.1
## beta : FALSE
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 10393.24SSE <- data1.ses$SSE
SSE## [1] 5954375MSE <- SSE/110
MSE## [1] 54130.69hasil_forecast.sess <- forecast(data1.ses)
hasil_forecast.ses <- as.data.frame(hasil_forecast.sess)
hasil_forecast.ses## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 111 10393.24 10093.50 10692.98 9934.830 10851.65
## 112 10393.24 10092.01 10694.48 9932.543 10853.94
## 113 10393.24 10090.52 10695.96 9930.268 10856.22
## 114 10393.24 10089.04 10697.44 9928.004 10858.48
## 115 10393.24 10087.57 10698.92 9925.751 10860.73
## 116 10393.24 10086.10 10700.38 9923.509 10862.97
## 117 10393.24 10084.64 10701.84 9921.278 10865.21
## 118 10393.24 10083.19 10703.29 9919.057 10867.43
## 119 10393.24 10081.74 10704.74 9916.846 10869.64
## 120 10393.24 10080.30 10706.18 9914.645 10871.84plot(hasil_forecast.sess)
prediksi titik, area berbayang biru tua menunjukkan interval prediksi
80% dan area berbayang biru muda menunjukkan interval prediksi 95% untuk
prediksi titik
data_forecast2 <- hasil_forecast.ses$`Point Forecast`
data_forecast2 <- as.data.frame(data_forecast2)
head(data_forecast2)## data_forecast2
## 1 10393.24
## 2 10393.24
## 3 10393.24
## 4 10393.24
## 5 10393.24
## 6 10393.24dataaktual <- as.data.frame(data2)
dataaktual## Sales (in thousands)
## 1 10210.1
## 2 10352.5
## 3 10423.8
## 4 10519.3
## 5 10596.7
## 6 10650.0
## 7 10741.6
## 8 10246.0
## 9 10354.4
## 10 10155.4akurasi.ses <- accuracy(data_forecast2$data_forecast2,dataaktual$`Sales (in thousands)`)
akurasi.ses## ME RMSE MAE MPE MAPE
## Test set 31.73833 189.8569 161.3 0.2723997 1.540696Perbandingan Model ARIMA dan Exponential Smoothing Forecasts.
akurasi.arima## ME RMSE MAE MPE MAPE
## Test set 53.57804 193.8149 161.3917 0.4821229 1.538676akurasi.ses## ME RMSE MAE MPE MAPE
## Test set 31.73833 189.8569 161.3 0.2723997 1.540696fit <- Arima(data.ts, order=c(3,0,1), method="ML")
fit## Series: data.ts
## ARIMA(3,0,1) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ma1 mean
## 0.6498 -0.0437 -0.2197 -0.6094 10374.0433
## s.e. 0.2868 0.1149 0.1105 0.2950 12.8663
##
## sigma^2 = 45555: log likelihood = -743.64
## AIC=1499.28 AICc=1500.1 BIC=1515.49hasil_forecast## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 111 10371.02 10097.49 10644.55 9952.692 10789.35
## 112 10361.97 10088.22 10635.72 9943.305 10780.64
## 113 10362.54 10088.75 10636.33 9943.808 10781.27
## 114 10367.76 10086.66 10648.86 9937.849 10797.67
## 115 10373.12 10088.65 10657.58 9938.063 10808.17
## 116 10376.24 10090.72 10661.76 9939.579 10812.91
## 117 10376.89 10091.37 10662.41 9940.230 10813.56
## 118 10376.00 10090.29 10661.72 9939.037 10812.97
## 119 10374.71 10088.73 10660.69 9937.338 10812.08
## 120 10373.76 10087.68 10659.84 9936.242 10811.29interval.prediction <- hasil_forecast[,-1]
interval.prediction## Lo 80 Hi 80 Lo 95 Hi 95
## 111 10097.49 10644.55 9952.692 10789.35
## 112 10088.22 10635.72 9943.305 10780.64
## 113 10088.75 10636.33 9943.808 10781.27
## 114 10086.66 10648.86 9937.849 10797.67
## 115 10088.65 10657.58 9938.063 10808.17
## 116 10090.72 10661.76 9939.579 10812.91
## 117 10091.37 10662.41 9940.230 10813.56
## 118 10090.29 10661.72 9939.037 10812.97
## 119 10088.73 10660.69 9937.338 10812.08
## 120 10087.68 10659.84 9936.242 10811.29plot(hasil_forecast$`Point Forecast`, type="n", ylim=range(hasil_forecast$`Lo 80`,hasil_forecast$`Hi 80`))
polygon(c(time(hasil_forecast$`Point Forecast`),rev(time(hasil_forecast$`Point Forecast`))), c(hasil_forecast$`Hi 80`,rev(hasil_forecast$`Lo 80`)),
col=rgb(0,0,0.6,0.2), border=FALSE)
lines(hasil_forecast$`Point Forecast`)
lines(fitted(fit),col='green')
out <- (hasil_forecast$`Point Forecast`< hasil_forecast$`Lo 80` | hasil_forecast$`Point Forecast` > hasil_forecast$`Hi 80`)plot(hasil_forecast$`Point Forecast`, type="n", ylim=range(hasil_forecast$`Lo 95`,hasil_forecast$`Hi 95`))
polygon(c(time(hasil_forecast$`Point Forecast`),rev(time(hasil_forecast$`Point Forecast`))), c(hasil_forecast$`Hi 95`,rev(hasil_forecast$`Lo 95`)),
col=rgb(0,0,0.6,0.2), border=FALSE)
lines(hasil_forecast$`Point Forecast`)
lines(fitted(fit),col='yellow')
out <- (hasil_forecast$`Point Forecast`< hasil_forecast$`Lo 95` | hasil_forecast$`Point Forecast` > hasil_forecast$`Hi 95`)Berdasarkan kedua plot yaitu plot interval prediksi 80% dan plot interval prediksi 95% menunjukkan hasil yang hampir sama hanya berbeda di rentang selangnya. Interval prediksi 95% memiliki rentang yang lebih luas dari pada interval prediksi 80%.