Eksplorasi dan Pemulusan Data Time Series
Kelompok 5
P2-MPDW 2022
Anggota:
Hanung Safrizal (G1401201050)
Maulana Ahsan Fadillah (G1401201062)
Muhammad Irsyad Robbani (G1401201064)
Muhammad Raziv Zulfikar (G1401201083)
Dhea Puspita Adinda (G1401201090)
Pendahuluan
Latar Belakang
Beras merupakan makanan pokok bagi sebagian besar masyarakat Indonesia. Konsumsi beras di Indonesia meningkat setiap tahunnya seiring dengan meningkatnya jumlah penduduk di Indonesia. Ketergantungan masyarakat Indonesia terhadap beras bisa menimbulkan masalah jika ketersediaan beras sudah tidak dapat tercukupi. Pemenuhan kebutuhan dan stabilitas harga beras merupakan isu yang tetap relevan dari waktu ke waktu. Jumlah penduduk yang besar menjadikan permasalahan beras bukan hanya pada jumlah ketersediaan dan harga saja, namun termasuk keberagaman kualitas dan jenis produk beras (BPS 2009). Pada awal tahun 2018 harga beras mengalami peningkatan. Kenaikan harga beras ini jika terus dibiarkan akan menyebabkan terjadinya inflasi yang berdampak pada melambatnya pertumbuhan ekonomi nasional serta dampak negatif lainnya. Dalam rangka perumusan kebijakan pengendalian inflasi maka data dan informasi terkait proyeksi keadaan pasar sangat dibutuhkan. Oleh karena itu, pemodelan harga beras di Indonesia sangat perlu dilakukan.
Tujuan
- Mengeksplorasi dan forecasting data time series harga beras premium di Indonesia per bulan pada tahun 2013-2022.
Tinjauan Pustaka
Metode Single Moving Average
Metode single moving average adalah metode peramalan yang menggunakan sejumlah data aktual permintaan yang baru untuk membangkitkan nilai ramalan untuk permintaan dimasa yang akan datang.(Naufal dan Andrean, 2017).
Metode Double Moving Average
Metode Double Moving Average adalah merupakan suatu metode peramalan yang dilakukan dengan menghitung nilai rata-rata bergerak sebanyak dua kali, lalu dilanjutkan dengan meramal menggunakan suatu persamaan. Metode DMA digunakan untuk mengatasi galat sistematis yang terjadi apabila rata-rata bergerak data dipakai pada data yang memiliki tren. Metode ini merupakan pengembangan dari metode SMA (Muzaki 2022).
Metode Single Exponential Smoothing
Metode Single Exponential Smoothing adalah metode yang menunjukan pembobotan menurun secara eksponensial terhadap nilai observasi yang lebih tua. Yaitu nilai yang lebih baru diberikan bobot yang relatif lebih besar dibanding nilai observasi yang lebih lama. Metode ini memberikan sebuah pembobotan eksponensial rata-rata bergerak dari semua nilai observasi sebelumnya (Hartono et.al, 2012).
Metode Double Exponential Smoothing
Metode Double Exponential Smoothing merupakan metode yang digunakan untuk meramalkan data yang mengalami trend kenaikan. Jika data yang digunakan semakin banyak dalam perhitungan peramalan, maka persentase error peramalannya akan semakin kecil dan sebaliknya (Hanief dan Purwanto 2017).
Winter Smoothing
Metode pemulusan eksponensial linear dari Winter’s digunakan untuk peramalan jika data memiliki komponen musiman. Metode winter didasarkan pada tiga persamaan pemulusan, yaitu stasioneritas, trend, dan musiman. Holt Winter sendiri memiliki dua metode yang berbeda, bergantung pada sifat musiman itu sendiri, yaitu aditif atau multiplikatif.
SSE
Sum Square Error (SSE) menjumlahkan kuadrat dari kesalahan.
MSE
Mean Squared Error (MSE) merata-ratakan kuadrat dari kesalahan menggunakan penyebut n tanpa memperhatikan derajat bebas model. MSE banyak dipakai dalam optimalisasi pembobotan karena memberikan ketelitian yang lebih baik.
RMSE
Root Mean Square Error (RMSE) diperoleh dengan mengakarkan hasil MSE.
MAD
Simpangan absolut rata-rata atau mean absolute deviation (MAD) merata-ratakan nilai absolut kesalahan peramalan dalam unit ukuran yang sama seperti data aktual.
MAPE
Mean absolute percentage error (MAPE) merata-ratakan persentase absolut kesalahan dan memberikan gambaran seberapa kesalahan peramalan dibanding dengan nilai sebenarnya. MAPE banyak dipakai dalam perbandingan data-data yang mempunyai skala interval waktu yang berbeda.
Package dan Data
1. Memanggil Library
Packages yang digunakan dalam analisis eksplorasi dan pemulusan data deret waktu, sebagai berikut:
library("tseries")## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library("forecast")
library("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':
##
## tar
library("graphics")
library("astsa")##
## Attaching package: 'astsa'
## The following object is masked from 'package:forecast':
##
## gas
library("car")## Loading required package: carData
library("readxl")
library("dplyr")##
## Attaching package: 'dplyr'
## The following object is masked from 'package:car':
##
## recode
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
2. Memanggil Data
Data yang digunakan pada pengerjaan tugas analisis eksplorasi dan pemulusan deret waktua adalah data beras premium di Indonesia per bulan dari tahun 2013 hingga 2022. Data diperoleh dari web BPS Indonesia ( https://www.bps.go.id/indicator/36/500/1/harga-beras-di-penggilingan-menurut-kualitas.html). Pada data harga beras terdapat 3 peubah yang digunakan yaitu beras premium, medium, dan luar kualitas tetapi peubah yang digunakan dalam eksplorasi kali ini adalah data beras premium dengan jumlah amatan sebanyak 115 amatan.
data <- read_xlsx("C:/0 SEM5/MPDW/p1/Data Kelompok 15.xlsx")
data$Premium <- as.numeric(data$Premium)3. Split Data
Spliting data dilakukan untuk memisahkan data menjadi dua bagian yaitu data training dan data testing. Data training digunakan untuk menentukan model yang sesuai dengan keseluruhan data dan membuat model, sedangkan data testing digunakan untuk menguji performa dari model yang telah didapatkan. Jumlah data training sebanyak 88 amatan (Januari 2013 hingga April 2020) dan data testing sebanyak 27 amatan (September 2020 hingga Juli 2022).
train <- data[1:88, 3]
test <- data[89:115, 3]
data.premium <- ts(data[,3])
train.ts <- ts(train)
test.ts <- ts(test)Eksplorasi Data Analisis
Plot Data Time Series
Single Moving Average (SMA)
Pemulusan SMA dengan n=4
data.sma<-SMA(train.ts, n=4)
data.sma## Time Series:
## Start = 1
## End = 88
## Frequency = 1
## [1] NA NA NA 7641.970 7578.912 7522.652 7584.505
## [8] 7669.647 7719.840 7794.297 7818.372 7872.230 7987.885 8102.180
## [15] 8170.402 8156.205 8106.365 8072.450 8081.307 8183.692 8258.962
## [22] 8316.285 8397.995 8570.225 8923.462 9270.422 9524.597 9551.942
## [29] 9298.530 9081.195 8924.512 8924.265 9107.827 9242.740 9397.527
## [36] 9531.725 9601.575 9683.997 9685.892 9551.942 9416.560 9308.970
## [43] 9259.602 9319.485 9301.647 9246.207 9216.955 9210.525 9290.722
## [50] 9359.675 9392.467 9388.222 9389.427 9398.400 9397.187 9425.222
## [57] 9433.822 9448.420 9487.187 9593.100 9812.930 10032.697 10121.150
## [64] 10037.292 9830.785 9604.862 9511.705 9494.982 9507.020 9548.832
## [71] 9611.592 9701.592 9836.355 9927.007 9937.897 9849.550 9687.292
## [78] 9564.192 9490.397 9506.727 9539.812 9575.727 9631.430 9708.495
## [85] 9818.062 9923.470 10008.537 10053.520
data.ramal<-c(NA,data.sma)
data.ramal #forecast 1 periode ke depan## [1] NA NA NA NA 7641.970 7578.912 7522.652
## [8] 7584.505 7669.647 7719.840 7794.297 7818.372 7872.230 7987.885
## [15] 8102.180 8170.402 8156.205 8106.365 8072.450 8081.307 8183.692
## [22] 8258.962 8316.285 8397.995 8570.225 8923.462 9270.422 9524.597
## [29] 9551.942 9298.530 9081.195 8924.512 8924.265 9107.827 9242.740
## [36] 9397.527 9531.725 9601.575 9683.997 9685.892 9551.942 9416.560
## [43] 9308.970 9259.602 9319.485 9301.647 9246.207 9216.955 9210.525
## [50] 9290.722 9359.675 9392.467 9388.222 9389.427 9398.400 9397.187
## [57] 9425.222 9433.822 9448.420 9487.187 9593.100 9812.930 10032.697
## [64] 10121.150 10037.292 9830.785 9604.862 9511.705 9494.982 9507.020
## [71] 9548.832 9611.592 9701.592 9836.355 9927.007 9937.897 9849.550
## [78] 9687.292 9564.192 9490.397 9506.727 9539.812 9575.727 9631.430
## [85] 9708.495 9818.062 9923.470 10008.537 10053.520
data.gab<-cbind(aktual=c(train.ts,rep(NA,5)),pemulusan=c(data.sma,rep(NA,5)),ramalan=c(data.ramal,rep(data.ramal[length(data.ramal)],4)))
data.gab #forecast 5 periode ke depan## aktual pemulusan ramalan
## [1,] 7797.63 NA NA
## [2,] 7773.26 NA NA
## [3,] 7576.27 NA NA
## [4,] 7420.72 7641.970 NA
## [5,] 7545.40 7578.912 7641.970
## [6,] 7548.22 7522.652 7578.912
## [7,] 7823.68 7584.505 7522.652
## [8,] 7761.29 7669.647 7584.505
## [9,] 7746.17 7719.840 7669.647
## [10,] 7846.05 7794.297 7719.840
## [11,] 7919.98 7818.372 7794.297
## [12,] 7976.72 7872.230 7818.372
## [13,] 8208.79 7987.885 7872.230
## [14,] 8303.23 8102.180 7987.885
## [15,] 8192.87 8170.402 8102.180
## [16,] 7919.93 8156.205 8170.402
## [17,] 8009.43 8106.365 8156.205
## [18,] 8167.57 8072.450 8106.365
## [19,] 8228.30 8081.307 8072.450
## [20,] 8329.47 8183.692 8081.307
## [21,] 8310.51 8258.962 8183.692
## [22,] 8396.86 8316.285 8258.962
## [23,] 8555.14 8397.995 8316.285
## [24,] 9018.39 8570.225 8397.995
## [25,] 9723.46 8923.462 8570.225
## [26,] 9784.70 9270.422 8923.462
## [27,] 9571.84 9524.597 9270.422
## [28,] 9127.77 9551.942 9524.597
## [29,] 8709.81 9298.530 9551.942
## [30,] 8915.36 9081.195 9298.530
## [31,] 8945.11 8924.512 9081.195
## [32,] 9126.78 8924.265 8924.512
## [33,] 9444.06 9107.827 8924.265
## [34,] 9455.01 9242.740 9107.827
## [35,] 9564.26 9397.527 9242.740
## [36,] 9663.57 9531.725 9397.527
## [37,] 9723.46 9601.575 9531.725
## [38,] 9784.70 9683.997 9601.575
## [39,] 9571.84 9685.892 9683.997
## [40,] 9127.77 9551.942 9685.892
## [41,] 9181.93 9416.560 9551.942
## [42,] 9354.34 9308.970 9416.560
## [43,] 9374.37 9259.602 9308.970
## [44,] 9367.30 9319.485 9259.602
## [45,] 9110.58 9301.647 9319.485
## [46,] 9132.58 9246.207 9301.647
## [47,] 9257.36 9216.955 9246.207
## [48,] 9341.58 9210.525 9216.955
## [49,] 9431.37 9290.722 9210.525
## [50,] 9408.39 9359.675 9290.722
## [51,] 9388.53 9392.467 9359.675
## [52,] 9324.60 9388.222 9392.467
## [53,] 9436.19 9389.427 9388.222
## [54,] 9444.28 9398.400 9389.427
## [55,] 9383.68 9397.187 9398.400
## [56,] 9436.74 9425.222 9397.187
## [57,] 9470.59 9433.822 9425.222
## [58,] 9502.67 9448.420 9433.822
## [59,] 9538.75 9487.187 9448.420
## [60,] 9860.39 9593.100 9487.187
## [61,] 10349.91 9812.930 9593.100
## [62,] 10381.74 10032.697 9812.930
## [63,] 9892.56 10121.150 10032.697
## [64,] 9524.96 10037.292 10121.150
## [65,] 9523.88 9830.785 10037.292
## [66,] 9478.05 9604.862 9830.785
## [67,] 9519.93 9511.705 9604.862
## [68,] 9458.07 9494.982 9511.705
## [69,] 9572.03 9507.020 9494.982
## [70,] 9645.30 9548.832 9507.020
## [71,] 9770.97 9611.592 9548.832
## [72,] 9818.07 9701.592 9611.592
## [73,] 10111.08 9836.355 9701.592
## [74,] 10007.91 9927.007 9836.355
## [75,] 9814.53 9937.897 9927.007
## [76,] 9464.68 9849.550 9937.897
## [77,] 9462.05 9687.292 9849.550
## [78,] 9515.51 9564.192 9687.292
## [79,] 9519.35 9490.397 9564.192
## [80,] 9530.00 9506.727 9490.397
## [81,] 9594.39 9539.812 9506.727
## [82,] 9659.17 9575.727 9539.812
## [83,] 9742.16 9631.430 9575.727
## [84,] 9838.26 9708.495 9631.430
## [85,] 10032.66 9818.062 9708.495
## [86,] 10080.80 9923.470 9818.062
## [87,] 10082.43 10008.537 9923.470
## [88,] 10018.19 10053.520 10008.537
## [89,] NA NA 10053.520
## [90,] NA NA 10053.520
## [91,] NA NA 10053.520
## [92,] NA NA 10053.520
## [93,] NA NA 10053.520
Hasil peramalan dengan menggunakan N=5 pada metode SMA adalah sebesar 10053.520
Plot time series
ts.plot(data.gab[,1], xlab="Time Period ", ylab="Sales", main= "SMA N=4 Data Premium")
points(data.gab[,1])
lines(data.gab[,2],col="green",lwd=2)
lines(data.gab[,3],col="red",lwd=2)
legend("topleft",c("data aktual","data pemulusan","data peramalan"), lty=8, col=c("black","green","red"), cex=0.8)Menghitung nilai keakuratan
data.premium=as.numeric(data.premium)error.sma = test.ts-data.ramal[1:length(test.ts)]
SSE.sma = sum(error.sma[5:length(test.ts)]^2)
MSE.sma = mean(error.sma[5:length(test.ts)]^2)
MAPE.sma = mean(abs((error.sma[5:length(test.ts)]/test.ts[5:length(test.ts)])*100))
akurasi.sma <- matrix(c(SSE.sma, MSE.sma, MAPE.sma))
row.names(akurasi.sma)<- c("SSE", "MSE", "MAPE")
colnames(akurasi.sma) <- c("Akurasi m = 4")
akurasi.sma## Akurasi m = 4
## SSE 6.172562e+07
## MSE 2.683723e+06
## MAPE 1.616197e+01
Double Moving Average (DMA)
Pemulusan DMA dengan n=4
dma <- SMA(data.sma, n = 4)
At <- 2*data.sma - dma
Bt <- 2/(4-1)*(data.sma - dma)
data.dma<- At+Bt
data.ramal2<- c(NA, data.dma)
t = 1:5
f = c()
for (i in t) {
f[i] = At[length(At)] + Bt[length(Bt)]*(i)
}
data.gab2 <- cbind(aktual = c(train.ts,rep(NA,5)), pemulusan1 = c(data.sma,rep(NA,5)),pemulusan2 = c(data.dma, rep(NA,5)),At = c(At, rep(NA,5)), Bt = c(Bt,rep(NA,5)),ramalan = c(data.ramal2, f[-1]))
data.gab2## aktual pemulusan1 pemulusan2 At Bt ramalan
## [1,] 7797.63 NA NA NA NA NA
## [2,] 7773.26 NA NA NA NA NA
## [3,] 7576.27 NA NA NA NA NA
## [4,] 7420.72 7641.970 NA NA NA NA
## [5,] 7545.40 7578.912 NA NA NA NA
## [6,] 7548.22 7522.652 NA NA NA NA
## [7,] 7823.68 7584.505 7588.663 7587.000 1.663333 NA
## [8,] 7761.29 7669.647 7804.178 7750.366 53.812083 7588.663
## [9,] 7746.17 7719.840 7879.305 7815.519 63.785833 7804.178
## [10,] 7846.05 7794.297 7964.672 7896.522 68.150000 7879.305
## [11,] 7919.98 7818.372 7931.428 7886.206 45.222083 7964.672
## [12,] 7976.72 7872.230 7990.638 7943.275 47.363333 7931.428
## [13,] 8208.79 7987.885 8187.366 8107.574 79.792500 7990.638
## [14,] 8303.23 8102.180 8363.869 8259.193 104.675417 8187.366
## [15,] 8192.87 8170.402 8399.116 8307.631 91.485417 8363.869
## [16,] 7919.93 8156.205 8242.933 8208.242 34.691250 8399.116
## [17,] 8009.43 8106.365 8060.660 8078.942 -18.282083 8242.933
## [18,] 8167.57 8072.450 7982.607 8018.544 -35.937083 8060.660
## [19,] 8228.30 8081.307 8043.350 8058.533 -15.182917 7982.607
## [20,] 8329.47 8183.692 8304.924 8256.431 48.492500 8043.350
## [21,] 8310.51 8258.962 8442.061 8368.822 73.239583 8304.924
## [22,] 8396.86 8316.285 8493.324 8422.508 70.815417 8442.061
## [23,] 8555.14 8397.995 8579.264 8506.756 72.507500 8493.324
## [24,] 9018.39 8570.225 8877.489 8754.583 122.905417 8579.264
## [25,] 9723.46 8923.462 9542.580 9294.933 247.647083 8877.489
## [26,] 9784.70 9270.422 10070.250 9750.319 319.930833 9542.580
## [27,] 9571.84 9524.597 10278.632 9977.018 301.613750 10070.250
## [28,] 9127.77 9551.942 9942.503 9786.279 156.224167 10278.632
## [29,] 8709.81 9298.530 9110.458 9185.687 -75.228750 9942.503
## [30,] 8915.36 9081.195 8609.743 8798.324 -188.580833 9110.458
## [31,] 8945.11 8924.512 8441.958 8634.980 -193.021667 8609.743
## [32,] 9126.78 8924.265 8702.831 8791.404 -88.573750 8441.958
## [33,] 9444.06 9107.827 9271.790 9206.205 65.585000 8702.831
## [34,] 9455.01 9242.740 9564.246 9435.644 128.602500 9271.790
## [35,] 9564.26 9397.527 9779.923 9626.965 152.958333 9564.246
## [36,] 9663.57 9531.725 9884.675 9743.495 141.180000 9779.923
## [37,] 9723.46 9601.575 9865.214 9759.758 105.455417 9884.675
## [38,] 9784.70 9683.997 9901.150 9814.289 86.860833 9865.214
## [39,] 9571.84 9685.892 9786.051 9745.987 40.063333 9901.150
## [40,] 9127.77 9551.942 9420.427 9473.033 -52.606250 9786.051
## [41,] 9181.93 9416.560 9136.496 9248.522 -112.025417 9420.427
## [42,] 9354.34 9308.970 9005.851 9127.099 -121.247500 9136.496
## [43,] 9374.37 9259.602 9051.825 9134.936 -83.110833 9005.851
## [44,] 9367.30 9319.485 9308.369 9312.816 -4.446250 9051.825
## [45,] 9110.58 9301.647 9308.683 9305.869 2.814167 9308.369
## [46,] 9132.58 9246.207 9186.994 9210.679 -23.685417 9308.683
## [47,] 9257.36 9216.955 9126.757 9162.836 -36.079167 9186.994
## [48,] 9341.58 9210.525 9155.010 9177.216 -22.205833 9126.757
## [49,] 9431.37 9290.722 9373.422 9340.342 33.080000 9155.010
## [50,] 9408.39 9359.675 9510.018 9449.881 60.137083 9373.422
## [51,] 9388.53 9392.467 9524.334 9471.587 52.746667 9510.018
## [52,] 9324.60 9388.222 9438.974 9418.673 20.300417 9524.334
## [53,] 9436.19 9389.427 9401.060 9396.407 4.652917 9438.974
## [54,] 9444.28 9398.400 9408.851 9404.671 4.180417 9401.060
## [55,] 9383.68 9397.187 9403.651 9401.066 2.585417 9408.851
## [56,] 9436.74 9425.222 9462.994 9447.886 15.108750 9403.651
## [57,] 9470.59 9433.822 9467.430 9453.987 13.442917 9462.994
## [58,] 9502.67 9448.420 9485.515 9470.677 14.837917 9467.430
## [59,] 9538.75 9487.187 9551.395 9525.712 25.682917 9485.515
## [60,] 9860.39 9593.100 9763.879 9695.567 68.311667 9551.395
## [61,] 10349.91 9812.930 10192.131 10040.451 151.680417 9763.879
## [62,] 10381.74 10032.697 10534.729 10333.916 200.812500 10192.131
## [63,] 9892.56 10121.150 10506.451 10352.331 154.120417 10534.729
## [64,] 9524.96 10037.292 10097.751 10073.567 24.183333 10506.451
## [65,] 9523.88 9830.785 9539.625 9656.089 -116.464167 10097.751
## [66,] 9478.05 9604.862 9115.429 9311.203 -195.773333 9539.625
## [67,] 9519.93 9511.705 9120.945 9277.249 -156.304167 9115.429
## [68,] 9458.07 9494.982 9302.314 9379.381 -77.067500 9120.945
## [69,] 9572.03 9507.020 9469.316 9484.398 -15.081667 9302.314
## [70,] 9645.30 9548.832 9604.162 9582.030 22.131667 9469.316
## [71,] 9770.97 9611.592 9729.902 9682.578 47.323750 9604.162
## [72,] 9818.07 9701.592 9883.814 9810.926 72.888750 9729.902
## [73,] 10111.08 9836.355 10105.958 9998.117 107.841250 9883.814
## [74,] 10007.91 9927.007 10190.125 10084.878 105.247083 10105.958
## [75,] 9814.53 9937.897 10083.205 10025.082 58.122917 10190.125
## [76,] 9464.68 9849.550 9785.962 9811.397 -25.435000 10083.205
## [77,] 9462.05 9687.292 9415.385 9524.148 -108.762917 9785.962
## [78,] 9515.51 9564.192 9238.291 9368.652 -130.360417 9415.385
## [79,] 9519.35 9490.397 9227.963 9332.937 -104.973750 9238.291
## [80,] 9530.00 9506.727 9414.353 9451.302 -36.950000 9227.963
## [81,] 9594.39 9539.812 9564.029 9554.342 9.686667 9414.353
## [82,] 9659.17 9575.727 9654.996 9623.289 31.707500 9564.029
## [83,] 9742.16 9631.430 9744.773 9699.436 45.337083 9654.996
## [84,] 9838.26 9708.495 9866.210 9803.124 63.085833 9744.773
## [85,] 10032.66 9818.062 10042.452 9952.696 89.755833 9866.210
## [86,] 10080.80 9923.470 10178.646 10076.576 102.070417 10042.452
## [87,] 10082.43 10008.537 10248.365 10152.434 95.930833 10178.646
## [88,] 10018.19 10053.520 10224.558 10156.143 68.415000 10248.365
## [89,] NA NA NA NA NA 10224.558
## [90,] NA NA NA NA NA 10292.973
## [91,] NA NA NA NA NA 10361.388
## [92,] NA NA NA NA NA 10429.803
## [93,] NA NA NA NA NA 10498.218
Hasil peramalan dengan menggunakan N=5 selama 5 periode kedepan pada metode DMA terus naik dari 10224.558 sampai 10498.218
Plot time series
data.gab2 <- ts(data.gab2, start = 2013, frequen = 12)
ts.plot(data.gab2[,1], xlab="Time Period ", ylab="Sales", main= "DMA N=4 Data Premium")
points(data.gab2[,1])
lines(data.gab2[,3],col="green",lwd=2)
lines(data.gab2[,6],col="red",lwd=2)
legend("topleft",c("data aktual","data pemulusan","data peramalan"), lty=8, col=c("black","green","red"), cex=0.8)Menghitung nilai keakuratan
error.dma = test.ts-data.ramal2[1:length(test.ts)]
SSE.dma = sum(error.dma[8:length(test.ts)]^2)
MSE.dma = mean(error.dma[8:length(test.ts)]^2)
MAPE.dma = mean(abs((error.dma[8:length(test.ts)]/test.ts[8:length(test.ts)])*100))
akurasi.dma <- matrix(c(SSE.dma, MSE.dma, MAPE.dma))
row.names(akurasi.dma)<- c("SSE", "MSE", "MAPE")
colnames(akurasi.dma) <- c("Akurasi m = 4")
akurasi.dma## Akurasi m = 4
## SSE 4.039817e+07
## MSE 2.019909e+06
## MAPE 1.377331e+01
akurasi.sma## Akurasi m = 4
## SSE 6.172562e+07
## MSE 2.683723e+06
## MAPE 1.616197e+01
Dengan metode SMA, diketahui SSE, MSE, dan MAPE bernilai lebih kecil Metode DMA lebih baik digunakan
Looping for Best M value
SMA
m = 2:30
akurasi.full <- c()
data.premium <- as.numeric(train.ts)
for(i in m){
data.sma <- SMA(train.ts, n = i)
data.ramal <- c(NA, data.sma)
error.sma = train.ts - data.ramal[1:length(train.ts)]
SSE.sma = sum(error.sma[(i+1):length(train.ts)]^2)
MSE.sma = mean(error.sma[(i+1):length(train.ts)]^2)
MAPE.sma = mean(abs((error.sma[(i+1):length(train.ts)]/train.ts[(i+1):length(train.ts)])*100))
tabel <- matrix(c(SSE.sma, MSE.sma, MAPE.sma))
colnames(tabel) <- paste("M =", i)
rownames(tabel) <- c("SSE", "MSE", "MAPE")
akurasi.full <- cbind(akurasi.full, tabel)
}
akurasi.full## M = 2 M = 3 M = 4 M = 5 M = 6
## SSE 5.404590e+06 7.220289e+06 8.430645e+06 9.231573e+06 9.815792e+06
## MSE 6.284407e+04 8.494458e+04 1.003648e+05 1.112238e+05 1.197048e+05
## MAPE 1.958499e+00 2.288197e+00 2.534714e+00 2.699857e+00 2.805653e+00
## M = 7 M = 8 M = 9 M = 10 M = 11
## SSE 1.026683e+07 1.066054e+07 1.099992e+07 1.117805e+07 1.122272e+07
## MSE 1.267510e+05 1.332567e+05 1.392395e+05 1.433083e+05 1.457496e+05
## MAPE 2.898577e+00 2.940943e+00 2.962072e+00 2.939377e+00 2.917784e+00
## M = 12 M = 13 M = 14 M = 15 M = 16
## SSE 1.131593e+07 1.151234e+07 1.196400e+07 1.275476e+07 1.373600e+07
## MSE 1.488938e+05 1.534979e+05 1.616757e+05 1.747228e+05 1.907777e+05
## MAPE 2.917183e+00 2.912221e+00 2.964623e+00 3.068680e+00 3.212885e+00
## M = 17 M = 18 M = 19 M = 20 M = 21
## SSE 1.466445e+07 1.550042e+07 1.632269e+07 1.709597e+07 1.787718e+07
## MSE 2.065416e+05 2.214345e+05 2.365608e+05 2.514113e+05 2.668237e+05
## MAPE 3.331594e+00 3.416996e+00 3.513380e+00 3.596031e+00 3.668926e+00
## M = 22 M = 23 M = 24 M = 25 M = 26
## SSE 1.849690e+07 1.883418e+07 1.837787e+07 1.609895e+07 1.388503e+07
## MSE 2.802561e+05 2.897566e+05 2.871542e+05 2.555389e+05 2.239521e+05
## MAPE 3.731977e+00 3.757461e+00 3.714365e+00 3.577343e+00 3.454139e+00
## M = 27 M = 28 M = 29 M = 30
## SSE 1.258754e+07 1.248655e+07 1.304894e+07 1.341710e+07
## MSE 2.063530e+05 2.081092e+05 2.211685e+05 2.313292e+05
## MAPE 3.379226e+00 3.383222e+00 3.487242e+00 3.565869e+00
Nilai error yang terkecil merupakan n optimal, terlihat error yang terkecil ada pada N=2
DMA
m = 2:30
akurasi.full2 <- c()
for(i in m){
data.sma <- SMA(train.ts, n = i)
dma <- SMA(data.sma, n = i)
At <- 2*data.sma - dma
Bt <- 2/(i - 1) * (data.sma - dma)
data.dma <- At + Bt
data.ramal2 <- c(NA, data.dma)
error.dma <- test.ts - data.ramal2[1:length(test.ts)]
SSE.dma <- sum(error.dma[(i*2):length(test.ts)]^2)
MSE.dma <- mean(error.dma[(i*2):length(test.ts)]^2)
MAPE.dma <- mean(abs(error.dma[(i*2):115]/test.ts[(i*2):115]) * 100)
tabel2 <- matrix(c(SSE.dma, MSE.dma, MAPE.dma))
colnames(tabel2) <- paste("M =", i)
rownames <- c("SSE", "MSE", "MAPE")
akurasi.full2 <- cbind(akurasi.full2, tabel2)
}
akurasi.full2 <- as.data.frame(akurasi.full2)
akurasi.full2## M = 2 M = 3 M = 4 M = 5 M = 6 M = 7 M = 8 M = 9
## 1 62652793 50900812 40398173 32035511 25409070 20262264 16816485 13675465
## 2 2610533 2313673 2019909 1779751 1588067 1447305 1401374 1367546
## 3 NA NA NA NA NA NA NA NA
## M = 10 M = 11 M = 12 M = 13 M = 14 M = 15 M = 16 M = 17 M = 18 M = 19
## 1 10495589 6475419 2618966.1 719520.3 NA NA NA NA NA NA
## 2 1311949 1079236 654741.5 359760.2 NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA NA NA NA
## M = 20 M = 21 M = 22 M = 23 M = 24 M = 25 M = 26 M = 27 M = 28 M = 29 M = 30
## 1 NA NA NA NA NA NA NA NA NA NA NA
## 2 NA NA NA NA NA NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA NA NA NA NA
Nilai error yang terkecil merupakan n optimal, terlihat error yang terkecil ada pada N=11
Perbandingan SMA dan DMA
comp <- cbind(akurasi.full[,1], akurasi.full2[,10])
colnames(comp) <- c("SMA","DMA")
rownames(comp) <- c("SSE","MSE","MAPE")
comp <- as.data.frame(comp)
comp## SMA DMA
## SSE 5.404590e+06 6475419
## MSE 6.284407e+04 1079236
## MAPE 1.958499e+00 NA
Jika dibandingkan antara SMA dan DMA pada nilai n optimum, metode DMA menghasilkan error yang lebih kecil dibanding SMA pada ketiga ukuran kebaikan diatas.
Pemulusan Eksponensial
Import data
dataekspo<-read_xlsx("C:/0 SEM5/MPDW/p1/Data Kelompok 15.xlsx")
dataekspo$Premium <- as.numeric(dataekspo$Premium)
str(dataekspo)## tibble [115 × 7] (S3: tbl_df/tbl/data.frame)
## $ Tahun : num [1:115] 2013 2013 2013 2013 2013 ...
## $ Bulan : chr [1:115] "Januari" "Februari" "Maret" "April" ...
## $ Premium : num [1:115] 7798 7773 7576 7421 7545 ...
## $ Medium : chr [1:115] "7697.37" "7645.05" "7503.27" "7290.96" ...
## $ Luar Kualitas: chr [1:115] "7545.32" "7328.44" "7033.14" "6870.91" ...
## $ IHK : num [1:115] 1.03 1.79 2.43 2.32 2.3 3.35 6.75 7.94 7.57 7.66 ...
## $ Harga Gabah : chr [1:115] "4333.19" "4265.58" "3783.15" "3669.04" ...
Membentuk objek time series
training<-dataekspo[1:88,3]
testing<-dataekspo[89:115,3]
dataekspo.ts<-ts(dataekspo$Premium, start = 2013, frequency = 12)
training.ts<-ts(training, start = 2013, frequency = 12)
testing.ts<-ts(testing, start = 2020, frequency = 12)Eksplorasi
dataekspo.ts <- as.numeric(dataekspo.ts)
plot(dataekspo.ts, col="red",main="Plot seluruh data")
points(dataekspo.ts)plot(training.ts, col="blue",main="Plot data training")
points(training.ts)plot(testing.ts, col="blue",main="Plot data testing")
points(testing.ts)Single Exponential
Fungsi Holtwinter
ses1<- HoltWinters(training.ts, gamma = FALSE, beta = FALSE, alpha = 0.2)
plot(ses1)Fungsi Holtwinter optimal
sesopt<- HoltWinters(training.ts, gamma = FALSE, beta = FALSE)
sesopt## Holt-Winters exponential smoothing without trend and without seasonal component.
##
## Call:
## HoltWinters(x = training.ts, beta = FALSE, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.9999339
## beta : FALSE
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 10018.19
plot(sesopt)Fungsi Holtwinter
ramalan1<- forecast(ses1, h=10)
ramalan1$mean## Jan Feb Mar Apr May Jun Jul Aug Sep
## 2020 9899.547 9899.547 9899.547 9899.547 9899.547
## 2021 9899.547 9899.547
## Oct Nov Dec
## 2020 9899.547 9899.547 9899.547
## 2021
Fungsi Holtwinter optimal
ramalanopt<- forecast(sesopt, h=10)
ramalanopt## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## May 2020 10018.19 9771.757 10264.63 9641.301 10395.09
## Jun 2020 10018.19 9669.691 10366.70 9485.205 10551.18
## Jul 2020 10018.19 9591.372 10445.02 9365.425 10670.96
## Aug 2020 10018.19 9525.345 10511.04 9264.446 10771.94
## Sep 2020 10018.19 9467.174 10569.21 9175.481 10860.91
## Oct 2020 10018.19 9414.583 10621.81 9095.050 10941.34
## Nov 2020 10018.19 9366.220 10670.17 9021.086 11015.30
## Dec 2020 10018.19 9321.206 10715.18 8952.242 11084.15
## Jan 2021 10018.19 9278.927 10757.46 8887.582 11148.81
## Feb 2021 10018.19 9238.938 10797.45 8826.425 11209.96
Akurasi Data Training
SSE
ses1<- HoltWinters(training.ts, gamma = FALSE, beta = FALSE, alpha = 0.2)
sse1.train <- ses1$SSE
sse1.train## [1] 9033633
sseopt.train <- sesopt$SSE
sseopt.train## [1] 3236764
MSE
mse1.train <- sse1.train/length(training.ts)
mse1.train## [1] 102654.9
mseopt.train <- sseopt.train/length(training.ts)
mseopt.train## [1] 36781.41
RMSE
rmse1.train <- sqrt(mse1.train)
rmse1.train## [1] 320.3981
rmseopt.train <- sqrt(mseopt.train)
rmseopt.train## [1] 191.7848
akurasi.ses1 <- matrix(c(sse1.train, mse1.train, rmse1.train, sseopt.train, mseopt.train, rmseopt.train), nrow=3, ncol=2)
row.names(akurasi.ses1)<- c("SSE", "MSE", "RMSE")
colnames(akurasi.ses1) <- c("alpha 0.2", "alpha 0.99")
akurasi.ses1## alpha 0.2 alpha 0.99
## SSE 9033633.2419 3236764.3104
## MSE 102654.9232 36781.4126
## RMSE 320.3981 191.7848
Metode SES memiliki parameter pemulusan yaitu alpha yang bernilai antara 0 dan 1. Pada awal pemulusan, digunakan alpha 0.2. Lalu dilakukan pemulusan dengan alpha yang optimal secara otomatis dengan program, didapatkan alpha sebesar 0.99 . Parameter pemulusan yang dipilih adalah alpha 0.99 karena memiliki error yang lebih kecil.
Double Exponential
des1<- HoltWinters(training.ts, gamma = FALSE, beta = 0.2, alpha = 0.2)
plot(des1)desopt<- HoltWinters(training.ts, gamma = FALSE)
desopt## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = training.ts, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 1
## beta : 0.01901135
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 10018.19000
## b 12.97931
plot(desopt)Forecasting
ramalandes1<- forecast(des1, h=10)
ramalandes1## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## May 2020 9975.226 9555.017 10395.44 9332.571 10617.88
## Jun 2020 10018.687 9586.545 10450.83 9357.782 10679.59
## Jul 2020 10062.148 9614.275 10510.02 9377.185 10747.11
## Aug 2020 10105.609 9637.985 10573.23 9390.440 10820.78
## Sep 2020 10149.070 9657.586 10640.55 9397.410 10900.73
## Oct 2020 10192.531 9673.100 10711.96 9398.130 10986.93
## Nov 2020 10235.992 9684.635 10787.35 9392.765 11079.22
## Dec 2020 10279.453 9692.361 10866.54 9381.573 11177.33
## Jan 2021 10322.914 9696.477 10949.35 9364.861 11280.97
## Feb 2021 10366.375 9697.198 11035.55 9342.957 11389.79
ramalandesopt<- forecast(desopt, h=10)
ramalandesopt## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## May 2020 10031.17 9780.149 10282.19 9647.267 10415.07
## Jun 2020 10044.15 9685.762 10402.54 9496.043 10592.25
## Jul 2020 10057.13 9614.031 10500.23 9379.469 10734.79
## Aug 2020 10070.11 9553.640 10586.57 9280.238 10859.98
## Sep 2020 10083.09 9500.251 10665.92 9191.716 10974.46
## Oct 2020 10096.07 9451.661 10740.47 9110.534 11081.60
## Nov 2020 10109.05 9406.575 10811.51 9034.710 11183.38
## Dec 2020 10122.02 9364.158 10879.89 8962.967 11281.08
## Jan 2021 10135.00 9323.833 10946.17 8894.426 11375.58
## Feb 2021 10147.98 9285.188 11010.78 8828.451 11467.52
Akurasi data training
SSE
sse.des.train <- des1$SSE
sse.des.train## [1] 9172032
sseopt.des.train <- desopt$SSE
sseopt.des.train## [1] 3305982
MSE
mse.des.train <- sse1.train/length(training.ts)
mse.des.train## [1] 102654.9
mseopt.des.train <- sseopt.train/length(training.ts)
mseopt.des.train## [1] 36781.41
RMSE
rmse.des.train <- sqrt(mse.des.train)
rmse.des.train## [1] 320.3981
rmseopt.des.train <- sqrt(mseopt.des.train)
rmseopt.des.train## [1] 191.7848
akurasi.des1 <- matrix(c(sse.des.train,mse.des.train,rmse.des.train,sseopt.des.train,mseopt.des.train,rmseopt.des.train), nrow=3, ncol=2)
row.names(akurasi.des1)<- c("SSE", "MSE", "RMSE")
colnames(akurasi.des1) <- c("alpha=0.2, beta=0.2", "alpha=1, beta=0.019")
akurasi.des1## alpha=0.2, beta=0.2 alpha=1, beta=0.019
## SSE 9172032.2713 3305981.8294
## MSE 102654.9232 36781.4126
## RMSE 320.3981 191.7848
Metode DES memiliki dua parameter pemulusan yaitu alpha dan beta yang bernilai antara 0 dan 1. Pada awal pemulusan, digunakan alpha dan beta sebesar 0.2. Lalu dilakukan pemulusan dengan alpha dan beta yang optimal secara otomatis dengan program, didapatkan alpha sebesar 1 dan beta sebesar 0.019. Parameter pemulusan yang dipilih adalah alpha = 1 dan beta = 0.0181 karena memiliki error yang lebih kecil.
Perbandingan Akurasi Test Single Exponential dan Double Exponential
selisihses<-ramalanopt$mean-testing.ts
selisihses## Jan Feb Mar Apr May Jun Jul Aug Sep
## 2020 147.0742 205.2942 303.6942 230.3342 238.0042
## 2021 391.1142 481.0342
## Oct Nov Dec
## 2020 246.0942 411.2242 468.4342
## 2021
SSEtestingses<-sum(selisihses^2)
str(testing.ts)## Time-Series [1:27, 1] from 2020 to 2022: 9827 9919 9932 9963 9871 ...
## - attr(*, "dimnames")=List of 2
## ..$ : NULL
## ..$ : chr "Premium"
selisihdes<-ramalandesopt$mean-testing.ts
selisihdes## Jan Feb Mar Apr May Jun Jul Aug Sep
## 2020 160.0493 231.2486 342.6279 282.2472 302.8965
## 2021 507.9238 610.8231
## Oct Nov Dec
## 2020 323.9659 502.0752 572.2645
## 2021
SSEtestingdes<-sum(selisihdes^2)
akurasi <- matrix(c(SSEtestingses, SSEtestingdes), nrow=1, ncol=2)
row.names(akurasi)<- "SSE"
colnames(akurasi) <- c("SES", "DES")
akurasi## SES DES
## SSE 1099169 1683507
Dibandingkan degan metode moving average, metode exponential memiliki error yang jauh lebih kecil, hal ini menandakan bahwa adanya perbedaan pengaruh harga beras antar periode sehingga metode exponential lebih cocok digunakan.
Pemulusan Winter
Import data
winter <- data
str(winter)## tibble [115 × 7] (S3: tbl_df/tbl/data.frame)
## $ Tahun : num [1:115] 2013 2013 2013 2013 2013 ...
## $ Bulan : chr [1:115] "Januari" "Februari" "Maret" "April" ...
## $ Premium : num [1:115] 7798 7773 7576 7421 7545 ...
## $ Medium : chr [1:115] "7697.37" "7645.05" "7503.27" "7290.96" ...
## $ Luar Kualitas: chr [1:115] "7545.32" "7328.44" "7033.14" "6870.91" ...
## $ IHK : num [1:115] 1.03 1.79 2.43 2.32 2.3 3.35 6.75 7.94 7.57 7.66 ...
## $ Harga Gabah : chr [1:115] "4333.19" "4265.58" "3783.15" "3669.04" ...
Membagi data menjadi training dan testing
training<-winter[1:88,3]
testing<-winter[89:115,3]Membentuk objek time series
winter.ts<-ts(winter$Premium,start = 2013, frequency = 12,)
training.ts<-ts(training,start = 2013, frequency = 12)
testing.ts<-ts(testing, start = 2020, frequency = 12)Membuat plot time series
plot(winter.ts, col="red")
points(winter.ts)plot(training.ts, col="blue")
points(training.ts)Winter Aditif
Pemulusan
aditif <- HoltWinters(training.ts, seasonal = "additive")
aditif ## Holt-Winters exponential smoothing with trend and additive seasonal component.
##
## Call:
## HoltWinters(x = training.ts, seasonal = "additive")
##
## Smoothing parameters:
## alpha: 0.9514662
## beta : 0
## gamma: 1
##
## Coefficients:
## [,1]
## a 10245.978564
## b 45.393039
## s1 -212.784340
## s2 -111.756023
## s3 24.000143
## s4 -7.181151
## s5 -67.295108
## s6 -43.632850
## s7 4.759410
## s8 44.169309
## s9 234.497408
## s10 240.376417
## s11 76.866665
## s12 -227.788564
Forecasting
ramalan1 <- forecast(aditif, h=27)
ramalan1## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## May 2020 10078.59 9872.835 10284.34 9763.916 10393.26
## Jun 2020 10225.01 9941.004 10509.01 9790.661 10659.36
## Jul 2020 10406.16 10061.218 10751.10 9878.618 10933.70
## Aug 2020 10420.37 10023.749 10816.99 9813.790 11026.95
## Sep 2020 10405.65 9963.345 10847.95 9729.203 11082.09
## Oct 2020 10474.70 9991.013 10958.40 9734.962 11214.45
## Nov 2020 10568.49 10046.683 11090.30 9770.455 11366.52
## Dec 2020 10653.29 10095.971 11210.61 9800.943 11505.64
## Jan 2021 10889.01 10298.309 11479.72 9985.610 11792.42
## Feb 2021 10940.29 10317.987 11562.58 9988.562 11892.01
## Mar 2021 10822.17 10169.803 11474.53 9824.462 11819.88
## Apr 2021 10562.91 9881.801 11244.01 9521.245 11604.57
## May 2021 10623.30 9911.799 11334.81 9535.151 11711.46
## Jun 2021 10769.73 10031.779 11507.67 9641.135 11898.32
## Jul 2021 10950.87 10187.403 11714.35 9783.246 12118.50
## Aug 2021 10965.09 10176.915 11753.26 9759.683 12170.49
## Sep 2021 10950.37 10138.246 11762.48 9708.336 12192.39
## Oct 2021 11019.42 10184.039 11854.80 9741.815 12297.03
## Nov 2021 11113.21 10255.192 11971.22 9800.988 12425.42
## Dec 2021 11198.01 10317.945 12078.07 9852.068 12543.95
## Jan 2022 11433.73 10532.156 12335.30 10054.891 12812.57
## Feb 2022 11485.00 10562.418 12407.59 10074.032 12895.97
## Mar 2022 11366.89 10423.760 12310.01 9924.499 12809.27
## Apr 2022 11107.62 10144.394 12070.85 9634.491 12580.75
## May 2022 11168.02 10183.061 12152.98 9661.655 12674.39
## Jun 2022 11314.44 10310.216 12318.67 9778.611 12850.27
## Jul 2022 11495.59 10472.462 12518.72 9930.850 13060.33
Akurasi data training
SSE
sse1.train <- aditif$SSE
sse1.train## [1] 1944917
Winter Multiplikatif
Pemulusan
multi <- HoltWinters(training.ts,seasonal = "multiplicative")
multi## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
##
## Call:
## HoltWinters(x = training.ts, seasonal = "multiplicative")
##
## Smoothing parameters:
## alpha: 0.9318329
## beta : 0
## gamma: 1
##
## Coefficients:
## [,1]
## a 1.028539e+04
## b 4.539304e+01
## s1 9.732413e-01
## s2 9.851614e-01
## s3 1.000882e+00
## s4 9.992562e-01
## s5 9.922297e-01
## s6 9.941485e-01
## s7 1.000243e+00
## s8 1.005849e+00
## s9 1.028732e+00
## s10 1.029457e+00
## s11 1.010173e+00
## s12 9.740216e-01
Forecasting
ramalan2 <- forecast(multi, h=27)
ramalan2## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## May 2020 10054.34 9845.032 10263.65 9734.230 10374.45
## Jun 2020 10222.21 9934.474 10509.94 9782.157 10662.25
## Jul 2020 10430.76 10079.609 10781.92 9893.719 10967.81
## Aug 2020 10459.17 10057.965 10860.38 9845.578 11072.77
## Sep 2020 10430.67 9987.007 10874.33 9752.147 11109.19
## Oct 2020 10495.97 10010.565 10981.37 9753.608 11238.33
## Nov 2020 10605.72 10079.897 11131.53 9801.546 11409.88
## Dec 2020 10710.82 10147.285 11274.35 9848.969 11572.67
## Jan 2021 11001.19 10392.946 11609.43 10070.961 11931.41
## Feb 2021 11055.66 10416.513 11694.81 10078.167 12033.16
## Mar 2021 10894.43 10237.457 11551.39 9889.679 11899.17
## Apr 2021 10548.76 9919.530 11177.98 9586.438 11511.07
## May 2021 10584.48 9917.499 11251.47 9564.419 11604.55
## Jun 2021 10758.84 10056.280 11461.40 9684.368 11833.31
## Jul 2021 10975.96 10236.272 11715.65 9844.705 12107.21
## Aug 2021 11003.49 10239.653 11767.32 9835.304 12171.67
## Sep 2021 10971.15 10187.911 11754.39 9773.287 12169.02
## Oct 2021 11037.50 10228.894 11846.10 9800.845 12274.15
## Nov 2021 11150.56 10314.036 11987.09 9871.204 12429.92
## Dec 2021 11258.72 10395.263 12122.18 9938.177 12579.26
## Jan 2022 11561.56 10657.466 12465.64 10178.870 12944.24
## Feb 2022 11616.43 10690.923 12541.93 10200.992 13031.86
## Mar 2022 11444.68 10515.580 12373.79 10023.741 12865.63
## Apr 2022 11079.32 10189.462 11969.18 9718.399 12440.24
## May 2022 11114.62 10194.882 12034.37 9708.000 12521.25
## Jun 2022 11295.47 10344.479 12246.47 9841.054 12749.89
## Jul 2022 11521.16 10535.780 12506.53 10014.153 13028.16
Akurasi data training
SSE
sse2.train <- multi$SSE
sse2.train## [1] 2011779
Akurasi data testing
selisih1<-as.numeric(ramalan1$mean)-as.numeric(testing.ts)
selisih1## [1] 251.9573 306.0586 474.4478 457.7196 534.5287 661.8039 853.9892
## [8] 865.4322 1108.8233 1168.1854 1215.1987 1013.1465 996.2237 1232.5651
## [15] 1549.2643 1465.7260 1494.8051 1569.9504 1573.9157 1525.4687 1609.4998
## [22] 1658.1218 1580.2551 1530.8729 1655.3902 1817.0416 1867.0208
SSEtesting1<-sum(selisih1^2)
selisih2<-as.numeric(ramalan2$mean)-as.numeric(testing.ts)
selisih2## [1] 227.7125 303.2557 499.0526 496.5242 559.5484 683.0679 891.2152
## [8] 922.9575 1220.9981 1283.5635 1287.4558 998.9956 957.4031 1221.6794
## [15] 1574.3497 1504.1255 1515.5922 1588.0269 1611.2741 1586.1802 1737.3256
## [22] 1789.5455 1658.0538 1502.5712 1601.9936 1798.0730 1892.5868
SSEtesting2<-sum(selisih2^2)
akurasi <- matrix(c(SSEtesting1, SSEtesting2), nrow=1, ncol=2)
row.names(akurasi)<- "SSE"
colnames(akurasi) <- c("Aditif", "Multiplikatif")
akurasi## Aditif Multiplikatif
## SSE 44163120 46549188
Kesimpulan
- Metode Exponential lebih cocok digunakan dibandingkan dengan metode moving average. Adanya perbedaan pengaruh harga beras antar periode membuat metode exponential lebih cocok digunakan.
- Dalam Metode Eksponensial, metode Double Exponential Smoothing lebih cocok digunakan dibandingkan Single Exponential Smoothing.
- Pada metode Winter, multiplikative lebih cocok digunakan dibandingkan additive.
Kestasioneran Data
Library
library("tseries")
library("forecast")
library("TTR")
library("TSA")
library("graphics")
library("astsa")
library("car")Import Data
data.premium <- ts(data$Premium, start = 2013, frequency = 12)
ts.plot(data.premium, main = "Harga Beras Premium Tahun 2013 - 2022", ylab = "Harga", lwd = 1.5)
points(data.premium)
legend("topleft", c("Premium"), cex = 0.7,
col = c("black", "red", "blue"), lty = 1)Evaluasi periode musiman
#Plot tahunan
library(forecast)
seasonplot(data.premium,12,main="Xt", ylab="Xt",year.labels = TRUE, col=rainbow(18))monthplot(data.premium,ylab="Xt")boxplot(data$Premium~data$Bulan,ylab="Xt",xlab="Bulan")Berdasarkan Boxplot terlihat bahwa ragam data dari tahun ke tahun cenderung tidak homogen, sehingga untuk memastikan, selanjutnya dilakukan pengujian dengan fligner test.
Plot Transformasi Log
logprem <- log10(data[,"Premium"])
ts.plot(logprem, type="l", ylab="Log Premium")Transformasi logaritma belum berhasil mengatasi ketidakstasioneran dalam ragam dan rataan.
Uji Homogenitas
library(car)
fligner.test(`Premium` ~ Tahun, data=data)##
## Fligner-Killeen test of homogeneity of variances
##
## data: Premium by Tahun
## Fligner-Killeen:med chi-squared = 17.725, df = 9, p-value = 0.0385
fligner.test(log10(Premium) ~ Bulan, data=data)##
## Fligner-Killeen test of homogeneity of variances
##
## data: log10(Premium) by Bulan
## Fligner-Killeen:med chi-squared = 1.8348, df = 11, p-value = 0.999
Berdasarkan hasil uji homogenitas pada data yang belum dilakukan transformasi logaritma di atas, diperoleh nilai p-value=0.0385 < 0.05 sehingga dapat disimpulkan bahwa data harga beras premium tidak homogen.
Oleh karena itu, berdasarkan hasil plot data dan pengujian homogenitas dapat disimpulkan bahwa data harga beras premium tidak stasioner dalam ragam dan rataan.
Sedangkan berdasarkan hasil uji homogenitas pada data yang telah dilakukan transformasi logaritma di atas, diperoleh nilai p-value=0.999 > 0.05 sehingga dapat disimpulkan bahwa data harga beras premium sudah homogen.
Pembagian Data Menjadi Training-Testing
log.train <- subset(log10(data[,"Premium"]),start=1,end=88)
log.test <- subset(log10(data[,"Premium"]),start=89,end=115)Untuk analisis berikutnya, data yang digunakan adalah yang telah ditransformasi logaritma (sebut saja Yt)
ACF Yt
acf(log.train,main="ACF",lag.max=50,xaxt="n")
axis(1, at=0:50, labels=0:50)Plot ACF data setelah transformasi logaritma menunjukkan pola nonstasioner.
ACF Musiman Yt
acf0<-acf(log.train,main="ACF",lag.max=50,xaxt="n")acf0$lag <- acf0$lag * 12
acf0.1 <- as.data.frame(cbind(acf0$acf,acf0$lag))
acf0.2 <- acf0.1[which(acf0.1$V2%%12==0),]
barplot(height = acf0.2$V1, names.arg=acf0.2$V2, ylab="ACF", xlab="Lag")Plot data differencing non musiman
log.train <- ts(log.train)
d1 <- diff(log.train)
ts.plot(d1, type="l", ylab="d1 Beras Premium")differencing non musiman d = 1 berhasil mengatasi ketidakstasioneran dalam rataan untuk komponen non musimannya.
Plot ACF Data d1
acf1 <- acf(d1,lag.max=50,xaxt="n", main="ACF d1")
axis(1, at=0:50, labels=0:50)acf1$lag <- acf1$lag * 12
acf1.1 <- as.data.frame(cbind(acf1$acf,acf1$lag))
acf1.2 <- acf1.1[which(acf1.1$V2%%12==0),]
barplot(height = acf1.2$V1, names.arg = acf1.2$V2, ylab = "ACF", xlab = "Lag")Plot ACF data differencing non musiman d = 1 mengkonfirmasi kestasioneran komponen non musiman
Plot Data Differencing Musiman d1
D1 <- diff(log.train,12)
ts.plot(D1, type="l", ylab="D1 Beras Premium")Plot ACF Data D1
D1 <- diff(log.train,12)
acf2 <- acf(D1,lag.max=50,xaxt="n", main="ACF d1")ts.plot(D1, type="l", ylab="D1 Harga Beras Premium")acf2$lag <- acf2$lag * 12
acf2.1 <- as.data.frame(cbind(acf2$acf,acf2$lag))
acf2.2 <- acf2.1[which(acf2.1$V2%%12==0),]
barplot(height = acf2.2$V1, names.arg=acf2.2$V2, ylab="ACF", xlab="Lag")differencing musiman D = 1 berhasil mengatasi ketidakstasioneran dalam rataan untuk komponen musimannya
Plot Data d1D1
d1D1 <- diff(D1)
ts.plot(d1D1, type="l", ylab="d1 D1 Harga Beras Premium")
# Correlogram Data d1D1
acf3 <- acf(d1D1,lag.max=50,xaxt="n", main="ACF d1D1")
axis(1, at=0:50, labels=0:50)pacf3 <- pacf(d1D1,lag.max=50,xaxt="n", main="PACF d1D1")
axis(1, at=0:50, labels=0:50)Kedua komponen telah stasioner. ARIMA(0,1,4)
Correlogram Musiman Data d1D1
acf3$lag <- acf3$lag * 12
acf3.1 <- as.data.frame(cbind(acf3$acf,acf3$lag))
acf3.2 <- acf3.1[which(acf3.1$V2%%12==0),]
barplot(height = acf3.2$V1, names.arg=acf3.2$V2, ylab="ACF", xlab="Lag")pacf3$lag <- pacf3$lag * 12
pacf3.1 <- as.data.frame(cbind(pacf3$acf,pacf3$lag))
pacf3.2 <- pacf3.1[which(pacf3.1$V2%%12==0),]
barplot(height = pacf3.2$V1, names.arg = pacf3.2$V2, ylab = "ACF", xlab = "Lag")Identifikasi komponen musiman adalah ARIMA()12, sehingga model tentatif adalah ARIMA(0,1,4)?(0,1,1)12.
Pengepasan ARIMA (0,1,4)(0,1,1)12
tmodel <- Arima(log.train,order=c(0,1,4),seasonal=c(0,1,1))
summary(tmodel)## Series: log.train
## ARIMA(0,1,4)
##
## Coefficients:
## ma1 ma2 ma3 ma4
## 0.4381 0.0439 -0.2887 -0.3129
## s.e. 0.0881 0.0920 0.0887 0.0823
##
## sigma^2 = 5.066e-05: log likelihood = 403.85
## AIC=-797.7 AICc=-797.15 BIC=-784.02
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0009743268 0.006960857 0.005191812 0.0245862 0.1310721 0.90942
## ACF1
## Training set -0.01831919
Diagnosa Sisaan
checkresiduals(tmodel)##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,4)
## Q* = 6.0052, df = 6, p-value = 0.4226
##
## Model df: 4. Total lags used: 10
Diagnosa Sisaan
tsdisplay(residuals(tmodel), lag.max=45, main='ARIMA(0,1,2)(0,1,1)12 Model Residuals')Diagnosa Sisaan
library(tseries)
jarque.bera.test(residuals(tmodel))##
## Jarque Bera Test
##
## data: residuals(tmodel)
## X-squared = 9.5284, df = 2, p-value = 0.00853
Validasi Model
ramalan_arima = forecast(tmodel, 12)plot(ramalan_arima)Validasi Model
Plot ACF dan Proses Differencing
ACF
acf(train.ts[,1], lag.max = 50)Berdasarkan plot ACF di atas terlihat bahwa nilai autokorelasi beras premium turun secara perlahan yang mengindikasikan bahwa pola data tidak stasioner.
ACF Musiman
acf0.Zt<-acf(train.ts, main="ACF", lag.max=50, xaxt="n")acf0.Zt$lag <- acf0.Zt$lag * 12
acf0.Zt.1 <- as.data.frame(cbind(acf0.Zt$acf,acf0.Zt$lag))
acf0.Zt.2 <- acf0.Zt.1[which(acf0.Zt.1$V2%%12==0),]
barplot(height = acf0.Zt.2$V1, names.arg = acf0.Zt.2$V2, ylab = "ACF", xlab = "Lag")Berdasarkan plot ACF di atas terlihat bahwa…
Differencing
# Plot data differencing non musiman
d1.Zt <- diff(train.ts)
ts.plot(d1.Zt, type="l", ylab="d1 Zt")Plot data hasil differencing untuk komponen non musiman di atas menunjukkan bahwa, dengan differencing ordo 1 berhasil mengatasi ketidakstasioneran dalam rataan untuk komponen non musiman.
Plot ACF Data Differencing
acf(d1.Zt[,1], lag.max = 50)acf1.Zt <- acf(d1.Zt,lag.max=50,xaxt="n", main="ACF d1.Zt")acf1.Zt$lag <- acf1.Zt$lag * 12
acf1.Zt.1 <- as.data.frame(cbind(acf1.Zt$acf,acf1.Zt$lag))
acf1.Zt.2 <- acf1.Zt.1[which(acf1.Zt.1$V2%%12==0),]
barplot(height = acf1.Zt.2$V1, names.arg=acf1.Zt.2$V2, ylab="ACF", xlab="Lag")ACF untuk data differencing beras premium untuk komponen non musiman menunjukkan bahwa komponen non musiman telah stasioner. Sehingga, selanjutnya dapat dilakukan identifikasi model ARIMA Musiman.
Identifikasi Model (Komponen Non Musiman)
# Correlogram Data d1.Zt
acf3.Zt<- acf(d1.Zt,lag.max=50,xaxt="n", main="ACF d1.Zt")
axis(1, at=0:50, labels=0:50)#PACF
pacf3.Zt <- pacf(d1.Zt,lag.max=50,xaxt="n", main="PACF d1.Zt")
axis(1, at=0:50, labels=0:50)#EACF
eacf(d1.Zt)## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x o x x o o o o o o x x x o
## 1 x o x x o o o o o o o x x o
## 2 x o o x o o o o o o o o o o
## 3 x o o x o o o o o o o o o o
## 4 o o o x o o o o o o o o o o
## 5 o x x o o o o o o o o o o o
## 6 x o x o o o o o o o o o o o
## 7 x o x o o o o x o o o o o o
Berdasarkan plot ACF (cut off setelah lag ke-1), PACF (cut off setelah lag ke-2), dan EACF, diperoleh kandidat model untuk komponen non musiman:
ARIMA (0,1,4)
ARIMA (5,1,3)
ARIMA (7,1,8)
ARIMA (0,1,1)
ARIMA (2,1,0)
Identifikasi Model (Komponen Musiman)
# Correlogram Data d1.Zt
acf3.Zt<- acf(d1.Zt,lag.max=120,xaxt="n", main="ACF d1.Zt")
axis(1, at=0:120, labels=0:120)#PACF
pacf3.Zt <- pacf(d1.Zt,lag.max=120,xaxt="n", main="PACF d1.Zt")
axis(1, at=0:120, labels=0:120)# Correlogram Musiman Data d1.Zt
acf3.Zt$lag <- acf3.Zt$lag * 12
acf3.Zt.1 <- as.data.frame(cbind(acf3.Zt$acf,acf3.Zt$lag))
acf3.Zt.2 <- acf3.Zt.1[which(acf3.Zt.1$V2%%12==0),]
barplot(height = acf3.Zt.2$V1, names.arg=acf3.Zt.2$V2, ylab="ACF", xlab="Lag")# PACF
pacf3.Zt$lag <- pacf3.Zt$lag * 12
pacf3.Zt.1 <- as.data.frame(cbind(pacf3.Zt$acf,pacf3.Zt$lag))
pacf3.Zt.2 <- pacf3.Zt.1[which(pacf3.Zt.1$V2%%12==0),]
barplot(height = pacf3.Zt.2$V1, names.arg=pacf3.Zt.2$V2, ylab="PACF", xlab="Lag")ARIMA (0,0,2)12
ARIMA (1,0,1)12
ARIMA (1,0,2)12
Pendugaan Parameter Model ARIMA Musiman
Model tentatif yang diperoleh berdasarkan hasil identifikasi menggunakan plot ACF, PACF, dan EACF adalah: ARIMA (0,1,4)x(0,0,2)12
ARIMA (5,1,3)x(0,0,2)12
ARIMA (7,1,8)x(0,0,2)12
ARIMA (0,1,1)x(0,0,2)12
ARIMA (2,1,0)x(0,0,2)12
ARIMA (0,1,4)x(1,0,1)12
ARIMA (5,1,3)x(1,0,1)12
ARIMA (7,1,8)x(1,0,1)12
ARIMA (0,1,1)x(1,0,1)12
ARIMA (2,1,0)x(1,0,1)12
ARIMA (0,1,4)x(1,0,2)12
ARIMA (5,1,3)x(1,0,2)12
ARIMA (7,1,8)x(1,0,2)12
ARIMA (0,1,1)x(1,0,2)12
ARIMA (2,1,0)x(1,0,2)12
Model untuk komponen musiman ARIMA (0,0,2)12
#ARIMA(0,1,4)x(0,0,2)12
modela <- Arima(train.ts,order=c(0,1,4),seasonal=list(order=c(1,0,1),period=12))
summary(modela)## Series: train.ts
## ARIMA(0,1,4)(1,0,1)[12]
##
## Coefficients:
## ma1 ma2 ma3 ma4 sar1 sma1
## 0.3433 0.0174 -0.1953 -0.3043 0.9996 -0.9779
## s.e. 0.1035 0.1071 0.1028 0.0973 0.0023 0.0660
##
## sigma^2 = 17338: log likelihood = -554.73
## AIC=1123.46 AICc=1124.88 BIC=1140.72
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.661292 126.3282 96.74171 0.04381646 1.065633 0.7244686
## ACF1
## Training set 0.005728288
#ARIMA(5,1,3)x(0,0,2)12
modelb <- Arima(train.ts,order=c(5,1,3),seasonal=list(order=c(0,0,2),period=12))
summary(modelb)## Series: train.ts
## ARIMA(5,1,3)(0,0,2)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 ma2 ma3 sma1
## 0.1842 -0.2464 0.4393 -0.4857 0.1750 0.2591 0.2453 -0.6750 0.2934
## s.e. 0.2899 0.1567 0.1820 0.1037 0.1693 0.2648 0.2472 0.2324 0.1215
## sma2
## 0.1234
## s.e. 0.1308
##
## sigma^2 = 23283: log likelihood = -556.96
## AIC=1135.93 AICc=1139.45 BIC=1163.05
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 20.66688 142.7323 111.3862 0.2239564 1.22949 0.8341367
## ACF1
## Training set 0.0002868045
#ARIMA(7,1,8)x(0,0,2)12
modelc <- Arima(train.ts,order=c(7,1,8),seasonal=list(order=c(0,0,2),period=12))
summary(modelc)## Series: train.ts
## ARIMA(7,1,8)(0,0,2)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ma1
## 0.7363 -0.0724 -0.3276 -0.1691 0.2922 0.2848 -0.5584 -0.3756
## s.e. 0.2579 0.3184 0.3217 0.2773 0.3295 0.3135 0.1914 0.2787
## ma2 ma3 ma4 ma5 ma6 ma7 ma8 sma1 sma2
## -0.2374 0.0397 0.1948 0.0060 -0.5676 0.3311 0.6113 0.1667 0.0071
## s.e. 0.2656 0.2632 0.3166 0.2448 0.2440 0.2031 0.1444 0.1444 0.1569
##
## sigma^2 = 20391: log likelihood = -551.09
## AIC=1138.17 AICc=1148.23 BIC=1182.56
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 18.65106 127.3572 101.1576 0.2021296 1.116991 0.7575376
## ACF1
## Training set -0.02963039
#ARIMA(0,1,1)x(0,0,2)12
modeld <- Arima(train.ts,order=c(0,1,1),seasonal=list(order=c(0,0,2),period=12))
summary(modeld)## Series: train.ts
## ARIMA(0,1,1)(0,0,2)[12]
##
## Coefficients:
## ma1 sma1 sma2
## 0.4459 0.2772 0.0889
## s.e. 0.0970 0.1103 0.1147
##
## sigma^2 = 26905: log likelihood = -566.23
## AIC=1140.46 AICc=1140.95 BIC=1150.33
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 13.99828 160.2572 120.8831 0.151143 1.329644 0.905256 0.003017726
#ARIMA(2,1,0)x(0,0,2)12
modele <- Arima(train.ts,order=c(2,1,0),seasonal=list(order=c(0,0,2),period=12))
summary(modele)## Series: train.ts
## ARIMA(2,1,0)(0,0,2)[12]
##
## Coefficients:
## ar1 ar2 sma1 sma2
## 0.4694 -0.2523 0.2467 0.0855
## s.e. 0.1064 0.1048 0.1121 0.1149
##
## sigma^2 = 26550: log likelihood = -565.07
## AIC=1140.14 AICc=1140.89 BIC=1152.47
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 15.8629 158.2453 118.8865 0.1699197 1.309519 0.8903041 -0.05123798
printstatarima.Zt <- function (x, digits = 4,se=TRUE){
if (length(x$coef) > 0) {
cat("\nCoefficients:\n")
coef <- round(x$coef, digits = digits)
if (se && nrow(x$var.coef)) {
ses <- rep(0, length(coef))
ses[x$mask] <- round(sqrt(diag(x$var.coef)), digits = digits)
coef <- matrix(coef, 1, dimnames = list(NULL, names(coef)))
coef <- rbind(coef, s.e. = ses)
statt <- coef[1,]/ses
pval <- 2*pt(abs(statt), df=length(x$residuals)-1, lower.tail = FALSE)
coef <- rbind(coef, t=round(statt,digits=digits),sign.=round(pval,digits=digits))
coef <- t(coef)
}
print.default(coef, print.gap = 2)
}
}printstatarima.Zt(modela)##
## Coefficients:
## s.e. t sign.
## ma1 0.3433 0.1035 3.3169 0.0013
## ma2 0.0174 0.1071 0.1625 0.8713
## ma3 -0.1953 0.1028 -1.8998 0.0608
## ma4 -0.3043 0.0973 -3.1274 0.0024
## sar1 0.9996 0.0023 434.6087 0.0000
## sma1 -0.9779 0.0660 -14.8167 0.0000
printstatarima.Zt(modelb)##
## Coefficients:
## s.e. t sign.
## ar1 0.1842 0.2899 0.6354 0.5268
## ar2 -0.2464 0.1567 -1.5724 0.1195
## ar3 0.4393 0.1820 2.4137 0.0179
## ar4 -0.4857 0.1037 -4.6837 0.0000
## ar5 0.1750 0.1693 1.0337 0.3042
## ma1 0.2591 0.2648 0.9785 0.3306
## ma2 0.2453 0.2472 0.9923 0.3238
## ma3 -0.6750 0.2324 -2.9045 0.0047
## sma1 0.2934 0.1215 2.4148 0.0178
## sma2 0.1234 0.1308 0.9434 0.3481
printstatarima.Zt(modelc)##
## Coefficients:
## s.e. t sign.
## ar1 0.7363 0.2579 2.8550 0.0054
## ar2 -0.0724 0.3184 -0.2274 0.8207
## ar3 -0.3276 0.3217 -1.0183 0.3113
## ar4 -0.1691 0.2773 -0.6098 0.5436
## ar5 0.2922 0.3295 0.8868 0.3776
## ar6 0.2848 0.3135 0.9085 0.3661
## ar7 -0.5584 0.1914 -2.9175 0.0045
## ma1 -0.3756 0.2787 -1.3477 0.1813
## ma2 -0.2374 0.2656 -0.8938 0.3739
## ma3 0.0397 0.2632 0.1508 0.8805
## ma4 0.1948 0.3166 0.6153 0.5400
## ma5 0.0060 0.2448 0.0245 0.9805
## ma6 -0.5676 0.2440 -2.3262 0.0223
## ma7 0.3311 0.2031 1.6302 0.1067
## ma8 0.6113 0.1444 4.2334 0.0001
## sma1 0.1667 0.1444 1.1544 0.2515
## sma2 0.0071 0.1569 0.0453 0.9640
printstatarima.Zt(modeld)##
## Coefficients:
## s.e. t sign.
## ma1 0.4459 0.0970 4.5969 0.0000
## sma1 0.2772 0.1103 2.5131 0.0138
## sma2 0.0889 0.1147 0.7751 0.4404
printstatarima.Zt(modele)##
## Coefficients:
## s.e. t sign.
## ar1 0.4694 0.1064 4.4117 0.0000
## ar2 -0.2523 0.1048 -2.4074 0.0182
## sma1 0.2467 0.1121 2.2007 0.0304
## sma2 0.0855 0.1149 0.7441 0.4588
Berdasarkan hasil di atas, dengan mempertimbangkan uji signifikansi penduga parameter dan ukuran keakuratan ramalan dapat disimpulkan bahwa dengan menggunakan model musiman (0,0,2)12 model terbaik adalah ARIMA model C
Model untuk komponen musiman ARIMA (1,0,1)12
modela2 <- Arima(train.ts,order=c(0,1,4),seasonal=list(order=c(1,0,1),period=12))
summary(modela2)## Series: train.ts
## ARIMA(0,1,4)(1,0,1)[12]
##
## Coefficients:
## ma1 ma2 ma3 ma4 sar1 sma1
## 0.3433 0.0174 -0.1953 -0.3043 0.9996 -0.9779
## s.e. 0.1035 0.1071 0.1028 0.0973 0.0023 0.0660
##
## sigma^2 = 17338: log likelihood = -554.73
## AIC=1123.46 AICc=1124.88 BIC=1140.72
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.661292 126.3282 96.74171 0.04381646 1.065633 0.7244686
## ACF1
## Training set 0.005728288
#ARIMA(5,1,3)x(1,0,1)12
modelb2 <- Arima(train.ts,order=c(5,1,3),seasonal=list(order=c(1,0,1),period=12))
summary(modelb2)## Series: train.ts
## ARIMA(5,1,3)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 ma2 ma3
## -0.0590 -0.3089 0.4069 -0.4669 0.0581 0.4473 0.3928 -0.5644
## s.e. 0.5162 0.1932 0.2426 0.0987 0.2726 0.6217 0.4465 0.3030
## sar1 sma1
## 0.9987 -0.9536
## s.e. 0.0085 0.1399
##
## sigma^2 = 15847: log likelihood = -549.66
## AIC=1121.32 AICc=1124.84 BIC=1148.44
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 1.155707 117.7543 94.3834 0.01732897 1.042683 0.706808 0.01425039
#ARIMA(7,1,8)x(1,0,1)12
modelc2 <- Arima(train.ts,order=c(7,1,8),seasonal=list(order=c(1,0,1),period=12))
summary(modelc2)## Series: train.ts
## ARIMA(7,1,8)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ma1
## 0.6304 -0.1539 -0.0171 -0.2284 0.3283 0.2347 -0.4719 -0.3145
## s.e. 0.3098 0.4593 0.4599 0.3770 0.3928 0.4077 0.2431 0.3246
## ma2 ma3 ma4 ma5 ma6 ma7 ma8 sar1
## -0.0820 -0.1953 0.1221 -0.1179 -0.4536 0.4224 0.5803 0.9996
## s.e. 0.3789 0.3742 0.2980 0.3627 0.3852 0.1494 0.2243 0.0022
## sma1
## -0.9801
## s.e. 0.0593
##
## sigma^2 = 15314: log likelihood = -546.45
## AIC=1128.89 AICc=1138.95 BIC=1173.28
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.771469 110.3702 87.65015 0.04371762 0.9696608 0.6563848
## ACF1
## Training set 0.003824925
#ARIMA(0,1,1)x(1,0,1)12
modeld2 <- Arima(train.ts,order=c(0,1,1),seasonal=list(order=c(1,0,1),period=12))
summary(modeld2)## Series: train.ts
## ARIMA(0,1,1)(1,0,1)[12]
##
## Coefficients:
## ma1 sar1 sma1
## 0.3757 0.9994 -0.9720
## s.e. 0.1048 0.0035 0.0817
##
## sigma^2 = 18440: log likelihood = -558.85
## AIC=1125.7 AICc=1126.19 BIC=1135.57
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.105948 132.6735 102.3998 0.03541785 1.127001 0.7668407
## ACF1
## Training set -0.0007157275
#ARIMA(2,1,0)x(1,0,1)12
modele2 <- Arima(train.ts,order=c(2,1,0),seasonal=list(order=c(1,0,1),period=12))
summary(modele2)## Series: train.ts
## ARIMA(2,1,0)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 sar1 sma1
## 0.3750 -0.1831 0.9996 -0.9775
## s.e. 0.1088 0.1072 0.0024 0.0653
##
## sigma^2 = 18502: log likelihood = -558.34
## AIC=1126.68 AICc=1127.42 BIC=1139.01
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.210448 132.1029 100.8782 0.03617671 1.111796 0.7554458
## ACF1
## Training set -0.0209211
printstatarima.Zt(modela2)##
## Coefficients:
## s.e. t sign.
## ma1 0.3433 0.1035 3.3169 0.0013
## ma2 0.0174 0.1071 0.1625 0.8713
## ma3 -0.1953 0.1028 -1.8998 0.0608
## ma4 -0.3043 0.0973 -3.1274 0.0024
## sar1 0.9996 0.0023 434.6087 0.0000
## sma1 -0.9779 0.0660 -14.8167 0.0000
printstatarima.Zt(modelb2)##
## Coefficients:
## s.e. t sign.
## ar1 -0.0590 0.5162 -0.1143 0.9093
## ar2 -0.3089 0.1932 -1.5989 0.1135
## ar3 0.4069 0.2426 1.6772 0.0971
## ar4 -0.4669 0.0987 -4.7305 0.0000
## ar5 0.0581 0.2726 0.2131 0.8317
## ma1 0.4473 0.6217 0.7195 0.4738
## ma2 0.3928 0.4465 0.8797 0.3814
## ma3 -0.5644 0.3030 -1.8627 0.0659
## sar1 0.9987 0.0085 117.4941 0.0000
## sma1 -0.9536 0.1399 -6.8163 0.0000
printstatarima.Zt(modelc2)##
## Coefficients:
## s.e. t sign.
## ar1 0.6304 0.3098 2.0349 0.0449
## ar2 -0.1539 0.4593 -0.3351 0.7384
## ar3 -0.0171 0.4599 -0.0372 0.9704
## ar4 -0.2284 0.3770 -0.6058 0.5462
## ar5 0.3283 0.3928 0.8358 0.4056
## ar6 0.2347 0.4077 0.5757 0.5663
## ar7 -0.4719 0.2431 -1.9412 0.0555
## ma1 -0.3145 0.3246 -0.9689 0.3353
## ma2 -0.0820 0.3789 -0.2164 0.8292
## ma3 -0.1953 0.3742 -0.5219 0.6031
## ma4 0.1221 0.2980 0.4097 0.6830
## ma5 -0.1179 0.3627 -0.3251 0.7459
## ma6 -0.4536 0.3852 -1.1776 0.2422
## ma7 0.4224 0.1494 2.8273 0.0058
## ma8 0.5803 0.2243 2.5872 0.0113
## sar1 0.9996 0.0022 454.3636 0.0000
## sma1 -0.9801 0.0593 -16.5278 0.0000
printstatarima.Zt(modeld2)##
## Coefficients:
## s.e. t sign.
## ma1 0.3757 0.1048 3.5849 6e-04
## sar1 0.9994 0.0035 285.5429 0e+00
## sma1 -0.9720 0.0817 -11.8972 0e+00
printstatarima.Zt(modele2)##
## Coefficients:
## s.e. t sign.
## ar1 0.3750 0.1088 3.4467 0.0009
## ar2 -0.1831 0.1072 -1.7080 0.0912
## sar1 0.9996 0.0024 416.5000 0.0000
## sma1 -0.9775 0.0653 -14.9694 0.0000
Berdasarkan hasil di atas, dengan mempertimbangkan uji signifikansi penduga parameter dan ukuran keakuratan ramalan dapat disimpulkan bahwa dengan menggunakan model musiman (1,0,1)12 model terbaik adalah ARIMA model C
Model untuk komponen musiman ARIMA (1,0,2)12
modela3 <- Arima(train.ts,order=c(0,1,4),seasonal=list(order=c(1,0,2),period=12))
summary(modela3)## Series: train.ts
## ARIMA(0,1,4)(1,0,2)[12]
##
## Coefficients:
## ma1 ma2 ma3 ma4 sar1 sma1 sma2
## 0.3398 0.0347 -0.1807 -0.2830 0.9994 -1.1543 0.1840
## s.e. 0.1040 0.1106 0.1042 0.1029 0.0035 0.2284 0.2063
##
## sigma^2 = 16456: log likelihood = -554.31
## AIC=1124.62 AICc=1126.46 BIC=1144.34
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.7546349 122.3097 95.18991 -0.00285626 1.047509 0.7128477
## ACF1
## Training set 0.005553062
#ARIMA(5,1,3)x(1,0,1)12
modelb3 <- Arima(train.ts,order=c(5,1,3),seasonal=list(order=c(1,0,1),period=12))
summary(modelb3)## Series: train.ts
## ARIMA(5,1,3)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 ma2 ma3
## -0.0590 -0.3089 0.4069 -0.4669 0.0581 0.4473 0.3928 -0.5644
## s.e. 0.5162 0.1932 0.2426 0.0987 0.2726 0.6217 0.4465 0.3030
## sar1 sma1
## 0.9987 -0.9536
## s.e. 0.0085 0.1399
##
## sigma^2 = 15847: log likelihood = -549.66
## AIC=1121.32 AICc=1124.84 BIC=1148.44
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 1.155707 117.7543 94.3834 0.01732897 1.042683 0.706808 0.01425039
tmodel <- Arima(train.ts, order=c(7,1,8),seasonal=c(1,0,2))
summary(tmodel)## Series: train.ts
## ARIMA(7,1,8)
##
## Coefficients:
## Warning in sqrt(diag(x$var.coef)): NaNs produced
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ma1
## 0.6893 0.0214 -0.4828 -0.0071 0.0587 0.4618 -0.6756 -0.336
## s.e. NaN 0.1405 NaN 0.1265 0.1068 0.1233 0.0705 NaN
## ma2 ma3 ma4 ma5 ma6 ma7 ma8
## -0.3149 0.1561 0.0727 0.1744 -0.7132 0.4139 0.5982
## s.e. NaN NaN NaN NaN NaN NaN 0.1218
##
## sigma^2 = 19971: log likelihood = -552.07
## AIC=1136.14 AICc=1143.91 BIC=1175.6
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 23.11811 127.829 99.53925 0.2505431 1.099021 0.7454186 -0.02677237
#ARIMA(7,1,8)x(1,0,1)12
modelc3 <- Arima(train.ts,order=c(7,1,8),seasonal=list(order=c(1,0,1),period=12))
summary(modelc3)## Series: train.ts
## ARIMA(7,1,8)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ma1
## 0.6304 -0.1539 -0.0171 -0.2284 0.3283 0.2347 -0.4719 -0.3145
## s.e. 0.3098 0.4593 0.4599 0.3770 0.3928 0.4077 0.2431 0.3246
## ma2 ma3 ma4 ma5 ma6 ma7 ma8 sar1
## -0.0820 -0.1953 0.1221 -0.1179 -0.4536 0.4224 0.5803 0.9996
## s.e. 0.3789 0.3742 0.2980 0.3627 0.3852 0.1494 0.2243 0.0022
## sma1
## -0.9801
## s.e. 0.0593
##
## sigma^2 = 15314: log likelihood = -546.45
## AIC=1128.89 AICc=1138.95 BIC=1173.28
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.771469 110.3702 87.65015 0.04371762 0.9696608 0.6563848
## ACF1
## Training set 0.003824925
#ARIMA(0,1,1)x(1,0,1)12
modeld3 <-Arima(train.ts,order=c(0,1,1),seasonal=list(order=c(1,0,1),period=12))
summary(modeld3)## Series: train.ts
## ARIMA(0,1,1)(1,0,1)[12]
##
## Coefficients:
## ma1 sar1 sma1
## 0.3757 0.9994 -0.9720
## s.e. 0.1048 0.0035 0.0817
##
## sigma^2 = 18440: log likelihood = -558.85
## AIC=1125.7 AICc=1126.19 BIC=1135.57
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.105948 132.6735 102.3998 0.03541785 1.127001 0.7668407
## ACF1
## Training set -0.0007157275
#ARIMA(2,1,0)x(1,0,1)12
modele3 <- Arima(train.ts,order=c(2,1,0),seasonal=list(order=c(1,0,1),period=12))
summary(modele3)## Series: train.ts
## ARIMA(2,1,0)(1,0,1)[12]
##
## Coefficients:
## ar1 ar2 sar1 sma1
## 0.3750 -0.1831 0.9996 -0.9775
## s.e. 0.1088 0.1072 0.0024 0.0653
##
## sigma^2 = 18502: log likelihood = -558.34
## AIC=1126.68 AICc=1127.42 BIC=1139.01
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 3.210448 132.1029 100.8782 0.03617671 1.111796 0.7554458
## ACF1
## Training set -0.0209211
printstatarima.Zt(modela3)##
## Coefficients:
## s.e. t sign.
## ma1 0.3398 0.1040 3.2673 0.0016
## ma2 0.0347 0.1106 0.3137 0.7545
## ma3 -0.1807 0.1042 -1.7342 0.0864
## ma4 -0.2830 0.1029 -2.7502 0.0072
## sar1 0.9994 0.0035 285.5429 0.0000
## sma1 -1.1543 0.2284 -5.0539 0.0000
## sma2 0.1840 0.2063 0.8919 0.3749
printstatarima.Zt(modelb3)##
## Coefficients:
## s.e. t sign.
## ar1 -0.0590 0.5162 -0.1143 0.9093
## ar2 -0.3089 0.1932 -1.5989 0.1135
## ar3 0.4069 0.2426 1.6772 0.0971
## ar4 -0.4669 0.0987 -4.7305 0.0000
## ar5 0.0581 0.2726 0.2131 0.8317
## ma1 0.4473 0.6217 0.7195 0.4738
## ma2 0.3928 0.4465 0.8797 0.3814
## ma3 -0.5644 0.3030 -1.8627 0.0659
## sar1 0.9987 0.0085 117.4941 0.0000
## sma1 -0.9536 0.1399 -6.8163 0.0000
printstatarima.Zt(modelc3)##
## Coefficients:
## s.e. t sign.
## ar1 0.6304 0.3098 2.0349 0.0449
## ar2 -0.1539 0.4593 -0.3351 0.7384
## ar3 -0.0171 0.4599 -0.0372 0.9704
## ar4 -0.2284 0.3770 -0.6058 0.5462
## ar5 0.3283 0.3928 0.8358 0.4056
## ar6 0.2347 0.4077 0.5757 0.5663
## ar7 -0.4719 0.2431 -1.9412 0.0555
## ma1 -0.3145 0.3246 -0.9689 0.3353
## ma2 -0.0820 0.3789 -0.2164 0.8292
## ma3 -0.1953 0.3742 -0.5219 0.6031
## ma4 0.1221 0.2980 0.4097 0.6830
## ma5 -0.1179 0.3627 -0.3251 0.7459
## ma6 -0.4536 0.3852 -1.1776 0.2422
## ma7 0.4224 0.1494 2.8273 0.0058
## ma8 0.5803 0.2243 2.5872 0.0113
## sar1 0.9996 0.0022 454.3636 0.0000
## sma1 -0.9801 0.0593 -16.5278 0.0000
printstatarima.Zt(modeld3)##
## Coefficients:
## s.e. t sign.
## ma1 0.3757 0.1048 3.5849 6e-04
## sar1 0.9994 0.0035 285.5429 0e+00
## sma1 -0.9720 0.0817 -11.8972 0e+00
printstatarima.Zt(modele3)##
## Coefficients:
## s.e. t sign.
## ar1 0.3750 0.1088 3.4467 0.0009
## ar2 -0.1831 0.1072 -1.7080 0.0912
## sar1 0.9996 0.0024 416.5000 0.0000
## sma1 -0.9775 0.0653 -14.9694 0.0000
Berdasarkan hasil di atas, dengan mempertimbangkan uji signifikansi penduga parameter dan ukuran keakuratan ramalan dapat disimpulkan bahwa dengan menggunakan model musiman (1,0,2)12 model terbaik adalah ARIMA model C
Diagnostik Model
#diagnostik model
checkresiduals(modelc)##
## Ljung-Box test
##
## data: Residuals from ARIMA(7,1,8)(0,0,2)[12]
## Q* = 6.2366, df = 3, p-value = 0.1006
##
## Model df: 17. Total lags used: 20
tsdisplay(residuals(modelc), lag.max=50, main='ARIMA C Model Residuals')library(tseries)
jarque.bera.test(residuals(modelc))##
## Jarque Bera Test
##
## data: residuals(modelc)
## X-squared = 3.1281, df = 2, p-value = 0.2093
Selanjutnya, apabila dilihat melalui histogram sisaan, terlihat bahwa sisaaan menyebar normal, selain itu, berdasarkan pengujian menggunakan Jarque bera Test diperoleh p-value = 0.1658 > 0.05 yang berarti bahwa H0 tidak ditolak sehingga dapat disimpulkan sisaan model menyebar normal.
#validasi model
ramalan_arima3 = forecast(modelc, 120)
ramalan_arima3## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 89 10021.384 9834.756 10208.01 9735.961 10306.81
## 90 10075.253 9755.998 10394.51 9586.995 10563.51
## 91 10115.971 9708.210 10523.73 9492.354 10739.59
## 92 10072.796 9623.828 10521.76 9386.159 10759.43
## 93 9984.161 9519.209 10449.11 9273.079 10695.24
## 94 9940.803 9459.413 10422.19 9204.581 10677.03
## 95 9997.291 9500.126 10494.46 9236.943 10757.64
## 96 10123.022 9609.808 10636.24 9338.130 10907.91
## 97 10235.421 9688.371 10782.47 9398.780 11072.06
## 98 10229.826 9647.732 10811.92 9339.590 11120.06
## 99 10154.331 9532.320 10776.34 9203.047 11105.61
## 100 10099.290 9428.589 10769.99 9073.541 11125.04
## 101 10133.449 9402.998 10863.90 9016.321 11250.58
## 102 10228.984 9433.287 11024.68 9012.071 11445.90
## 103 10286.067 9425.052 11147.08 8969.259 11602.87
## 104 10231.780 9322.016 11141.54 8840.417 11623.14
## 105 10115.064 9166.743 11063.38 8664.732 11565.40
## 106 10043.774 9062.557 11024.99 8543.133 11544.41
## 107 10083.114 9074.132 11092.10 8540.009 11626.22
## 108 10191.618 9155.073 11228.16 8606.359 11776.88
## 109 10259.814 9194.201 11325.43 8630.100 11889.53
## 110 10219.648 9127.906 11311.39 8549.972 11889.32
## 111 10119.268 9001.616 11236.92 8409.967 11828.57
## 112 10064.194 8920.484 11207.90 8315.041 11813.35
## 113 10115.594 8944.294 11286.89 8324.245 11906.94
## 114 10226.114 9023.317 11428.91 8386.594 12065.63
## 115 10285.855 9049.799 11521.91 8395.470 12176.24
## 116 10235.363 8967.964 11502.76 8297.044 12173.68
## 117 10126.641 8829.844 11423.44 8143.362 12109.92
## 118 10067.754 8743.855 11391.65 8043.025 12092.48
## 119 10116.733 8766.194 11467.27 8051.262 12182.20
## 120 10221.513 8843.372 11599.65 8113.828 12329.20
## 121 10273.327 8868.112 11678.54 8124.235 12422.42
## 122 10218.292 8788.119 11648.46 8031.032 12405.55
## 123 10111.441 8658.551 11564.33 7889.438 12333.44
## 124 10060.322 8586.301 11534.34 7806.001 12314.64
## 125 10117.139 8621.443 11612.84 7829.670 12404.61
## 126 10224.609 8705.693 11743.52 7901.628 12547.59
## 127 10274.594 8731.937 11817.25 7915.304 12633.89
## 128 10217.817 8652.264 11783.37 7823.510 12612.12
## 129 10112.950 8526.013 11699.89 7685.940 12539.96
## 130 10067.019 8459.551 11674.49 7608.608 12525.43
## 131 10127.072 8498.264 11755.88 7636.026 12618.12
## 132 10232.048 8580.594 11883.50 7706.367 12757.73
## 133 10275.397 8601.020 11949.78 7714.658 12836.14
## 134 10213.082 8516.952 11909.21 7619.075 12807.09
## 135 10107.941 8391.885 11824.00 7483.460 12732.42
## 136 10066.120 8331.096 11801.14 7412.630 12719.61
## 137 10129.451 8374.920 11883.98 7446.127 12812.78
## 138 10233.113 8458.046 12008.18 7518.382 12947.84
## 139 10271.841 8476.017 12067.66 7525.366 13018.32
## 140 10206.498 8390.983 12022.01 7429.907 12983.09
## 141 10103.555 8269.878 11937.23 7299.188 12907.92
## 142 10067.589 8216.362 11918.82 7236.382 12898.80
## 143 10135.099 8265.617 12004.58 7275.972 12994.23
## 144 10237.650 8348.788 12126.51 7348.885 13126.42
## 145 10271.502 8362.958 12180.05 7352.636 13190.37
## 146 10202.651 8275.404 12129.90 7255.182 13150.12
## 147 10101.163 8156.599 12045.73 7127.209 13075.12
## 148 10069.968 8108.606 12031.33 7070.323 13069.61
## 149 10140.455 8161.633 12119.28 7114.108 13166.80
## 150 10240.895 8243.607 12238.18 7186.307 13295.48
## 151 10269.374 8253.446 12285.30 7186.278 13352.47
## 152 10197.087 8163.539 12230.64 7087.043 13307.13
## 153 10097.414 8047.554 12147.27 6962.424 13232.40
## 154 10071.338 8005.612 12137.06 6912.083 13230.59
## 155 10145.063 8062.805 12227.32 6960.525 13329.60
## 156 10243.668 8143.914 12343.42 7032.371 13454.97
## 157 10267.223 8149.826 12384.62 7028.943 13505.50
## 158 10192.071 8057.988 12326.15 6928.273 13455.87
## 159 10094.628 7945.037 12244.22 6807.112 13382.14
## 160 10073.715 7908.952 12238.48 6762.995 13384.44
## 161 10150.379 7969.747 12331.01 6815.390 13485.37
## 162 10246.759 8049.334 12444.18 6886.087 13607.43
## 163 10265.180 8050.865 12479.49 6878.678 13651.68
## 164 10187.157 7956.895 12417.42 6776.265 13598.05
## 165 10091.950 7846.842 12337.06 6658.354 13525.55
## 166 10076.032 7816.356 12335.71 6620.156 13531.91
## 167 10155.324 7880.387 12430.26 6676.109 13634.54
## 168 10249.254 7958.196 12540.31 6745.383 13753.13
## 169 10262.609 7955.389 12569.83 6734.021 13791.20
## 170 10182.031 7859.572 12504.49 6630.136 13733.93
## 171 10089.367 7752.693 12426.04 6515.732 13663.00
## 172 10078.522 7727.840 12429.20 6483.464 13673.58
## 173 10160.301 7794.909 12525.69 6542.746 13777.86
## 174 10251.650 7870.730 12632.57 6610.347 13892.95
## 175 10260.010 7863.560 12656.46 6594.957 13925.06
## 176 10177.114 7766.032 12588.19 6489.683 13864.54
## 177 10087.162 7662.400 12511.92 6378.809 13795.52
## 178 10081.307 7643.010 12519.61 6352.253 13810.36
## 179 10165.289 7712.750 12617.83 6414.455 13916.12
## 180 10253.817 7786.271 12721.36 6480.031 14027.60
## 181 10257.190 7774.677 12739.70 6460.514 14053.87
## 182 10172.182 7675.585 12668.78 6353.967 13990.40
## 183 10085.115 7575.322 12594.91 6246.718 13923.51
## 184 10084.196 7561.300 12607.09 6225.760 13942.63
## 185 10170.154 7633.451 12706.86 6290.602 14049.71
## 186 10255.703 7704.473 12806.93 6353.933 14157.47
## 187 10254.172 7688.494 12819.85 6330.306 14178.04
## 188 10167.314 7588.053 12746.57 6222.675 14111.95
## 189 10083.323 7491.307 12675.34 6119.176 14047.47
## 190 10087.287 7482.555 12692.02 6103.694 14070.88
## 191 10174.989 7556.838 12793.14 6170.873 14179.10
## 192 10257.386 7625.142 12889.63 6231.715 14283.06
## 193 10251.016 7604.792 12897.24 6203.966 14298.07
## 194 10162.538 7503.183 12821.89 6095.405 14229.67
## 195 10081.778 7410.061 12753.49 5995.740 14167.82
## 196 10090.537 7406.455 12774.62 5985.587 14195.49
## 197 10179.734 7482.586 12876.88 6054.801 14304.67
## 198 10258.817 7547.970 12969.66 6112.935 14404.70
## 199 10247.688 7523.294 12972.08 6081.087 14414.29
## 200 10157.830 7420.721 12894.94 5971.783 14343.88
## 201 10080.459 7331.350 12829.57 5876.060 14284.86
## 202 10093.931 7332.775 12855.09 5871.108 14316.75
## 203 10184.387 7410.491 12958.28 5942.078 14426.70
## 204 10260.010 7472.782 13047.24 5997.313 14522.71
## 205 10244.219 7443.843 13044.59 5961.413 14527.02
## 206 10153.224 7340.514 12965.93 5851.556 14454.89
## 207 10079.386 7255.007 12903.76 5759.871 14398.90
## 208 10097.472 7261.338 12933.61 5759.979 14434.96
plot(ramalan_arima3)Berdasarkan ukuran keakuratan ramalan di atas, dapat disimpulkan bahwa model ARIMA()x()12 dapat digunakan untuk meramalkan data harga beras premium dengan baik. Selain itu, apabila dilihat berdasarkan plot ramalan, nilai ramalan yang dihasilkan memiliki nilai yang .. dengan nilai data sebelumnya.
Peramalan (5 bulan ke depan)
Berikut merupakan hasil ramalan data harga beras premium untuk 5 bulan ke depan:
modelZt<-Arima(data$Premium,order=c(7,1,8),seasonal=c(1,0,2))
ramalan_arimaZt = forecast(modelZt, 5)
ramalan_arimaZt## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 116 9726.153 9554.557 9897.749 9463.719 9988.587
## 117 9809.721 9514.223 10105.218 9357.796 10261.645
## 118 9907.472 9524.601 10290.343 9321.921 10493.023
## 119 9930.126 9506.461 10353.792 9282.186 10578.067
## 120 9903.078 9463.719 10342.437 9231.137 10575.020
plot(ramalan_arimaZt)Daftar Pustaka
Fani E, Widjajati FA, Soehardjoepri. 2017. Perbandingan Metode Winter Eksponensial Smoothing dan Metode Event Based untuk Menentukan Penjualan Produk Terbaik di Perusahaan X. Jurnal Sains dan Seni ITS. 6(1):2337-3529.
H. R. Naufal and R Adrean, “R. naufal Hayâ and R. Adrean, ‘Sistem Informasi Inventory Berdasarkan Prediksi Data Penjualan Barang Menggunakan Metode Single Moving Average Pada CV. Agung Youanda,’ ProTekInfo (Pengembangan Ris. dan Obs.Tek. Inform., vol. 4, pp. 29??”33, 2017,” ProTekInfo (Pengembangan Ris. dan Obs. Tek.Inf., vol. 4, pp. 29??“33, 2017.
Hartono A, Dwijana D, Handiwidjojo. 2012. Perbandingan metode single exponential smoothing dan exponential smoothing adjusted for trend untuk meramalkan penjualan. Studi kasus: Toko onderdil mobil “Prodi, Purwodadi”. Jurnal EKSIS. 5(1):8-18.
Hapsari, V. 2013. Perbandingan Metode Dekomposisi Klasik dengan Metode Pemulusan Eksponensial Holt-Winter dalam meramalkan Tingkat Pencemaran Udara di Kota Bandung Periode 2003-2012. Universitas Pendidikan Indonesia.
Muzaki L. 2022. Contoh studi kasus metode double moving average, kelebihan kekurangan, dan cara menghitung [internet]. [diacu pada 2022 Oktober 25]. Tersedia dari: https://www.pengadaanbarang.co.id/2022/05/metode-double-moving-average.html.
M. Layakana, S. Iskandar. “Penerapan Metode Double Moving Average dan Double Eksponensial Smoothing Dalam Meramalkan Jumlah Produksi Crude Palm Oil(CPO) Pada PT.Perkebunan Nusantara IV Unit Dolok Sinumbah. 6(1): 47.
Purwanto A, Hanief S. 2017. Teknik peramalan dengan double exponential smoothing pada distributor gula. Jurnal Teknologi Informasi dan Komputer. 3(1):362-367