Implementasi Model Autoregressive Integrated Moving Average (ARIMA) untuk Peramalan Email Ham
Novi Fitianingrum
10/01/2020
Project ini merupakan Project Kelompok Mata Kuliah Analisis Runtun Waktu
Pendahuluan
Peramalan pada dasarnya merupakan suatu dugaan atau perkiraan atas kejadian di waktu mendatang. Peramalan merupakan alat bantu yang penting dalam perencanaan yang efektif dan efisien (Makridakis, dkk., 1999: 3). Peramalan sering digunakan pada bidang ekonomi, perencanaan produksi, peramalan penjualan, dan kontrol stok. Meramal juga dapat didasarkan pada keahlian penilaian, yang ada pada gilirannya didasarkan pada data historis dan pengalaman (Makridakis, dkk., 1999: 519). Keakuratan dari hasil peramalan dalam pemilihan metode peramalan harus dilakukan dengan teliti agar bisa dipertanggungjawabkan. Salah satu metode peramalan yang sering digunakan adalah runtun waktu ( time series).
Runtun waktu ( time series ) adalah analisis yang mempertimbangkan pengaruh waktu secara beruntun. Data-data yang dikumpulkan berdasarkan urutan waktu seperti, jam, hari, minggu, bulan, kuartal, semester, dan tahun dapat dianalisis menggunakan metode runtun waktu. Data runtun waktu dapat dijadikan dasar dalam pengambilan ketentuan untuk meramalkan kejadian yang terjadi di era yang akan datang.
Model Autoregresif Integrated Moving Average (ARIMA) adalah model yang secara penuh mengabaikan independen variabel dalam membuat peramalan. ARIMA menggunakan nilai masa lalu dan sekarang dari variabel dependen untuk menghasilkan peramalan jangka pendek yang akurat.Model ARIMA merupakan model .Metode ARIMA sering digunakan untuk peramalan di berbagai bidang, salah satunya bidang komunikasi.
Komunikasi adalah salah satu bidang yang perkembangannya sudah pesat. Saat ini sudah ada internet yang akan memudahan kita untuk berkomunikasi secara mudah. Email atau surat elektronik merupakan salah satu penerapan internet dalam bidang komunikasi yang paling digemari karena mudah digunakan, cepat dan berbiaya murah. Penggunaan email yang terbilang intens ini dapat membawa dampak positif dan negatif karena pada kenyataannya tidak semua orang menggunakan email dengan baik dan bahkan ada banyak sekali penyalahgunaan email sehingga berpotensi untuk merugikan orang lain. Email yang disalah gunakan ini biasa kita kenal sebagai spam atau junkmail (email sampah) yang mana email tersebut berisikan iklan, penipuan dan bahkan virus. Selain email spam ada juga email ham yang merupakan lawan dari email spam atau email yang diharapkan datang.
Email ini jumlahnya beragam perharinya. Dikarenakan jumlahnya yang memungkinkan akan meningkat setiap harinya. Maka perlu dilakuan pengendalian terhadap jumlah email ham yang masuk agar email ham tidak menumpuk. Untuk melakukan pengendalian, bisa direncanakan dengan meramalkan migrasi di waktu yang akan dating. Maka akan dilakukan permalan mengenai data runtun waktu banyaknya jumlah email ham yang masuk ke suatu akun dalam periode harian dimulai dari Januari 2017 sampai dengan Juni 2018, dilakukan peramalan ini memiliki tujuan yaitu menentukan model ARIMA terbaik yang mampu menirukan dan meramalkan perilaku data historis, meramalkan jumlah email ham menggunakan model ARIMA terbaik dan mengetahui tingkat kebaikan model ARIMA terbaik yang diperoleh.
Data
Data yang akan diteliti adalah data runtun waktu mengenai jumlah email Ham yang masuk secara harian sebanyak 506 hari mulai dari 13 Januari 2000 sampai dengan 1 Juni 2001. Data akan diteliti untuk dilakukan peramalan jumlah email ham yang masuk menggunakan model ARIMA.
Persiapan
Input packages yang akan digunakan dan Input Data
library(TSA)##
## Attaching package: 'TSA'## The following objects are masked from 'package:stats':
##
## acf, arima## The following object is masked from 'package:utils':
##
## tarlibrary(tseries)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(forecast)## Registered S3 methods overwritten by 'forecast':
## method from
## fitted.Arima TSA
## plot.Arima TSAlibrary(EnvStats)##
## Attaching package: 'EnvStats'## The following objects are masked from 'package:TSA':
##
## kurtosis, skewness## The following objects are masked from 'package:stats':
##
## predict, predict.lm## The following object is masked from 'package:base':
##
## print.defaultHasil Analisis
Plot Email Ham
Data Email Ham memiliki plot time series sebagai berikut:
plot.ts(data, ylab='Count',main='Plot Email Ham')
Uji Stasioner
Selanjutnya, data akan di uji apakah memenuhi data stasioner atau tidak. Untuk mengecek apakah variansi sudah stasioner, dilakukan uji BoxCox. Untuk pengujian BoxCox dalam software R, data yang di inputkan harus bernilai positif dan bukan nol. Karna data yang dimiliki terdapat nilai nol. Maka perlu di tambah konstanta c agar nilai tidak terdapat nol. Dan nilai c yang digunakan adalah satu.
Selanjutnya dilakukan uji stasioner terhadap variansinya.
boxcox<-boxcox(Dataku, optimize = TRUE)
boxcox$lambda## [1] 0.308434Setelah diuji, didapat data tidak stasioner karna memiliki lambda mendekati nol atau tidak sama dengan satu (𝜆 ≠ 1) . Maka dilakukan transformasi untuk membuat data menjadi stasioner.
data1.trans <- BoxCox(data, 0.308434)Setelah data stasioner terhadap varian, data akan diuji apakah data stasioner terhadap mean menggunakan uji KPSS dan dilihat secara visual dengan Plot ACF dan PACF. Hasil Uji KPSS sebagai berikut:
kpss.test(data1.trans)## Warning in kpss.test(data1.trans): p-value smaller than printed p-value##
## KPSS Test for Level Stationarity
##
## data: data1.trans
## KPSS Level = 1.3685, Truncation lag parameter = 5, p-value = 0.01Berdasarkan hasil yang diperoleh, P-Value < 𝛼 = 0,05 menunjukan data tidak stasioner. Dan juga berdasarkan plot ACF dan PACF di bawah,
par(mfrow=c(1,2))
acf(data1.trans)
pacf(data1.trans)
dapat dilihat, masih terdapat seasonal dalam data. Maka perlu dilakukannya differencing agar data stasioner terhadap mean.
Data.diff=diff(data1.trans, d=1)Setelah di differencing, didapatkan hasil Uji KPSS sebagai berikut.
kpss.test(Data.diff)## Warning in kpss.test(Data.diff): p-value greater than printed p-value##
## KPSS Test for Level Stationarity
##
## data: Data.diff
## KPSS Level = 0.0089994, Truncation lag parameter = 5, p-value = 0.1Karna P-Value > 𝛼 = 0.05 maka didapat data sudah stasioner terhadap mean.
Identifikasi Model
Ketika data dikatakan telah stasioner, maka langkah yang dilakukan untuk membentuk model ARIMA adalah dengan identifikasi model. Untuk mengidentifikasi model ARIMA dapat dilihat dari Plot ACF, PACF, EACF, dan Grafik BIC dari data “Data.diff” .
Plot ACF dan PACF.
par(mfrow=c(1,2))
acf(Data.diff)
pacf(Data.diff)
Berdasarkan Plot ACF dan PACF diatas didapat model ARIMA (6,1,2) yang merupakan kandidat model pertama.
Plot EACF
Selanjutnya untuk mengidentifikasi model dapat juga dilihat dari plot EACF berikut
eacf(Data.diff, ar.max = 10, ma.max = 10)## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10
## 0 x x o o x o x o x o o
## 1 x x o o x o x o x o o
## 2 x x o o x o x o x o o
## 3 x o x o x o x o x x x
## 4 x x x x o x x x o o o
## 5 x x x x o x x x o o o
## 6 o x x x o x x o o o o
## 7 x x o x o o x x o o o
## 8 x o o o o o x o o o o
## 9 x o x o o o x x x x o
## 10 x x o o o o x x o x oDari plot EACF didapatkan tiga kandidat model ARIMA lainnya yaitu ARIMA(5,1,8), ARIMA(4,1,8) dan ARIMA(6,1,7).
Grafik BIC
x <- armasubsets(Data.diff, nar = 10, nma = 10)
plot(x)
Berdasarkan grafik BIC diatas, didapatkan kandidat model yaitu model ARIMA (7,1,8).
Dari ketiga proses identifikasi model ARIMA didapat lima kandidat model ARIMA(6,1,2), ARIMA(5,1,8), ARIMA(4,1,8), ARIMA(6,1,7), dan ARIMA(7,1,8)
Setelah identifikasi model, langkah selanjutnya yaitu akan dilakukan estimasi pada kandidat model ARIMA terpilih.
Estimasi Parameter Model ARIMA
Berikut adalah estimasi parameter model ARIMA yang terpilih sebagai kandidat.
model.arima1 <- arima(Data.diff, order = c(6,1,2))
model.arima2 <- arima(Data.diff, order = c(5,1,8))## Warning in log(s2): NaNs produced## Warning in stats::arima(x = x, order = order, seasonal = seasonal, xreg =
## xreg, : possible convergence problem: optim gave code = 1model.arima3 <- arima(Data.diff, order = c(4,1,8))## Warning in stats::arima(x = x, order = order, seasonal = seasonal, xreg =
## xreg, : possible convergence problem: optim gave code = 1model.arima4 <- arima(Data.diff, order = c(6,1,7))## Warning in stats::arima(x = x, order = order, seasonal = seasonal, xreg =
## xreg, : possible convergence problem: optim gave code = 1model.arima5 <- arima(Data.diff, order = c(7,1,8))
model.arima1##
## Call:
## arima(x = Data.diff, order = c(6, 1, 2))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ma1 ma2
## -0.7375 -0.7782 -0.7666 -0.6934 -0.7614 -0.5732 -0.9612 -0.0388
## s.e. 0.0605 0.0407 0.0474 0.0451 0.0394 0.0439 0.0710 0.0707
##
## sigma^2 estimated as 3.183: log likelihood = -1013.84, aic = 2043.67model.arima2##
## Call:
## arima(x = Data.diff, order = c(5, 1, 8))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 ma2 ma3
## -0.1962 -0.6408 -0.6359 -0.1951 -0.9955 -1.5442 0.9280 -0.2836
## s.e. 0.0022 0.0035 0.0053 0.0027 0.0046 0.0453 0.0747 0.0644
## ma4 ma5 ma6 ma7 ma8
## -0.4755 1.1636 -1.5201 0.6641 0.0677
## s.e. 0.0340 0.0332 0.0613 0.0740 0.0472
##
## sigma^2 estimated as 2.669: log likelihood = -974.63, aic = 1975.27model.arima3##
## Call:
## arima(x = Data.diff, order = c(4, 1, 8))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ma1 ma2 ma3 ma4
## -0.2303 -0.0761 -0.3260 -0.9187 -1.5221 0.1267 0.7078 0.4968
## s.e. 0.0218 0.0185 0.0169 0.0195 0.0491 0.0932 0.0841 0.0434
## ma5 ma6 ma7 ma8
## -1.6292 0.6268 0.4319 -0.2387
## s.e. 0.0455 0.0837 0.0880 0.0484
##
## sigma^2 estimated as 3.388: log likelihood = -1032.55, aic = 2089.1model.arima4##
## Call:
## arima(x = Data.diff, order = c(6, 1, 7))
##
## Coefficients:## Warning in sqrt(diag(x$var.coef)): NaNs produced## ar1 ar2 ar3 ar4 ar5 ar6 ma1 ma2
## -0.0841 -0.6230 -0.5616 -0.1302 -0.9711 0.1074 -1.6491 1.0937
## s.e. 0.0602 0.0071 0.0313 0.0403 NaN 0.0643 0.0419 0.0245
## ma3 ma4 ma5 ma6 ma7
## -0.3901 -0.4322 1.1956 -1.6318 0.814
## s.e. NaN NaN NaN 0.0689 0.055
##
## sigma^2 estimated as 2.67: log likelihood = -974.57, aic = 1975.14model.arima5##
## Call:
## arima(x = Data.diff, order = c(7, 1, 8))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ma1
## -0.7570 -0.7709 -0.7603 -0.7637 -0.7725 -0.7596 0.2225 -0.9450
## s.e. 0.0806 0.0652 0.0708 0.0752 0.0620 0.0824 0.0601 0.0747
## ma2 ma3 ma4 ma5 ma6 ma7 ma8
## 0.0737 -0.0407 0.0132 0.0322 -0.0573 -0.9236 0.8496
## s.e. 0.0961 0.0566 0.0373 0.0735 0.1008 0.1169 0.0635
##
## sigma^2 estimated as 2.54: log likelihood = -965.73, aic = 1961.47Untuk mengetahui parameter signifikan atau tidak, nilai dari koefisien parameter dibagi dengan nilai standar error nya. Jika hasilnya kurang dari -1.96 atau lebih besar dari 1.96 maka parameter dapat dikatakan signifikan.
Berdasarkan output diatas, tidak di dapatkan model yang seluruh parameternya signifikan. Maka untuk mencari model terbaiknya, dapat dilihat dari nilai AIC terendah. Di dapat nilai AIC terendah terdapat pada Model ARIMA (7,1,8). Sehingga model ARIMA (7,1,8) lah yang akan dilanjutkan ke tahap diagnosis model.
fitARIMA <- stats::arima(Data.diff, order=c(7,1,8))Diagnosis Model
Pada tahap ini model terpilih yaitu ARIMA (7,1,8) akan diuji kenormalan, autokorelasi, dan heteroskedastisitas dari residualnya. Uji ini bertujuan untuk melihat apakah data “fitARIMA” banyak email ham yang masuk telah memenuhi asumsi tersebut sebelum diprediksi.
Cek Heteroskedastisitas
Berikut adalah plot dari residual model ARIMA (7,1,8).
plot(rstandard(fitARIMA)) 
Berdasarkan residual plot residual diatas didapat residual cenderung stabil maka dikatakan terdapat homokedastisitas.
Uji Normalitas
Untuk menguji kenormalan dapat menggunakan uji Jarque Bera
jarque.bera.test(rstandard(fitARIMA))##
## Jarque Bera Test
##
## data: rstandard(fitARIMA)
## X-squared = 5.3288, df = 2, p-value = 0.06964Karena P-value > 0.05, maka residual berdistribusi normal.
Uji Autokorelasi
Uji autokorelasi dapat dilakukan dengan Ljung-Box
Box.test(rstandard(fitARIMA), type=c("Ljung-Box"))##
## Box-Ljung test
##
## data: rstandard(fitARIMA)
## X-squared = 0.0010326, df = 1, p-value = 0.9744Karena P-value > 0.05, maka data residual tidak terdapat autokorelasi.
Peramalan ( Forecasting )
Pada langkah forecasting, karna data yang dipakai merupakan data hasil transformasi, maka data hasil forecasting harus dikembalikan ke data awal atau di transformasi balik ke data awal dan juga dikurangkan dengan konstanta c=1. Karna akan mencari forcasting untuk data awal email ham count yang di berikan, bukan data hasil transformasi.
forecasting <- forecast(fitARIMA , h = 20)
plot(forecast(fitARIMA , h = 20))
forecasting=data.frame(forecasting)
forecasting## Point.Forecast Lo.80 Hi.80 Lo.95 Hi.95
## 506 -4.9131524937 -6.9621681 -2.864137 -8.0468514 -1.7794536
## 507 -0.4691793047 -2.9704389 2.032080 -4.2945257 3.3561671
## 508 2.8209012577 0.3085731 5.333229 -1.0213732 6.6631757
## 509 0.5655325051 -1.9485415 3.079606 -3.2794118 4.4104769
## 510 0.0005427919 -2.5135798 2.514665 -3.8444760 3.8455616
## 511 1.0721176307 -1.4423309 3.586566 -2.7733996 4.9176348
## 512 -0.3244549082 -2.8395816 2.190672 -4.1710092 3.5220994
## 513 -3.9771100098 -6.4941416 -1.460078 -7.8265777 -0.1276423
## 514 -0.2868924149 -2.8064216 2.232637 -4.1401797 3.5663949
## 515 2.8604838623 0.3403078 5.380660 -0.9937928 6.7147605
## 516 0.5761976632 -1.9440492 3.096445 -3.2781873 4.4305826
## 517 0.0021931467 -2.5180534 2.522440 -3.8521913 3.8565776
## 518 1.1077778518 -1.4128302 3.628386 -2.7471595 4.9627152
## 519 -0.3403357510 -2.8615252 2.180854 -4.1961623 3.5154908
## 520 -3.9627833984 -6.4857048 -1.439862 -7.8212587 -0.1043081
## 521 -0.3116199444 -2.8368096 2.213570 -4.1735642 3.5503243
## 522 2.8524291124 0.3266903 5.378168 -1.0103551 6.7152133
## 523 0.5761642659 -1.9496191 3.101948 -3.2866880 4.4390165
## 524 0.0018330922 -2.5239507 2.527617 -3.8610199 3.8646861
## 525 1.1420127756 -1.3841077 3.668133 -2.7213550 5.0053806forecast.asli=(exp(forecasting$Point.Forecast))-1
forecast.asli## [1] -0.9926507168 -0.3744845848 15.7919777545 0.7603849471 0.0005429392
## [6] 1.9215597386 -0.2770786953 -0.9812602812 -0.2494075230 16.4699779555
## [11] 0.7792602059 0.0021955534 2.0276230866 -0.2884686147 -0.9809898722
## [16] -0.2677402252 16.3298268506 0.7792007843 0.0018347733 2.1330681862Setelah hasil forecastingnya di transformasi, didapat hasil forecasting sudah cukup baik untuk meramalkan data Ham yang masuk dalam 20 langkah kedepan atau 20 hari kedepan.
Kesimpulan
Berdasarkan hasil penelitian dan pembahasan, diperoleh kesimpulan bahwa model terbaik dari metode ARIMA ( Autoregressive Integrated Moving Average ) yang digunakan untuk melakukan peramalan/prediksi jumlah email Ham yang masuk secara harian sebanyak 506 hari mulai dari 13 Januari 2000 sampai dengan 1 Juni 2001 yaitu ARIMA (7,1,8).
Saran
Perlu dilakukan persiapan yang cermat terhadap data pengamatan yang akan digunakan dalam peramalan karena kesalahan dalam mengidentifikasi karakteristik/ pola data yang diinput akan menyebabkan kekeliruan saat menentukan model peramalan. Data historis yang digunakan harus bersifat representatif dan dalam jumlah yang cukup agar diperoleh model yang layak untuk peramalan.
Model peramalan ini hanya efektif digunakan untuk peramalan jangka pendek dan kurang efektif untuk peramalan jangka panjang. Jika ingin meramalkan untuk waktu yang berbeda, sebaiknya model yang telah ada (ARIMA (7,1,8)) diperbaharui sesuai dengan data baru yang diperoleh.