TUGAS PRAKTIKUM 8

Install Packagaes

library(readxl)
library(forecast)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(TTR)
library(imputeTS)
library(tseries)
## 
## Attaching package: 'tseries'
## The following object is masked from 'package:imputeTS':
## 
##     na.remove
library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(graphics)
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(tidyverse)
## ── Attaching packages
## ───────────────────────────────────────
## tidyverse 1.3.2 ──
## ✔ tibble  3.1.8     ✔ purrr   0.3.4
## ✔ tidyr   1.2.0     ✔ stringr 1.4.0
## ✔ readr   2.1.2     ✔ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ✖ readr::spec()   masks TSA::spec()
library(lubridate)
## 
## Attaching package: 'lubridate'
## 
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     combine
library(ggfortify)
## Registered S3 methods overwritten by 'ggfortify':
##   method                 from    
##   autoplot.Arima         forecast
##   autoplot.acf           forecast
##   autoplot.ar            forecast
##   autoplot.bats          forecast
##   autoplot.decomposed.ts forecast
##   autoplot.ets           forecast
##   autoplot.forecast      forecast
##   autoplot.stl           forecast
##   autoplot.ts            forecast
##   fitted.ar              forecast
##   fortify.ts             forecast
##   residuals.ar           forecast
library(cowplot)
## 
## Attaching package: 'cowplot'
## 
## The following object is masked from 'package:lubridate':
## 
##     stamp
library(graphics)
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## 
## The following object is masked from 'package:imputeTS':
## 
##     na.locf
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(stats)
library(MASS)
## 
## Attaching package: 'MASS'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
library(FinTS)
## 
## Attaching package: 'FinTS'
## 
## The following object is masked from 'package:forecast':
## 
##     Acf
library(orcutt)

Latar Belakang

Data yang digunakan pada proses forecasting adalah Kurs Transaksi Mata Uang yang dilakukan oleh Bank Indonesia pada rentang periode 1 tahun dengan perhitungan data harian, 1 Januari hingga 31 Desember 2021 terhadap mata uang Britania Raya, Pound sterling.

Import Data

data <- read_excel("~/Downloads/kurs transaksi GBP.xlsx")
ts.data <- ts(data)
head(data)
## # A tibble: 6 × 2
##   Tanggal             Kurs_Jual
##   <dttm>                  <dbl>
## 1 2021-01-04 00:00:00    19128.
## 2 2021-01-05 00:00:00    19038.
## 3 2021-01-06 00:00:00    19044.
## 4 2021-01-07 00:00:00    19024.
## 5 2021-01-08 00:00:00    19158.
## 6 2021-01-11 00:00:00    19225.
str(data)
## tibble [256 × 2] (S3: tbl_df/tbl/data.frame)
##  $ Tanggal  : POSIXct[1:256], format: "2021-01-04" "2021-01-05" ...
##  $ Kurs_Jual: num [1:256] 19128 19038 19044 19024 19158 ...

Eksplorasi Data

data$Tanggal <- as.Date(data$Tanggal)
ggplot(data, aes(x=Tanggal, y= Kurs_Jual))+
  geom_line()+
  scale_x_date(date_breaks = "2 month", date_labels = "%Y %b %d") +
  labs(title = "Plot Time Series Data Kurs Transaksi GBP",
       subtitle = "(Januari 2021 sd Desember 2021)",
       y="Kurs Jual GBP") +
  theme(plot.title =  element_text(face = "bold", hjust=.5),
        plot.subtitle = element_text(hjust=.5))

p=0.8
freq_train=as.integer(p*nrow(data))
data$Tanggal <- as.Date(data$Tanggal)
ggplot(data, aes(x=Tanggal, y=Kurs_Jual)) +
  geom_line() + 
  scale_x_date(date_breaks = "2 month", date_labels = "%Y %b %d") +
  labs(title = "Plot Time Series Data Kurs Transaksi GBP",
       subtitle = "(Januari 2021 sd Desember 2021)",
       y="Kurs Jual GBP") +
  geom_vline(aes(xintercept = Tanggal[freq_train], 
                 col="Frequency_Train_Data"), lty=2, lwd=.7) +
  theme(plot.title =  element_text(face = "bold", hjust=.5),
        plot.subtitle = element_text(hjust=.5),
        legend.position = "bottom") +
  scale_color_manual(name = "", values = c(Frequency_Train_Data="red"))

data[freq_train,1]
## # A tibble: 1 × 1
##   Tanggal   
##   <date>    
## 1 2021-10-19

Berdasarkan hasil eksplorasi, terlihat bahwa sebaran data cenderung mengikuti pola siklik. Sekilas terlihat mirip namun pola siklik berbeda dengan pola musiman, pada pola musiman mempunyai panjang gelombang yang tetap dan terjadi pada waktu yang tetap, sedangkan pola siklik memiliki jarak waktu yang bervariasi dari satu siklus ke siklus lainnya (Ilham 2016). Pola data siklik tersebut menunjukkan bahwa siklus transaksi Kurs Mata Uang terus berfluktuasi dalam jangka panjang. Walaupun di bagian akhir periode, data cenderung mengalami penurunan, tetapi datanya tetap berfluktuasi sehingga lebih tepat dikatakan siklik daripada tren menurun.

#SPlitting Data

training <- data[1:232,]
testing<-data[233:256,]

training.ts<-ts(training)
testing.ts<-ts(testing)

Dari keseluruhan dataset dilakukan pengaturan komposisi menjadi 2 sub bagian yaitu data training dan testing dengan presentase 80 : 20 dengan tujuan menghindari adanya overfitting.

Training

plot(training.ts, col="blue",main="Plot Data Training")
points(training.ts)

Testing

plot(testing.ts, col="blue",main="Plot dataTtraining")
points(testing.ts)

Time Series Plot Data Training dan Testing

p=0.8
freq_train=as.integer(p*nrow(data))
data$Tanggal <- as.Date(data$Tanggal)
ggplot(data, aes(x=Tanggal, y=Kurs_Jual)) +
  geom_line() + 
  scale_x_date(date_breaks = "2 month", date_labels = "%Y %b %d") +
  labs(title = "Plot Time Series Data Kurs Transaksi GBP",
       subtitle = "(Januari 2021 sd Desember 2021)",
       y="Kurs Jual GBP") +
  geom_vline(aes(xintercept = Tanggal[freq_train], 
                 col="Frequency_Train_Data"), lty=2, lwd=.7) +
  theme(plot.title =  element_text(face = "bold", hjust=.5),
        plot.subtitle = element_text(hjust=.5),
        legend.position = "bottom") +
  scale_color_manual(name = "", values = c(Frequency_Train_Data="red"))

Pemulusan

pemulusan <- function(x,min,max,method=c("SMA","DMA"),
                      metrik=c("SSE","MSE","RMSE","MAD","MAPE")){
  df.master <- data.frame()
  z <- max-min+1
  temp <- sqrt(z)
  a <- round(temp) 
  b=a
  if(a*b<z){
    b=b+1
  }
  par(mfrow=c(a,b))
  if(method=="SMA"){
    for(i in min:max) {
      sma <- TTR::SMA(x,i)
      ramal <- c(NA,sma)
      df <- cbind(aktual=c(x,NA),mulus=c(sma,NA),ramal)
      error <- df[,1]-df[,3]
      sse <- sum(error^2,na.rm=T)
      mse <- mean(error^2,na.rm=T)
      rmse <- sqrt(mse)
      mad <- mean(abs(error),na.rm=T)
      rerror <- error/df[,1]*100
      mape <- mean(abs(rerror),na.rm=T)
      ak <- data.frame(n=i,SSE=sse,MSE=mse,RMSE=rmse,MAD=mad,MAPE=mape)
      df.master <- rbind(df.master,ak)
      ts.plot(x,col="gray",main=paste("Single Moving Average n =",i))
      lines(sma,col="green",lwd=2)
      lines(ramal,col="red",lwd=1)
    }
  }
  else if(method=="DMA") {
    for(i in min:max) {
      df_ts_sma <- TTR::SMA(x,  n=i)
      df_ts_dma <- TTR::SMA(df_ts_sma, n=i)
      At <- 2*df_ts_sma-df_ts_dma
      Bt <- df_ts_sma-df_ts_dma
      pemulusan_dma <- At+Bt
      ramal_dma <- c(NA, pemulusan_dma)
      df_dma <- cbind(df_aktual=c(x,NA), pemulusan_dma=c(pemulusan_dma,NA), ramal_dma)
      error.dma <- df_dma[, 1] - df_dma[, 3]
      SSE.dma <- sum(error.dma^2, na.rm = T)
      MSE.dma <- mean(error.dma^2, na.rm = T)
      RMSE.dma <- sqrt(mean(error.dma^2, na.rm = T))
      MAD.dma <- mean(abs(error.dma), na.rm = T)
      r.error.dma <- (error.dma/df_dma[, 1])*100 
      MAPE.dma <- mean(abs(r.error.dma), na.rm = T)
      ak <- data.frame(n=i,SSE=SSE.dma,MSE=MSE.dma,RMSE=RMSE.dma,
                       MAD=MAD.dma,MAPE=MAPE.dma)
      df.master <- rbind(df.master,ak)
      ts.plot(x,main=paste("Double Moving Average n =",i))
      lines(pemulusan_dma,col="green",lwd=2)
      lines(ramal_dma,col="red",lwd= 2)
    }
  }
  opt <- df.master$n[which.min(df.master[,metrik])]
  return(list(Akurasi=df.master,n.optimum=paste("n optimum adalah",opt)))
}

pemulusan2 <- function(x,alfa=NULL,beta=NULL,method=c("SES","DES"),
                       metrik=c("SSE","MSE","RMSE")){
  df.master <- data.frame()
  if(method=="SES"){
    for(i in alfa){
      df_ses <- HoltWinters(x, alpha = i, beta=F, gamma=F)
      sse <- df_ses$SSE 
      mse <- sse/length(x)
      rmse <-sqrt(mse)
      ak <- data.frame(Alpha=i,SSE=sse,MSE=mse,RMSE=rmse)
      df.master <- rbind(df.master,ak)
      datases <- data.frame(x, c(NA, df_ses$fitted[,1]))
      colnames(datases) = c("y","yhat")
      ts.plot(x, xlab="Periode  waktu", ylab="Yt", col="blue", lty=3,main=paste("Single Exponential Smoothing Alpha=",i))
      points(x)
      lines(datases[,2], col="red",lwd=2) #nilai dugaan
    }
  }
  else if(method=="DES"){
    for (i in alfa) {
      for (j in beta) {
        df_des <- HoltWinters(x, alpha = i, beta=j, gamma=F)
        sse <- df_des$SSE 
        mse <- sse/length(x)
        rmse <-sqrt(mse)
        ak <- data.frame(Alpha=i,Beta=j,SSE=sse,MSE=mse,RMSE=rmse)
        df.master <- rbind(df.master,ak)
        datades <- data.frame(x, c(NA, NA, df_des$fitted[,1]))
        colnames(datades) = c("y","yhat")
        ts.plot(x,xlab="periode  waktu", ylab="Yt",  col="blue", lty=3,main=paste("Double Exponential Smoothing (Alpha=",i,"Beta=",j))
        points(x)
        lines (datades[,2], col="red",lwd=2)
      }
    }
  }
  if(method=="SES"){
    opt <- df.master$Alpha[which.min(df.master[,metrik])]
    return(list(Akurasi=df.master,n.optimum=paste("Alpha optimum adalah",opt)))
  }
  else if(method=="DES"){
    opt <- df.master$Alpha[which.min(df.master[,metrik])]
    opt2 <- df.master$Beta[which.min(df.master[,metrik])]
    return(list(Akurasi=df.master,n.optimum=paste("Alpha optimum adalah",opt,
                                                  ",Beta optimum adalah",opt2)))
  }
}

SAS

par(mar=c(1,1,1,1))
pemulusan(x=data$Kurs_Jual,min=2,max=10,method = "SMA",metrik = "MSE")

## $Akurasi
##    n     SSE       MSE      RMSE      MAD      MAPE
## 1  2 1571354  6186.431  78.65387 63.05423 0.3185231
## 2  3 1956408  7732.837  87.93655 69.81639 0.3528312
## 3  4 2331726  9252.882  96.19190 75.33272 0.3805974
## 4  5 2671985 10645.357 103.17634 80.28536 0.4054753
## 5  6 2992723 11970.892 109.41157 84.10787 0.4246626
## 6  7 3246460 13037.991 114.18402 87.28620 0.4403455
## 7  8 3455174 13932.153 118.03454 90.66235 0.4571334
## 8  9 3703719 14994.814 122.45331 94.56226 0.4766684
## 9 10 3967532 16128.178 126.99677 98.11385 0.4945390
## 
## $n.optimum
## [1] "n optimum adalah 2"

SES

pemulusan2(x=data$Kurs_Jual,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           method="SES",metrik="MSE")

## $Akurasi
##   Alpha     SSE       MSE      RMSE
## 1   0.1 5823358 22747.493 150.82272
## 2   0.2 3146003 12289.075 110.85610
## 3   0.3 2316206  9047.681  95.11930
## 4   0.4 1901653  7428.330  86.18776
## 5   0.5 1653383  6458.529  80.36497
## 6   0.6 1491548  5826.358  76.33058
## 7   0.7 1383478  5404.211  73.51334
## 8   0.8 1313546  5131.040  71.63128
## 9   0.9 1273686  4975.335  70.53605
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9"

DMA

par(mar=c(1,1,1,1))
pemulusan(x=data$Kurs_Jual,min=2,max=10,method = "DMA",metrik = "MSE")

## $Akurasi
##    n     SSE       MSE      RMSE       MAD      MAPE
## 1  2 1565904  6189.344  78.67239  62.89745 0.3173407
## 2  3 2140060  8526.136  92.33708  75.58943 0.3817482
## 3  4 2790207 11205.650 105.85674  83.83842 0.4233595
## 4  5 3435454 13908.720 117.93524  94.27952 0.4760195
## 5  6 3917089 15988.117 126.44413 101.74941 0.5133435
## 6  7 4310675 17739.402 133.18935 103.90173 0.5232927
## 7  8 4585508 19027.003 137.93840 108.43651 0.5456519
## 8  9 4882720 20429.791 142.93282 115.36217 0.5802503
## 9 10 5295199 22342.611 149.47445 122.07013 0.6142200
## 
## $n.optimum
## [1] "n optimum adalah 2"

DES

pemulusan2(x=data$Kurs_Jual,alfa=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),
           beta=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9),method="DES",metrik="MSE")

## $Akurasi
##    Alpha Beta      SSE       MSE      RMSE
## 1    0.1  0.1 10030614 39182.085 197.94465
## 2    0.1  0.2  7575060 29590.076 172.01766
## 3    0.1  0.3  8331141 32543.519 180.39822
## 4    0.1  0.4  9207969 35968.627 189.65397
## 5    0.1  0.5  9032041 35281.409 187.83346
## 6    0.1  0.6  8110353 31681.068 177.99176
## 7    0.1  0.7  7115765 27795.958 166.72120
## 8    0.1  0.8  6271333 24497.393 156.51643
## 9    0.1  0.9  5608390 21907.772 148.01274
## 10   0.2  0.1  4520596 17658.579 132.88559
## 11   0.2  0.2  4214841 16464.221 128.31298
## 12   0.2  0.3  4016562 15689.694 125.25851
## 13   0.2  0.4  3731997 14578.112 120.73985
## 14   0.2  0.5  3578138 13977.102 118.22479
## 15   0.2  0.6  3602363 14071.730 118.62432
## 16   0.2  0.7  3749067 14644.792 121.01567
## 17   0.2  0.8  3942928 15402.061 124.10504
## 18   0.2  0.9  4138958 16167.803 127.15268
## 19   0.3  0.1  3023630 11811.056 108.67868
## 20   0.3  0.2  2833785 11069.473 105.21156
## 21   0.3  0.3  2750940 10745.859 103.66223
## 22   0.3  0.4  2757961 10773.287 103.79444
## 23   0.3  0.5  2845066 11113.540 105.42078
## 24   0.3  0.6  2967313 11591.066 107.66181
## 25   0.3  0.7  3094939 12089.605 109.95274
## 26   0.3  0.8  3211953 12546.690 112.01201
## 27   0.3  0.9  3305363 12911.574 113.62911
## 28   0.4  0.1  2336714  9127.791  95.53947
## 29   0.4  0.2  2244999  8769.527  93.64575
## 30   0.4  0.3  2250494  8790.992  93.76029
## 31   0.4  0.4  2306974  9011.619  94.92955
## 32   0.4  0.5  2382546  9306.821  96.47186
## 33   0.4  0.6  2451600  9576.564  97.85992
## 34   0.4  0.7  2499875  9765.137  98.81871
## 35   0.4  0.8  2521743  9850.557  99.24997
## 36   0.4  0.9  2519386  9841.352  99.20359
## 37   0.5  0.1  1954271  7633.869  87.37202
## 38   0.5  0.2  1915223  7481.338  86.49473
## 39   0.5  0.3  1942530  7588.009  87.10918
## 40   0.5  0.4  1990099  7773.823  88.16929
## 41   0.5  0.5  2035040  7949.376  89.15927
## 42   0.5  0.6  2066151  8070.902  89.83820
## 43   0.5  0.7  2082379  8134.294  90.19032
## 44   0.5  0.8  2088894  8159.743  90.33130
## 45   0.5  0.9  2092873  8175.285  90.41728
## 46   0.6  0.1  1716218  6703.975  81.87781
## 47   0.6  0.2  1703786  6655.414  81.58072
## 48   0.6  0.3  1734906  6776.977  82.32240
## 49   0.6  0.4  1773371  6927.230  83.22998
## 50   0.6  0.5  1806687  7057.371  84.00816
## 51   0.6  0.6  1832691  7158.951  84.61058
## 52   0.6  0.7  1854890  7245.664  85.12146
## 53   0.6  0.8  1878045  7336.113  85.65111
## 54   0.6  0.9  1905670  7444.024  86.27876
## 55   0.7  0.1  1561913  6101.222  78.11032
## 56   0.7  0.2  1566137  6117.721  78.21586
## 57   0.7  0.3  1600856  6253.344  79.07809
## 58   0.7  0.4  1639730  6405.195  80.03246
## 59   0.7  0.5  1676814  6550.053  80.93240
## 60   0.7  0.6  1713065  6691.661  81.80258
## 61   0.7  0.7  1751643  6842.357  82.71854
## 62   0.7  0.8  1795229  7012.615  83.74136
## 63   0.7  0.9  1845315  7208.263  84.90149
## 64   0.8  0.1  1464018  5718.822  75.62289
## 65   0.8  0.2  1482373  5790.520  76.09547
## 66   0.8  0.3  1524564  5955.330  77.17078
## 67   0.8  0.4  1571604  6139.077  78.35226
## 68   0.8  0.5  1620831  6331.369  79.56990
## 69   0.8  0.6  1673834  6538.413  80.86045
## 70   0.8  0.7  1732933  6769.269  82.27557
## 71   0.8  0.8  1799857  7030.691  83.84921
## 72   0.8  0.9  1875691  7326.920  85.59743
## 73   0.9  0.1  1408769  5503.003  74.18223
## 74   0.9  0.2  1442046  5632.990  75.05325
## 75   0.9  0.3  1496387  5845.262  76.45431
## 76   0.9  0.4  1557983  6085.871  78.01199
## 77   0.9  0.5  1625956  6351.389  79.69560
## 78   0.9  0.6  1702214  6649.274  81.54308
## 79   0.9  0.7  1789074  6988.568  83.59766
## 80   0.9  0.8  1888697  7377.723  85.89367
## 81   0.9  0.9  2003336  7825.530  88.46202
## 
## $n.optimum
## [1] "Alpha optimum adalah 0.9 ,Beta optimum adalah 0.1"

Pemulusan Winter

winter.ts<-ts(data$Kurs_Jual, frequency = 23,)
training.ts<-ts(training$Kurs_Jual, frequency = 23)
testing.ts<-ts(testing$Kurs_Jual, start=233, frequency = 23)

Aditif

aditif <- HoltWinters(training.ts)
aditif
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = training.ts)
## 
## Smoothing parameters:
##  alpha: 0.914396
##  beta : 0.01882048
##  gamma: 1
## 
## Coefficients:
##             [,1]
## a   19179.034901
## b      -8.445389
## s1    -84.524442
## s2    -61.605144
## s3    -54.149268
## s4    -43.053158
## s5    -66.094801
## s6    -61.588248
## s7    -76.063215
## s8    -96.647790
## s9    -88.215760
## s10   -18.239487
## s11    14.819391
## s12    20.584515
## s13    74.776688
## s14    89.340018
## s15    95.684838
## s16   112.470237
## s17    99.193606
## s18    94.617583
## s19    74.384848
## s20    45.264356
## s21   -38.707249
## s22   -81.133289
## s23   -88.804901

Forcasting

ramalan1 <- forecast(aditif, h=23)
ramalan1
##          Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 11.08696       19086.07 18980.03 19192.10 18923.90 19248.23
## 11.13043       19100.54 18955.62 19245.46 18878.90 19322.18
## 11.17391       19099.55 18923.13 19275.97 18829.74 19369.36
## 11.21739       19102.20 18898.20 19306.20 18790.21 19414.19
## 11.26087       19070.71 18841.61 19299.81 18720.34 19421.09
## 11.30435       19066.77 18814.31 19319.24 18680.66 19452.89
## 11.34783       19043.85 18769.29 19318.42 18623.95 19463.76
## 11.39130       19014.82 18719.15 19310.50 18562.62 19467.03
## 11.43478       19014.81 18698.79 19330.83 18531.49 19498.13
## 11.47826       19076.34 18740.59 19412.09 18562.86 19589.82
## 11.52174       19100.96 18745.99 19455.92 18558.08 19643.83
## 11.56522       19098.27 18724.52 19472.03 18526.66 19669.89
## 11.60870       19144.02 18751.82 19536.22 18544.20 19743.84
## 11.65217       19150.14 18739.79 19560.49 18522.56 19777.72
## 11.69565       19148.04 18719.79 19576.29 18493.09 19802.98
## 11.73913       19156.38 18710.45 19602.31 18474.39 19838.37
## 11.78261       19134.66 18671.23 19598.09 18425.90 19843.41
## 11.82609       19121.64 18640.86 19602.41 18386.35 19856.92
## 11.86957       19092.96 18594.96 19590.95 18331.34 19854.57
## 11.91304       19055.39 18540.30 19570.49 18267.62 19843.16
## 11.95652       18962.97 18430.87 19495.08 18149.19 19776.76
## 12.00000       18912.10 18363.07 19461.13 18072.43 19751.77
## 12.04348       18895.99 18330.09 19461.88 18030.53 19761.44

Akurasi Training (SSE)

sse1.train <- aditif$SSE
sse1.train
## [1] 1436348

Multikatif

multi <- HoltWinters(training.ts)
multi
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = training.ts)
## 
## Smoothing parameters:
##  alpha: 0.914396
##  beta : 0.01882048
##  gamma: 1
## 
## Coefficients:
##             [,1]
## a   19179.034901
## b      -8.445389
## s1    -84.524442
## s2    -61.605144
## s3    -54.149268
## s4    -43.053158
## s5    -66.094801
## s6    -61.588248
## s7    -76.063215
## s8    -96.647790
## s9    -88.215760
## s10   -18.239487
## s11    14.819391
## s12    20.584515
## s13    74.776688
## s14    89.340018
## s15    95.684838
## s16   112.470237
## s17    99.193606
## s18    94.617583
## s19    74.384848
## s20    45.264356
## s21   -38.707249
## s22   -81.133289
## s23   -88.804901

Forcasting

ramalan2 <- forecast(multi, h=23)
ramalan2
##          Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 11.08696       19086.07 18980.03 19192.10 18923.90 19248.23
## 11.13043       19100.54 18955.62 19245.46 18878.90 19322.18
## 11.17391       19099.55 18923.13 19275.97 18829.74 19369.36
## 11.21739       19102.20 18898.20 19306.20 18790.21 19414.19
## 11.26087       19070.71 18841.61 19299.81 18720.34 19421.09
## 11.30435       19066.77 18814.31 19319.24 18680.66 19452.89
## 11.34783       19043.85 18769.29 19318.42 18623.95 19463.76
## 11.39130       19014.82 18719.15 19310.50 18562.62 19467.03
## 11.43478       19014.81 18698.79 19330.83 18531.49 19498.13
## 11.47826       19076.34 18740.59 19412.09 18562.86 19589.82
## 11.52174       19100.96 18745.99 19455.92 18558.08 19643.83
## 11.56522       19098.27 18724.52 19472.03 18526.66 19669.89
## 11.60870       19144.02 18751.82 19536.22 18544.20 19743.84
## 11.65217       19150.14 18739.79 19560.49 18522.56 19777.72
## 11.69565       19148.04 18719.79 19576.29 18493.09 19802.98
## 11.73913       19156.38 18710.45 19602.31 18474.39 19838.37
## 11.78261       19134.66 18671.23 19598.09 18425.90 19843.41
## 11.82609       19121.64 18640.86 19602.41 18386.35 19856.92
## 11.86957       19092.96 18594.96 19590.95 18331.34 19854.57
## 11.91304       19055.39 18540.30 19570.49 18267.62 19843.16
## 11.95652       18962.97 18430.87 19495.08 18149.19 19776.76
## 12.00000       18912.10 18363.07 19461.13 18072.43 19751.77
## 12.04348       18895.99 18330.09 19461.88 18030.53 19761.44

Akurasi Testing (SSE)

selisih1<-as.numeric(ramalan1$mean)-as.numeric(testing.ts)
## Warning in as.numeric(ramalan1$mean) - as.numeric(testing.ts): longer object
## length is not a multiple of shorter object length
selisih1
##  [1] -139.144930 -116.561022  -63.760536 -129.159815 -151.536846 -171.935683
##  [7] -156.676039  -89.886003  -34.769362    2.491522   32.065010   35.374745
## [13]   13.601529  -21.440531  -23.541099   67.338910   43.026890  102.965477
## [19]  -68.792646  -89.058528 -195.265521 -320.966951 -353.403952 -213.174930
SSEtesting1<-sum(selisih1^2)
selisih2<-as.numeric(ramalan2$mean)-as.numeric(testing.ts)
## Warning in as.numeric(ramalan2$mean) - as.numeric(testing.ts): longer object
## length is not a multiple of shorter object length
selisih2
##  [1] -139.144930 -116.561022  -63.760536 -129.159815 -151.536846 -171.935683
##  [7] -156.676039  -89.886003  -34.769362    2.491522   32.065010   35.374745
## [13]   13.601529  -21.440531  -23.541099   67.338910   43.026890  102.965477
## [19]  -68.792646  -89.058528 -195.265521 -320.966951 -353.403952 -213.174930
SSEtesting2<-sum(selisih2^2)
akurasi <- matrix(c(SSEtesting1, SSEtesting2), nrow=1, ncol=2)
row.names(akurasi)<- "SSE"
colnames(akurasi) <- c("Aditif", "Multiplikatif")

Simpulan

akurasi
##     Aditif Multiplikatif
## SSE 484679        484679

Pemodelan

Kestasioneran Data

Eksploratif

ACF dan PACF Plot

par(mfrow=c(1,1))
acf(training.ts)

pacf(training.ts)

ACF tails off slowly, sedangkan PACF cuts off setelah lag ke-1. Karena ACF tails off slowly atau polanya menurun secara perlahan, maka datanya tidak stasioner.

Uji Formal

ADF Test

Hipotesis:

H0: Data tidak stasioner

H1: Data stasioner

tseries::adf.test(training.ts)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  training.ts
## Dickey-Fuller = -1.4054, Lag order = 6, p-value = 0.8257
## alternative hypothesis: stationary

Bedasarkan hasil Augmented Dickey-Fuller Test (ADF Test) diperoleh p-value (0.8257) > alpha (0.05) , maka tak tolak H0. Artinya, tidak cukup bukti untuk menyatakan data stasioner pada taraf nyata 5% atau dengan kata lain data deret waktu tersebut tidak stasioner. Setelah melalui uji ADF dan menggunakan plot ACF & PACF data tidak stasioner sehingga perlu dilakukan differencing.

Differencing

data.dif1<-diff(training.ts, differences = 1)

Kestasionearan Data Setelah Differencing

plot(x = data.dif1,
     col = "black",
     lwd = 1,
     type = "o",
     main = "Setelah Differencing")

Setelah dilakukan differencing satu kali ( d = 1). Pola data kurs jual rupiah terhadap british pound sudah stasioner dilihat dari time series plot di atas.

Hipotesis:

H0: Data tidak stasioner

H1: Data stasioner

tseries::adf.test(data.dif1)
## Warning in tseries::adf.test(data.dif1): p-value smaller than printed p-value
## 
##  Augmented Dickey-Fuller Test
## 
## data:  data.dif1
## Dickey-Fuller = -7.1647, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationary

Setelah dilakukan differencing satu kali ( d = 1), dengan menggunakan uji ADF p-value (0.01) < alpha (0.05), maka tolak H0. Artinya, dengan menggunakan uji ADF data deret waktu tersebut stasioner pada taraf nyata 5%.

Identifikasi Model ARIMA

ACF

acf(data.dif1)

PACF

pacf(data.dif1)

EACF

eacf(data.dif1)
## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o o o o o o o o o o  o  o  o 
## 1 o o o o o o o o o o o  o  o  o 
## 2 o x o o o o o o o o o  o  o  o 
## 3 x x o o o o o o o o o  o  o  o 
## 4 o x x o o o o o o o o  o  o  o 
## 5 o x x o x o o o o o o  o  o  o 
## 6 x x x x o x o o o o o  o  o  o 
## 7 x x x x x x x o o o o  o  o  o

Berdasarkan plot ACF, PACF, dan EACF, diperoleh … model sebagai berikut: ARIMA(0,1,0), ARIMA (0,1,1), ARIMA (2,1,2), ARIMA (3,1,2)

Perbandingan Kebaikan Model

model1 <- arima(training.ts,order = c(0,1,0))
model2 <- arima(training.ts,order = c(0,1,1))
model3 <- arima(training.ts,order = c(2,1,2))
model4 <- arima(training.ts,order = c(3,1,2))
Model <- c("ARIMA (0,1,0)","ARIMA (0,1,1)","ARIMA (2,1,2)","ARIMA (3,1,2)")
AIC <- c(model1$aic,model2$aic,model3$aic,model4$aic)
Akurasi <- data.frame(Model,AIC)
kableExtra::kable(Akurasi)
Model AIC
ARIMA (0,1,0) 2626.367
ARIMA (0,1,1) 2628.310
ARIMA (2,1,2) 2634.295
ARIMA (3,1,2) 2632.150

Simpulan

Setelah dibandingkan, model ARIMA yang paling baik berdasarkan log likelihood dan ragam terkecil adalah ARIMA(0,1,0).

Data kurs jual rupiah terhadap british pound merupakan data tak stasioner. Oleh karena itu, perlu dilakukan penstasioneran terlebih dahulu yaitu dengan differencing sebelum menentukan orde ARIMA yang cocok. Berdasarkan analisis yang telah dilakukan, pendugaan model terbaik pada data kurs jual rupiah terhadap british pound adalah ARIMA (0, 1, 0).

Diagnostik Model

Eksploratif

residual <- model1$residuals
par(mfrow = c(2,2))
qqnorm(residual)
qqline(residual, col = "blue", lwd = 2)
plot(residual, type="o", 
     ylab = "Sisaan", xlab = "Order", main = "Sisaan Model1 vs Order")
abline(h = 0, col='red')
acf(residual)
pacf(residual)

  1. Normal Q-Q Plot

Berdasarkan hasil eksplorasi di atas, terlihat bahwa banyak amatan yang di sepanjang garis qq-plot distribusi normal. Sehingga secara eksploratif, dapat disimpulkan bahwa sisaan menyebar normal.

Residual (SisaanModel1) vs Order

Berdasarkan hasil eksplorasi di atas, titik pada plot kebebasan sisaan kebanyakan di sekitar titik nol. Namun, terdapat beberapa titik amatan yang terletak cukup jauh dari titik nol. Sehingga, belum dapat disimpulkan apakah terdapat autokorelasi atau tidak.

Plot ACF dan PACF

Berdasarkan hasil eksplorasi di atas, baik dari plot ACF maupun plot PACF, pada keduanya terdapat garis vertikal di lag tertentu yang melebihi tinggi garis biru horizontal. Artinya, menurut kedua plot ini, terdapat autokorelasi pada model.

Uji Formal

1. Sisaan Menyebar Normal

Hipotesis:

H0: Sisaan menyebar normal

H1: Sisaan tidak menyebar normal

shapiro.test(residual)
## 
##  Shapiro-Wilk normality test
## 
## data:  residual
## W = 0.99418, p-value = 0.5108

2. Sisaan Saling Bebas/Tidak Terdapat Autokorelasi

Hipotesis:

H0: Tidak terdapat autokorelasi

H1: Terdapat autokorelasi

Box.test(residual, type = "Ljung")
## 
##  Box-Ljung test
## 
## data:  residual
## X-squared = 0.047757, df = 1, p-value = 0.827

3. Nilai Tengah Sisaan Sama Dengan Nol

Hipotesis:

H0: Nilai tengah sisaan sama dengan nol

H1: Nilai tengah sisaan tidak sama dengan nol

t.test(residual, mu = 0, conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  residual
## t = -0.017531, df = 231, p-value = 0.986
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -9.296728  9.132749
## sample estimates:
##   mean of x 
## -0.08198978

Overfitting

1. Model ARIMA (0,1,0)

model1.a <- Arima(training.ts, order=c(1,1,0), method="ML")
summary(model1.a)
## Series: training.ts 
## ARIMA(1,1,0) 
## 
## Coefficients:
##          ar1
##       0.0158
## s.e.  0.0660
## 
## sigma^2 = 5095:  log likelihood = -1313.15
## AIC=2630.31   AICc=2630.36   BIC=2637.19
## 
## Training set error measures:
##                       ME     RMSE      MAE         MPE      MAPE      MASE
## Training set -0.08286392 71.07315 56.00379 -0.00106266 0.2819512 0.2310819
##                      ACF1
## Training set -0.001552789
lmtest::coeftest(model1.a)
## 
## z test of coefficients:
## 
##     Estimate Std. Error z value Pr(>|z|)
## ar1 0.015818   0.065950  0.2398   0.8105

2. Model ARIMA (0,1,1)

model2.a <- Arima(training.ts,order = c(1,1,1), method="ML")
summary(model2.a)
## Series: training.ts 
## ARIMA(1,1,1) 
## 
## Coefficients:
##          ar1     ma1
##       0.0080  0.0078
## s.e.  1.7322  1.7265
## 
## sigma^2 = 5118:  log likelihood = -1313.15
## AIC=2632.31   AICc=2632.42   BIC=2642.64
## 
## Training set error measures:
##                       ME     RMSE      MAE          MPE     MAPE      MASE
## Training set -0.08284219 71.07318 56.00392 -0.001062646 0.281952 0.2310824
##                      ACF1
## Training set -0.001531111
lmtest::coeftest(model2.a)
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)
## ar1 0.0079659  1.7321819  0.0046   0.9963
## ma1 0.0078251  1.7264693  0.0045   0.9964
model2.b <- Arima(training.ts,order = c(0,1,2), method="ML")
summary(model2.b)
## Series: training.ts 
## ARIMA(0,1,2) 
## 
## Coefficients:
##          ma1     ma2
##       0.0161  0.0038
## s.e.  0.0662  0.0695
## 
## sigma^2 = 5118:  log likelihood = -1313.15
## AIC=2632.31   AICc=2632.41   BIC=2642.63
## 
## Training set error measures:
##                       ME     RMSE      MAE          MPE     MAPE     MASE
## Training set -0.08342985 71.07274 56.00308 -0.001063169 0.281947 0.231079
##                      ACF1
## Training set -0.001664652
lmtest::coeftest(model2.b)
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)
## ma1 0.0160681  0.0661987  0.2427   0.8082
## ma2 0.0038052  0.0694669  0.0548   0.9563

3. Model ARIMA (2,1,2)

model3.a <- arima(training.ts,order = c(2,1,3), method = "ML")
## Warning in stats::arima(x = x, order = order, seasonal = seasonal, xreg =
## xreg, : possible convergence problem: optim gave code = 1
summary(model3.a)
## 
## Call:
## arima(x = training.ts, order = c(2, 1, 3), method = "ML")
## 
## Coefficients:
## Warning in sqrt(diag(x$var.coef)): NaNs produced
##          ar1      ar2      ma1     ma2      ma3
##       1.5692  -0.9494  -1.5817  0.9821  -0.0529
## s.e.     NaN      NaN      NaN     NaN   0.0701
## 
## sigma^2 estimated as 4958:  log likelihood = -1310.7,  aic = 2631.39
## 
## Training set error measures:
## Warning in trainingaccuracy(object, test, d, D): test elements must be within
## sample
##               ME RMSE MAE MPE MAPE
## Training set NaN  NaN NaN NaN  NaN
lmtest::coeftest(model3.a)
## Warning in sqrt(diag(se)): NaNs produced
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)
## ar1  1.569183        NaN     NaN      NaN
## ar2 -0.949354        NaN     NaN      NaN
## ma1 -1.581725        NaN     NaN      NaN
## ma2  0.982069        NaN     NaN      NaN
## ma3 -0.052892   0.070059  -0.755   0.4503

4. Model ARIMA (3,1,2)

model4.a <- arima(training.ts,order = c(3,1,3), method = "ML")
summary(model4.a)
## 
## Call:
## arima(x = training.ts, order = c(3, 1, 3), method = "ML")
## 
## Coefficients:
##          ar1      ar2     ar3      ma1     ma2      ma3
##       0.9088  -1.0590  0.7676  -0.8919  1.0348  -0.7876
## s.e.  0.6711   0.0785  0.6301   0.6253  0.0678   0.5918
## 
## sigma^2 estimated as 4966:  log likelihood = -1311.02,  aic = 2634.03
## 
## Training set error measures:
## Warning in trainingaccuracy(object, test, d, D): test elements must be within
## sample
##               ME RMSE MAE MPE MAPE
## Training set NaN  NaN NaN NaN  NaN
lmtest::coeftest(model4.a)
## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value Pr(>|z|)    
## ar1  0.908807   0.671088   1.3542   0.1757    
## ar2 -1.059041   0.078530 -13.4859   <2e-16 ***
## ar3  0.767587   0.630083   1.2182   0.2231    
## ma1 -0.891858   0.625269  -1.4264   0.1538    
## ma2  1.034787   0.067826  15.2565   <2e-16 ***
## ma3 -0.787567   0.591847  -1.3307   0.1833    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
model4.b <- arima(training.ts,order = c(4,1,2), method = "ML")
summary(model4.b)
## 
## Call:
## arima(x = training.ts, order = c(4, 1, 2), method = "ML")
## 
## Coefficients:
##          ar1     ar2      ar3      ar4      ma1      ma2
##       0.1055  0.2464  -0.0417  -0.0507  -0.0923  -0.2464
## s.e.  0.8786  0.5513   0.0671   0.0760   0.8780   0.5419
## 
## sigma^2 estimated as 5047:  log likelihood = -1312.6,  aic = 2637.21
## 
## Training set error measures:
## Warning in trainingaccuracy(object, test, d, D): test elements must be within
## sample
##               ME RMSE MAE MPE MAPE
## Training set NaN  NaN NaN NaN  NaN
lmtest::coeftest(model4.b)
## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)
## ar1  0.105508   0.878613  0.1201   0.9044
## ar2  0.246393   0.551328  0.4469   0.6549
## ar3 -0.041728   0.067065 -0.6222   0.5338
## ar4 -0.050742   0.075995 -0.6677   0.5043
## ma1 -0.092285   0.877976 -0.1051   0.9163
## ma2 -0.246411   0.541908 -0.4547   0.6493

Perbandingan Kebaikan Model dari Model Overfitting

Model <- c("ARIMA (0,1,0)", "ARIMA (1,1,0)","ARIMA (0,1,1)", "ARIMA (1,1,1)","ARIMA (0,1,2)","ARIMA (2,1,2)","ARIMA (2,1,3)","ARIMA (3,1,2)", "ARIMA (3,1,3)", "ARIMA (4,1,2)")
AIC <- c(model1$aic,model1.a$aic,model2$aic,model2.a$aic, model2.b$aic,model3$aic,model3.a$aic,model4$aic,model4.a$aic,model4.b$aic)
Akurasi <- data.frame(Model,AIC)
akurasioverfitting <- kableExtra::kable(Akurasi)

model.ov_aic <- data.frame(
  "Model" = c("ARIMA (0,1,0)", "ARIMA (1,1,0)","ARIMA (0,1,1)", "ARIMA (1,1,1)","ARIMA (0,1,2)","ARIMA (2,1,2)","ARIMA (2,1,3)","ARIMA (3,1,2)", "ARIMA (3,1,3)", "ARIMA (4,1,2)"),
  "AIC" = c(model1$aic,model1.a$aic,model2$aic,model2.a$aic, model2.b$aic,model3$aic,model3.a$aic,model4$aic,model4.a$aic,model4.b$aic)
)

model.ov_aic
##            Model      AIC
## 1  ARIMA (0,1,0) 2626.367
## 2  ARIMA (1,1,0) 2630.309
## 3  ARIMA (0,1,1) 2628.310
## 4  ARIMA (1,1,1) 2632.310
## 5  ARIMA (0,1,2) 2632.307
## 6  ARIMA (2,1,2) 2634.295
## 7  ARIMA (2,1,3) 2631.393
## 8  ARIMA (3,1,2) 2632.150
## 9  ARIMA (3,1,3) 2634.031
## 10 ARIMA (4,1,2) 2637.207
dplyr::arrange(.data=model.ov_aic, AIC)
##            Model      AIC
## 1  ARIMA (0,1,0) 2626.367
## 2  ARIMA (0,1,1) 2628.310
## 3  ARIMA (1,1,0) 2630.309
## 4  ARIMA (2,1,3) 2631.393
## 5  ARIMA (3,1,2) 2632.150
## 6  ARIMA (0,1,2) 2632.307
## 7  ARIMA (1,1,1) 2632.310
## 8  ARIMA (3,1,3) 2634.031
## 9  ARIMA (2,1,2) 2634.295
## 10 ARIMA (4,1,2) 2637.207

Setelah overfitting model yang paling baik ternyata tetap ARIMA (0,1,0)

Forecasting

1. Forecasting ARIMA (0,1,0)

ramalan1 <- forecast::forecast(training.ts,model=model1,h=24)
ramalan1
##          Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 11.08696       19090.23 18998.95 19181.51 18950.63 19229.83
## 11.13043       19090.23 18961.14 19219.32 18892.81 19287.65
## 11.17391       19090.23 18932.13 19248.33 18848.44 19332.02
## 11.21739       19090.23 18907.67 19272.79 18811.03 19369.43
## 11.26087       19090.23 18886.13 19294.33 18778.08 19402.38
## 11.30435       19090.23 18866.65 19313.81 18748.29 19432.17
## 11.34783       19090.23 18848.73 19331.73 18720.89 19459.57
## 11.39130       19090.23 18832.06 19348.40 18695.39 19485.07
## 11.43478       19090.23 18816.40 19364.06 18671.44 19509.02
## 11.47826       19090.23 18801.58 19378.88 18648.78 19531.68
## 11.52174       19090.23 18787.50 19392.96 18627.24 19553.22
## 11.56522       19090.23 18774.03 19406.43 18606.65 19573.81
## 11.60870       19090.23 18761.12 19419.34 18586.90 19593.56
## 11.65217       19090.23 18748.70 19431.76 18567.90 19612.56
## 11.69565       19090.23 18736.71 19443.75 18549.57 19630.89
## 11.73913       19090.23 18725.12 19455.34 18531.84 19648.62
## 11.78261       19090.23 18713.88 19466.58 18514.65 19665.81
## 11.82609       19090.23 18702.97 19477.49 18497.97 19682.49
## 11.86957       19090.23 18692.36 19488.10 18481.74 19698.72
## 11.91304       19090.23 18682.02 19498.44 18465.93 19714.53
## 11.95652       19090.23 18671.94 19508.52 18450.51 19729.95
## 12.00000       19090.23 18662.10 19518.36 18435.46 19745.00
## 12.04348       19090.23 18652.48 19527.98 18420.74 19759.72
## 12.08696       19090.23 18643.06 19537.40 18406.34 19774.12
data.ramalan1 <- ramalan1$mean
plot(ramalan1)
lines((length(training.ts)+1):(length(training.ts)+24), testing.ts,col="black")

2. Forecasting ARIMA (2,1,3)

ramalan2 <- forecast::forecast(training.ts,model=model3.a,h=24)
ramalan2
##          Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 11.08696       19098.35 19008.11 19188.59 18960.34 19236.36
## 11.13043       19112.41 18985.59 19239.24 18918.46 19306.37
## 11.17391       19129.38 18973.71 19285.06 18891.30 19367.47
## 11.21739       19142.66 18963.62 19321.71 18868.84 19416.49
## 11.26087       19147.39 18949.45 19345.33 18844.66 19450.12
## 11.30435       19142.20 18928.76 19355.64 18815.77 19468.63
## 11.34783       19129.57 18902.82 19356.31 18782.79 19476.35
## 11.39130       19114.67 18875.65 19353.70 18749.12 19480.23
## 11.43478       19103.30 18852.06 19354.53 18719.06 19487.53
## 11.47826       19099.58 18835.61 19363.55 18695.88 19503.28
## 11.52174       19104.55 18827.25 19381.85 18680.45 19528.65
## 11.56522       19115.88 18825.01 19406.75 18671.03 19560.73
## 11.60870       19128.94 18824.91 19432.96 18663.96 19593.91
## 11.65217       19138.67 18822.49 19454.85 18655.11 19622.23
## 11.69565       19141.55 18814.50 19468.60 18641.37 19641.73
## 11.73913       19136.83 18800.10 19473.55 18621.85 19651.80
## 11.78261       19126.68 18781.11 19472.26 18598.17 19655.19
## 11.82609       19115.25 18761.15 19469.35 18573.70 19656.79
## 11.86957       19106.93 18744.18 19469.69 18552.15 19661.72
## 11.91304       19104.74 18732.92 19476.57 18536.08 19673.40
## 11.95652       19109.20 18727.86 19490.55 18525.98 19692.42
## 12.00000       19118.28 18727.19 19509.36 18520.16 19716.39
## 12.04348       19128.28 18727.62 19528.94 18515.53 19741.03
## 12.08696       19135.37 18725.66 19545.07 18508.77 19761.96
data.ramalan2 <- ramalan2$mean
plot(ramalan2)
lines((length(training.ts)+1):(length(training.ts)+24), testing.ts,col="black")