UAS Analisis Runtun Waktu
~ Exponential Smoothing and Box Jenkins ~
| Kontak | : \(\downarrow\) |
| dsciencelabs@outlook.com | |
| https://www.instagram.com/dsciencelabs/ | |
| RPubs | https://rpubs.com/dsciencelabs/ |
getwd()## [1] "C:/Users/HP/OneDrive/Documents/kuliah/R/optimasi"
library(forecast)
library(tseries)
bipolar <- read.csv("gangguan_bipolar.csv", sep=";")
print(bipolar)## Minggu Gangguan.bipolar...Indonesia.
## 1 02/01/2022 37
## 2 09/01/2022 30
## 3 16/01/2022 43
## 4 23/01/2022 42
## 5 30/01/2022 80
## 6 06/02/2022 100
## 7 13/02/2022 50
## 8 20/02/2022 38
## 9 27/02/2022 38
## 10 06/03/2022 34
## 11 13/03/2022 43
## 12 20/03/2022 39
## 13 27/03/2022 54
## 14 03/04/2022 33
## 15 10/04/2022 36
## 16 17/04/2022 38
## 17 24/04/2022 28
## 18 01/05/2022 33
## 19 08/05/2022 53
## 20 15/05/2022 66
## 21 22/05/2022 33
## 22 29/05/2022 38
## 23 05/06/2022 37
## 24 12/06/2022 39
## 25 19/06/2022 45
## 26 26/06/2022 100
## 27 03/07/2022 53
## 28 10/07/2022 44
## 29 17/07/2022 38
## 30 24/07/2022 28
## 31 31/07/2022 36
## 32 07/08/2022 40
## 33 14/08/2022 27
## 34 21/08/2022 28
## 35 28/08/2022 25
## 36 04/09/2022 26
## 37 11/09/2022 24
## 38 18/09/2022 27
## 39 25/09/2022 28
## 40 02/10/2022 26
## 41 09/10/2022 69
## 42 16/10/2022 39
## 43 23/10/2022 30
## 44 30/10/2022 31
## 45 06/11/2022 29
## 46 13/11/2022 30
## 47 20/11/2022 29
## 48 27/11/2022 24
## 49 04/12/2022 27
## 50 11/12/2022 26
Exponential Smoothing
Plot Data Aktual
Mengubah data menjadi data time series
bip.ts <- ts(bipolar$Gangguan.bipolar...Indonesia., start = 1, end = 50)
plot(bip.ts, col = "blue", lwd = 1, type = "b", pch = 8, main = "Gangguan Bipolar 2022")
grid()Deskripsi data:
Data terdiri dari 50 observasi
yang berarti data tersebut diukur dalam rentang waktu 4 minggu dalam
sebulan yang artinya 52 minggu dalam setahun. Berdasarkan summary data,
terlihat bahwa data berada pada range 30-100. Selanjutnya, berdasarkan
plot deret waktu diats, terlihat bahwa data tidak stasioner.
menggunakan fungsi ses yang terdapat pada package forecast,
pemulusan single exponential dapat dilakukan menggunakan fungsi
HoltWinters yang terdapat di R. Berikut ilustrasinya.
# lambda = 0.5
ses1<- HoltWinters(bip.ts, gamma = FALSE, beta = FALSE, alpha = 0.5)
plot(ses1)
legend("topright",c("data aktual","data pemulusan"), lty=8, col=c("black","red"), cex=0.8)ramalan1<- forecast(ses1, h=3)
#menampilkan ramalan 3 minggu
ramalan1## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 51 26.48198 4.4674161 48.49655 -7.186392 60.15036
## 52 26.48198 1.8689489 51.09502 -11.160405 64.12437
## 53 26.48198 -0.4802451 53.44421 -14.753187 67.71715
Menghitung Nilai Keakuratan
Ukuran kebaikan metode pemulusan yang akan digunakan yaitu SSE, MSE, dan MAPE. Semakin kecil nilai SSE, MSE, dan RMSE berarti metode pemulusan semakin baik.
# Nilai akurasiSES dengan lambda = 0.5
# Berdasarkan hasil pemulusan menggunakan fungsi HoltWinters
SSE1 <- ses1$SSE
MSE1 <- ses1$SSE/length(bip.ts)
RMSE1 <- sqrt(MSE1)
akurasi1 <- matrix(c(SSE1,MSE1,RMSE1))
row.names(akurasi1)<- c("SSE", "MSE", "RMSE")
colnames(akurasi1) <- c("Akurasi lambda = 0.5")
akurasi1## Akurasi lambda = 0.5
## SSE 14173.14424
## MSE 283.46288
## RMSE 16.83636
Selanjutnya, akan dilakukan mulusan menggunakan Double Exponential Smoothing.
Pemulusan Double Exponential Smoothing (DES)
Metode pemulusan Double Exponential Smoothing (DES) digunakan untuk data yang memiliki pola trend. Berbeda dengan metode sebelumnya, pemulusan menggunakan metode ini akan menghasilkan permalan tidak konstan untuk periode berikutnya. Berikut penerapan metode Double Exponential Smoothing pada R dengan kombinasi nilai parameter (λ = 0.5, γ = 0.5)
modeldes <- HoltWinters(bip.ts, gamma = FALSE)
modeldes## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = bip.ts, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.6953373
## beta : 0.0551066
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 25.7381907
## b -0.9460252
ramalan.1<- forecast(modeldes, h=3)
#menampilkan ramalan 10 periode
ramalan.1## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 51 24.79217 1.290530 48.29380 -11.15048 60.73482
## 52 23.84614 -5.302044 52.99432 -20.73216 68.42444
## 53 22.90012 -11.433124 57.23335 -29.60804 75.40827
Menghitung nilai keakuratan
# Akurasi DES lambda = 0.5 dan gamma = 0.5
SSE1= modeldes$SSE
MSE1= SSE1/length(bip.ts)
RMSE1 = sqrt(MSE1)
akurasi.des1 <- matrix(c(SSE1,MSE1,RMSE1))
row.names(akurasi.des1)<- c("SSE", "MSE", "RMSE")
colnames(akurasi.des1) <- c("Akurasi alpha: 0.6953373, beta : 0.0551066, gamma: FALSE")
akurasi.des1## Akurasi alpha: 0.6953373, beta : 0.0551066, gamma: FALSE
## SSE 16326.04203
## MSE 326.52084
## RMSE 18.06989
accuracy(modeldes$fitted[,1], bip.ts)## ME RMSE MAE MPE MAPE ACF1 Theil's U
## Test set 3.291548 18.4425 12.23246 1.56318 25.41557 0.03819848 1.089658
Mean Absolute Percentage Error adalah satuan yang digunakan untuk mengetahui persentase error hasil peramalan. Satuan MAPE lebih sering digunakan untuk reporting karena konsep yang mudah dipahami. Dari output diatas menjelaskan bahwa MAPE nya sebesar 25.41% berarti peramalan dikatakan layak.
Menentukan Metode Pemulusan Terbaik
akurasi1## Akurasi lambda = 0.5
## SSE 14173.14424
## MSE 283.46288
## RMSE 16.83636
akurasi.des1## Akurasi alpha: 0.6953373, beta : 0.0551066, gamma: FALSE
## SSE 16326.04203
## MSE 326.52084
## RMSE 18.06989
Berdasarkan hasil diatas, didapatkan model terbaik antara masing-masing metode. Berdasarkan hasil diatas terlihat bahwa nilai SSE, MSE, dan RMSE metode pemulusan SES dengan λ = 0.5 lebih kecil dibandingkan metode pemulusan DES dengan λ = 0.5 dan γ= 0.5. Hal tersebut menunjukkan bahwa metode SES dengan λ = 0.5 lebih baik digunakan.
ARIMA
getwd()## [1] "C:/Users/HP/OneDrive/Documents/kuliah/R/optimasi"
library(forecast)
library(tseries)
bipolar <- read.csv("gangguan_bipolar.csv", sep=";")
print(bipolar)## Minggu Gangguan.bipolar...Indonesia.
## 1 02/01/2022 37
## 2 09/01/2022 30
## 3 16/01/2022 43
## 4 23/01/2022 42
## 5 30/01/2022 80
## 6 06/02/2022 100
## 7 13/02/2022 50
## 8 20/02/2022 38
## 9 27/02/2022 38
## 10 06/03/2022 34
## 11 13/03/2022 43
## 12 20/03/2022 39
## 13 27/03/2022 54
## 14 03/04/2022 33
## 15 10/04/2022 36
## 16 17/04/2022 38
## 17 24/04/2022 28
## 18 01/05/2022 33
## 19 08/05/2022 53
## 20 15/05/2022 66
## 21 22/05/2022 33
## 22 29/05/2022 38
## 23 05/06/2022 37
## 24 12/06/2022 39
## 25 19/06/2022 45
## 26 26/06/2022 100
## 27 03/07/2022 53
## 28 10/07/2022 44
## 29 17/07/2022 38
## 30 24/07/2022 28
## 31 31/07/2022 36
## 32 07/08/2022 40
## 33 14/08/2022 27
## 34 21/08/2022 28
## 35 28/08/2022 25
## 36 04/09/2022 26
## 37 11/09/2022 24
## 38 18/09/2022 27
## 39 25/09/2022 28
## 40 02/10/2022 26
## 41 09/10/2022 69
## 42 16/10/2022 39
## 43 23/10/2022 30
## 44 30/10/2022 31
## 45 06/11/2022 29
## 46 13/11/2022 30
## 47 20/11/2022 29
## 48 27/11/2022 24
## 49 04/12/2022 27
## 50 11/12/2022 26
Identifikasi Model
bip.ts <- ts(bipolar$Gangguan.bipolar...Indonesia., start = 1, end = 50)
plot(bip.ts, col = "blue", lwd = 1, type = "b", pch = 8, main = "Gangguan Bipolar 2022")
grid()Data diatas tidak stasioner. Uji stasioner terhadap mean menggunakan Augmented Dickey Fuller(ADF) Test.
library(tseries)
adf.test(bip.ts)##
## Augmented Dickey-Fuller Test
##
## data: bip.ts
## Dickey-Fuller = -4.4774, Lag order = 3, p-value = 0.01
## alternative hypothesis: stationary
berdasarkan hasil uji ADF, diperoleh p-value sebesar 0.01. dalam penggunaan alpha = 0.05 maka dapat disimpulkan bahwa Hipotesis nol ditolak dan deret waktu telah stasioner. Pada kasus ini, diasumsikan bahwa data telat stasioner dalam varians. tahap selanjutnya adalah menampilkan ACF dan PACF dari deret waktu.
acf(bip.ts)pacf(bip.ts)estimasi parameter
arimafit <- auto.arima(bip.ts)
summary(arimafit)## Series: bip.ts
## ARIMA(0,1,0)
##
## sigma^2 = 325.5: log likelihood = -211.27
## AIC=424.54 AICc=424.62 BIC=426.43
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.21926 17.86001 11.10074 -6.129483 24.27557 0.9800653 -0.2235273
dugaan1 = arima(bip.ts, order = c(1,0,2));dugaan1##
## Call:
## arima(x = bip.ts, order = c(1, 0, 2))
##
## Coefficients:
## ar1 ma1 ma2 intercept
## 0.2864 0.1817 0.0004 39.6546
## s.e. 0.6620 0.6623 0.3151 3.5140
##
## sigma^2 estimated as 229.5: log likelihood = -206.95, aic = 423.91
summary(dugaan1)##
## Call:
## arima(x = bip.ts, order = c(1, 0, 2))
##
## Coefficients:
## ar1 ma1 ma2 intercept
## 0.2864 0.1817 0.0004 39.6546
## s.e. 0.6620 0.6623 0.3151 3.5140
##
## sigma^2 estimated as 229.5: log likelihood = -206.95, aic = 423.91
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0103458 15.14883 10.14072 -9.394886 23.67397 0.8953067
## ACF1
## Training set 0.009195142
dugaan2 = arima(bip.ts, order = c(1,0,1));dugaan2##
## Call:
## arima(x = bip.ts, order = c(1, 0, 1))
##
## Coefficients:
## ar1 ma1 intercept
## 0.2871 0.181 39.6545
## s.e. 0.3098 0.321 3.5135
##
## sigma^2 estimated as 229.5: log likelihood = -206.95, aic = 421.91
summary(dugaan2)##
## Call:
## arima(x = bip.ts, order = c(1, 0, 1))
##
## Coefficients:
## ar1 ma1 intercept
## 0.2871 0.181 39.6545
## s.e. 0.3098 0.321 3.5135
##
## sigma^2 estimated as 229.5: log likelihood = -206.95, aic = 421.91
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0103695 15.14883 10.14101 -9.394671 23.67444 0.8953327
## ACF1
## Training set 0.00919068
Diagnosa model
bip.jung <- dugaan1$residuals
bip.jungg <- dugaan2$residuals
Box.test(bip.jung, lag = 9)##
## Box-Pierce test
##
## data: bip.jung
## X-squared = 3.1139, df = 9, p-value = 0.9596
Box.test(bip.jungg, lag = 5)##
## Box-Pierce test
##
## data: bip.jungg
## X-squared = 0.63823, df = 5, p-value = 0.9862
Berdasarkan hasil uji Ljung-Box, diperoleh nilai p-value sebesar 0.9596 untuk residual dugaan1 dengan lag sebanyak 9 dan 0.9862 untuk residual dugaan2 dengan lag sebanyak 5. Hipotesis Nol yang digunakan adalah model layak digunakan (uncorrelated residual). Jika menggunakan nilai alpha 5% maka dapat disimpulkan bahwa Hipotesis Nol gagal ditolak dan kesimpulannya adalah model telah layak digunakan. Berdasarkan hasil tersebut, dapat diketahui bahwa ternyata model dugaan1 maupun dugaan2 telah sesuai dengan masuk ke kandidat model terbaik. Untuk memilih salah satu model terbaik, digunakan perbandingan akurasi prediksi berdarkan nilai Mean Absolute Percentage Error (MAPE) dengan perintah sebagai berikut:
library(forecast)
accuracy(dugaan1)## ME RMSE MAE MPE MAPE MASE
## Training set 0.0103458 15.14883 10.14072 -9.394886 23.67397 0.8953067
## ACF1
## Training set 0.009195142
accuracy(dugaan2)## ME RMSE MAE MPE MAPE MASE
## Training set 0.0103695 15.14883 10.14101 -9.394671 23.67444 0.8953327
## ACF1
## Training set 0.00919068
Berdasarkan nilai MAPE yang diperoleh, ternyata model dugaan1 memberikan nilai MAPE yang lebih kecil. Maka dapat disimpulkan bahwa model terbaik untuk studi kasus ini adalah model ARMA(1,2).
Forecasting
Tahap terakhir adalah melakukan prediksi untuk jangka waktu ke depan. Jika ingin diketahui nilai prediksi untuk 3 minggu ke depan, maka perintah yang digunakan adalah:
library(forecast)
prediksi = forecast(dugaan1, h = 3);prediksi## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 51 34.12140 14.70739 53.53540 4.430232 63.81256
## 52 38.06655 16.63066 59.50244 5.283185 70.84991
## 53 39.19983 17.60560 60.79405 6.174313 72.22534