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

  1. 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

Plot Time Series Harga Beras Premium

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)

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

  1. Metode Exponential lebih cocok digunakan dibandingkan dengan metode moving average. Adanya perbedaan pengaruh harga beras antar periode membuat metode exponential lebih cocok digunakan.
  2. Dalam Metode Eksponensial, metode Double Exponential Smoothing lebih cocok digunakan dibandingkan Single Exponential Smoothing.
  3. 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

https://www.bps.go.id/indicator/36/500/10/rata-rata-harga-beras-bulanan-di-tingkat-penggilingan-menurut-kualitas.html