Pertemuan 6 MPDW - Pendugaan Parameter, Diagnostik Model, dan Peramalan

library(ggplot2)

## Warning: package 'ggplot2' was built under R version 4.3.3

library(tsibble)

## Warning: package 'tsibble' was built under R version 4.3.3

## Registered S3 method overwritten by 'tsibble':
##   method               from 
##   as_tibble.grouped_df dplyr

## 
## Attaching package: 'tsibble'

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, union

library(tseries)

## Warning: package 'tseries' was built under R version 4.3.2

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

library(forecast)

## Warning: package 'forecast' was built under R version 4.3.3

library(MASS)

## Warning: package 'MASS' was built under R version 4.3.3

library(TSA)

## Warning: package 'TSA' was built under R version 4.3.3

## 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(TTR)

## Warning: package 'TTR' was built under R version 4.3.2

library(aTSA)

## Warning: package 'aTSA' was built under R version 4.3.2

## 
## Attaching package: 'aTSA'

## The following object is masked from 'package:forecast':
## 
##     forecast

## The following objects are masked from 'package:tseries':
## 
##     adf.test, kpss.test, pp.test

## The following object is masked from 'package:graphics':
## 
##     identify

library(graphics)

I. Data Keseluruhan

Data Time Series

data <- read.csv("C:\\Users\\user\\Downloads\\dt_pt2.csv")
data

##        Tanggal Penumpang.Berangkat
## 1   2023-04-26                7538
## 2   2023-04-27                6954
## 3   2023-04-28                7420
## 4   2023-04-29                8118
## 5   2023-04-30                8970
## 6   2023-05-01                9472
## 7   2023-05-02                7850
## 8   2023-05-03                7657
## 9   2023-05-04                7515
## 10  2023-05-05                7646
## 11  2023-05-06                7336
## 12  2023-05-07                7525
## 13  2023-05-08                6769
## 14  2023-05-09                6509
## 15  2023-05-10                7562
## 16  2023-05-11                6471
## 17  2023-05-12                6810
## 18  2023-05-13                6686
## 19  2023-05-14                6932
## 20  2023-05-15                6952
## 21  2023-05-16                6071
## 22  2023-05-17                8257
## 23  2023-05-18                7842
## 24  2023-05-19                5854
## 25  2023-05-20                5879
## 26  2023-05-21                7127
## 27  2023-05-22                6446
## 28  2023-05-23                6292
## 29  2023-05-24                6786
## 30  2023-05-25                6979
## 31  2023-05-26                7588
## 32  2023-05-27                6390
## 33  2023-05-28                7530
## 34  2023-05-29                6593
## 35  2023-05-30                7177
## 36  2023-05-31                8949
## 37  2023-06-01                8128
## 38  2023-06-02                5582
## 39  2023-06-03                5132
## 40  2023-06-04                6595
## 41  2023-06-05                6304
## 42  2023-06-06                5984
## 43  2023-06-07                6616
## 44  2023-06-08                6899
## 45  2023-06-09                6989
## 46  2023-06-10                6830
## 47  2023-06-11                7631
## 48  2023-06-12                7022
## 49  2023-06-13                7091
## 50  2023-06-14                8598
## 51  2023-06-15                9241
## 52  2023-06-16                8803
## 53  2023-06-17                9322
## 54  2023-06-18               10038
## 55  2023-06-19                9265
## 56  2023-06-20                9051
## 57  2023-06-21                9895
## 58  2023-06-22                9556
## 59  2023-06-23                9250
## 60  2023-06-24                9446
## 61  2023-06-25                9704
## 62  2023-06-26                8872
## 63  2023-06-27                3673
## 64  2023-06-28                8940
## 65  2023-06-29                5570
## 66  2023-06-30                6883
## 67  2023-07-01                6791
## 68  2023-07-02                7998
## 69  2023-07-03                7496
## 70  2023-07-04                7273
## 71  2023-07-05                8557
## 72  2023-07-06                8148
## 73  2023-07-07                7973
## 74  2023-07-08                7397
## 75  2023-07-09                7762
## 76  2023-07-10                7674
## 77  2023-07-11                7265
## 78  2023-07-12                8403
## 79  2023-07-13                7763
## 80  2023-07-14                8054
## 81  2023-07-15                7800
## 82  2023-07-16                7758
## 83  2023-07-17                6046
## 84  2023-07-18                6562
## 85  2023-07-19                7811
## 86  2023-07-20                7045
## 87  2023-07-21                7294
## 88  2023-07-22                7001
## 89  2023-07-23                7560
## 90  2023-07-24                7610
## 91  2023-07-25                6927
## 92  2023-07-26                7188
## 93  2023-07-27                6824
## 94  2023-07-28                7521
## 95  2023-07-29                6204
## 96  2023-07-30                5918
## 97  2023-07-31                6013
## 98  2023-08-01                6121
## 99  2023-08-02                7177
## 100 2023-08-03                6217
## 101 2023-08-04                6922
## 102 2023-08-05                6672
## 103 2023-08-06                6539
## 104 2023-08-07                6457
## 105 2023-08-08                5438
## 106 2023-08-09                7027
## 107 2023-08-10                7257
## 108 2023-08-11                7154
## 109 2023-08-12                7083
## 110 2023-08-13                7122
## 111 2023-08-14                6481
## 112 2023-08-15                5898
## 113 2023-08-16                6762
## 114 2023-08-17                6139
## 115 2023-08-18                5795

data$Tanggal <- as.Date(data$Tanggal, format ="%Y-%m-%d")
summary(data)

##     Tanggal           Penumpang.Berangkat
##  Min.   :2023-04-26   Min.   : 3673      
##  1st Qu.:2023-05-24   1st Qu.: 6594      
##  Median :2023-06-22   Median : 7154      
##  Mean   :2023-06-22   Mean   : 7298      
##  3rd Qu.:2023-07-20   3rd Qu.: 7806      
##  Max.   :2023-08-18   Max.   :10038

Plot Time Series

ggplot(data) + geom_line(aes(x=Tanggal, y=Penumpang.Berangkat), size=1)

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

data_ts <- ts(data$Penumpang.Berangkat)
plot(data_ts)

Plot deret waktu di atas menunjukkan bahwa data tidak stasioner dalam rataan, ditandai dengan adanya beberapa perubahan tren, terutama setelah titik sekitar lag 60, di mana terdapat lonjakan yang cukup besar diikuti oleh penurunan. Ini menandakan bahwa rataan data tidak konstan sepanjang waktu. Selain itu, lebar pita memiliki lebar yang sekilas hampir sama sehingga kemungkinan data stasioner dalam ragam.

Plot ACF

acf(data_ts)

Berdasarkan plot AFC tersebut, terlihat bahwa plot ACF tails off dan menurun secara perlahan seiring bertambahnya lag.

Uji ADF

tseries::adf.test(data_ts)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  data_ts
## Dickey-Fuller = -2.8291, Lag order = 4, p-value = 0.2324
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.2324 yang lebih besar dari taraf nyata 5% sehingga terima \(H_0\) dan menandakan bahwa data tidak stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

Plot Box-Cox

index <- seq(1:115)
bc = boxcox(data_ts~index, lambda = seq(-2,3,by=0.01))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda

## [1] 0.65

#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]

##   [1] -0.12 -0.11 -0.10 -0.09 -0.08 -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01
##  [13]  0.00  0.01  0.02  0.03  0.04  0.05  0.06  0.07  0.08  0.09  0.10  0.11
##  [25]  0.12  0.13  0.14  0.15  0.16  0.17  0.18  0.19  0.20  0.21  0.22  0.23
##  [37]  0.24  0.25  0.26  0.27  0.28  0.29  0.30  0.31  0.32  0.33  0.34  0.35
##  [49]  0.36  0.37  0.38  0.39  0.40  0.41  0.42  0.43  0.44  0.45  0.46  0.47
##  [61]  0.48  0.49  0.50  0.51  0.52  0.53  0.54  0.55  0.56  0.57  0.58  0.59
##  [73]  0.60  0.61  0.62  0.63  0.64  0.65  0.66  0.67  0.68  0.69  0.70  0.71
##  [85]  0.72  0.73  0.74  0.75  0.76  0.77  0.78  0.79  0.80  0.81  0.82  0.83
##  [97]  0.84  0.85  0.86  0.87  0.88  0.89  0.90  0.91  0.92  0.93  0.94  0.95
## [109]  0.96  0.97  0.98  0.99  1.00  1.01  1.02  1.03  1.04  1.05  1.06  1.07
## [121]  1.08  1.09  1.10  1.11  1.12  1.13  1.14  1.15  1.16  1.17  1.18  1.19
## [133]  1.20  1.21  1.22  1.23  1.24  1.25  1.26  1.27  1.28  1.29  1.30  1.31
## [145]  1.32  1.33  1.34  1.35  1.36  1.37  1.38  1.39  1.40  1.41  1.42  1.43
## [157]  1.44  1.45  1.46  1.47  1.48  1.49

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 0,65 dan pada selang kepercayaan 95% nilai memiliki batas bawah -0.12 dan batas atas 1.49. Selang tersebut memuat nilai satu sehingga dapat dikatakan bahwa data bangkitan stasioner dalam ragam.

II. Partisi Data - Data Awal Sebesar 50%

Data Time Series

dt_stas1 <- data_ts[1:58] |> ts()
head(dt_stas1)

## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
## [1] 7538 6954 7420 8118 8970 9472

mean(dt_stas1)

## [1] 7431.448

var(dt_stas1)

## [1] 1307058

Plot Time Series

dt_stas1 |> as_tsibble() |> 
  ggplot(aes(x = index, y = value)) +
  geom_line() + theme_bw() +
  xlab("Obs") + ylab("Nilai")

Plot deret waktu di atas menunjukkan bahwa data tidak stasioner dalam rataan, ditandai dengan adanya beberapa perubahan tren, terutama setelah titik sekitar lag 30 dan 40, di mana terdapat penurunan yang cukup besar diikuti oleh lonjakan. Ini menandakan bahwa rataan data tidak konstan sepanjang waktu. Selain itu, lebar pita memiliki lebar yang sekilas hampir sama sehingga kemungkinan data stasioner dalam ragam.

Plot ACF

acf(dt_stas1)

Berdasarkan plot ACF, terlihat bahwa plot ACF pada data tersebut cenderung tails off dan menurun secara perlahan seiring bertambahnya lag.

Uji ADF

tseries::adf.test(dt_stas1)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  dt_stas1
## Dickey-Fuller = -1.3173, Lag order = 3, p-value = 0.8503
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.8503 yang lebih besar dari taraf nyata 5% sehingga terima \(H_0\) dan menandakan bahwa data tidak stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

Plot Box-Cox

index <- seq(1:58)
bc = boxcox(dt_stas1~index, lambda = seq(-4,6,by=1))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda

## [1] -0.2626263

#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]

##  [1] -1.77777778 -1.67676768 -1.57575758 -1.47474747 -1.37373737 -1.27272727
##  [7] -1.17171717 -1.07070707 -0.96969697 -0.86868687 -0.76767677 -0.66666667
## [13] -0.56565657 -0.46464646 -0.36363636 -0.26262626 -0.16161616 -0.06060606
## [19]  0.04040404  0.14141414  0.24242424  0.34343434  0.44444444  0.54545455
## [25]  0.64646465  0.74747475  0.84848485  0.94949495  1.05050505  1.15151515
## [31]  1.25252525

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar -0.2626263 dan pada selang kepercayaan 95% nilai memiliki batas bawah -1.77777778 dan batas atas 1.25252525. Selang tersebut memuat nilai satu sehingga dapat dikatakan bahwa data bangkitan stasioner dalam ragam.

III. Partisi Data - Data Akhir Sebesar 50%

Data Time Series

dt_stas2 <- data_ts[59:115] |> ts()
head(dt_stas2)

## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
## [1] 9250 9446 9704 8872 3673 8940

mean(dt_stas2)

## [1] 7162.018

var(dt_stas2)

## [1] 1132565

Plot Time Series

dt_stas2 |> as_tsibble() |> 
  ggplot(aes(x = index, y = value)) +
  geom_line() + theme_bw() +
  xlab("Obs") + ylab("Nilai")

Plot deret waktu di atas menunjukkan bahwa data kemungkinan stasioner dalam rataan ditandai dengan data yang menyebar di sekitar nilai tengahnya (7162.018). Namun, plot juga mengindikasikan adanya beberapa perubahan tren. Ini menandakan bahwa rataan data mungkin saja tidak konstan sepanjang waktu. Selain itu, lebar pita memiliki lebar yang sekilas hampir sama sehingga kemungkinan data stasioner dalam ragam. Sehingga, perlu dilakukan analisis lanjut.

Plot ACF

acf(dt_stas2)

Berdasarkan plot ACF tersebut, terlihat bahwa data stasioner dalam rataan ditandai dengan plot ACF yang cuts off pada lag ke-2.

Uji ADF

tseries::adf.test(dt_stas2)

## Warning in tseries::adf.test(dt_stas2): p-value smaller than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  dt_stas2
## Dickey-Fuller = -4.8825, Lag order = 3, p-value = 0.01
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) dan menandakan bahwa data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF.

Plot Box-Cox

index <- seq(59:115)
bc = boxcox(dt_stas2~index, lambda = seq(0,4,by=1))

#Nilai Rounded Lambda
lambda <- bc$x[which.max(bc$y)]
lambda

## [1] 2.020202

#SK
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]

##  [1] 0.969697 1.010101 1.050505 1.090909 1.131313 1.171717 1.212121 1.252525
##  [9] 1.292929 1.333333 1.373737 1.414141 1.454545 1.494949 1.535354 1.575758
## [17] 1.616162 1.656566 1.696970 1.737374 1.777778 1.818182 1.858586 1.898990
## [25] 1.939394 1.979798 2.020202 2.060606 2.101010 2.141414 2.181818 2.222222
## [33] 2.262626 2.303030 2.343434 2.383838 2.424242 2.464646 2.505051 2.545455
## [41] 2.585859 2.626263 2.666667 2.707071 2.747475 2.787879 2.828283 2.868687
## [49] 2.909091 2.949495 2.989899 3.030303 3.070707 3.111111 3.151515

Gambar di atas menunjukkan nilai rounded value (\(\lambda\)) optimum sebesar 2.020202 dan pada selang kepercayaan 95% nilai memiliki batas bawah 0.969697 dan batas atas 3.151515. Selang tersebut memuat nilai satu sehingga dapat dikatakan bahwa data bangkitan stasioner dalam ragam.

IV. Penanganan Ketidakstasioneran Data

train.diff<-diff(dt_stas1,differences = 1) 
plot.ts(train.diff, lty=1, main="Plot Difference")

Berdasarkan plot data deret waktu, terlihat bahwa data sudah stasioner dalam rataan ditandai dengan data bergerak pada nilai tengah tertentu (tidak terdapat trend ataupun musiman pada data).

acf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot ACF cuts off pada lag ke 2. Hal ini menandakan data sudah stasioner dalam rataan dan ketidakstasioneran data telah berhasil tertangani.

tseries::adf.test(train.diff)

## Warning in tseries::adf.test(train.diff): p-value smaller than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  train.diff
## Dickey-Fuller = -5.0654, Lag order = 3, p-value = 0.01
## alternative hypothesis: stationary

\(H_0\) : Data tidak stasioner dalam rataan

\(H_1\) : Data stasioner dalam rataan

Berdasarkan uji ADF tersebut, didapat p-value sebesar 0.01 yang lebih kecil dari taraf nyata 5% sehingga tolak \(H_0\) atau data stasioner dalam rataan. Hal ini sesuai dengan hasil eksplorasi menggunakan plot time series dan plot ACF, sehingga dalam hal ini ketidakstasioneran data sudah berhasil ditangani dan dapat dilanjutkan ke pemodelan

V. Identifikasi Model

acf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot ACF cenderung cuts off pada lag ke 2, sehingga jika plot PACF dianggap tails of, maka model tentatifnya adalah ARIMA(0,1,2).

Plot PACF

pacf(train.diff)

Berdasarkan plot tersebut, terlihat bahwa plot PACF cenderung cuts off pada lag ke 2, sehingga jika plot ACF dianggap tails of, maka model tentatifnya adalah ARIMA(2,1,0).

Jika baik plot ACF maupun plot PACF keduanya dianggap tails off, maka model yang terbentuk adalah ARIMA(2,1,2)

Plot EACF

eacf(train.diff)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o x o o o o o o o o o  x  o  x 
## 1 x x o o o o o o o o o  x  o  x 
## 2 o o o o o o o o o o o  o  o  o 
## 3 x o o o o o o o o o o  o  o  o 
## 4 x o o o o o o o o o o  o  o  o 
## 5 x o o o o o o o o o o  o  o  o 
## 6 x o o o o o o o o o o  o  o  o 
## 7 x o o o o o o o o o o  o  o  o

Identifikasi model menggunakan plot EACF dilakukan dengan melihat ujung segitiga pada pola segitiga nol. Dalam hal ini model tentatif yang terbentuk adalah ARIMA(0,1,2), ARIMA(1,1,2), ARIMA(2,1,0), ARIMA(2,1,1), dan ARIMA(2,1,2).

VI. Pendugaan Parameter Model Tentatif

#ARIMA(0,1,2)
model1.da=Arima(train.diff, order=c(0,0,2),method="ML")
summary(model1.da) #AIC=929.41

## Series: train.diff 
## ARIMA(0,0,2) with non-zero mean 
## 
## Coefficients:
##           ma1      ma2     mean
##       -0.1695  -0.4076  30.6923
## s.e.   0.1352   0.1383  45.7179
## 
## sigma^2 = 642727:  log likelihood = -460.71
## AIC=929.41   AICc=930.18   BIC=937.59
## 
## Training set error measures:
##                    ME     RMSE    MAE      MPE     MAPE      MASE        ACF1
## Training set 3.765875 780.3197 632.82 58.06482 158.2342 0.5729933 -0.06310222

lmtest::coeftest(model1.da) #terdapat parameter yang tidak signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value Pr(>|z|)   
## ma1       -0.16951    0.13517 -1.2541 0.209812   
## ma2       -0.40756    0.13831 -2.9467 0.003212 **
## intercept 30.69229   45.71785  0.6713 0.502003   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#ARIMA(1,1,2)
model2.da=Arima(train.diff, order=c(1,0,2),method="ML")
summary(model2.da) #AIC=929.54

## Series: train.diff 
## ARIMA(1,0,2) with non-zero mean 
## 
## Coefficients:
##           ar1     ma1      ma2     mean
##       -0.4244  0.2099  -0.5325  34.0596
## s.e.   0.2204  0.1931   0.1107  49.7125
## 
## sigma^2 = 631075:  log likelihood = -459.77
## AIC=929.54   AICc=930.72   BIC=939.75
## 
## Training set error measures:
##                    ME     RMSE      MAE      MPE     MAPE      MASE
## Training set 1.765858 766.0215 617.7286 33.44867 193.6986 0.5593287
##                      ACF1
## Training set -0.007168734

lmtest::coeftest(model2.da) #terdapat parameter tidak signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value  Pr(>|z|)    
## ar1       -0.42439    0.22035 -1.9259   0.05411 .  
## ma1        0.20994    0.19312  1.0871   0.27700    
## ma2       -0.53254    0.11070 -4.8106 1.505e-06 ***
## intercept 34.05959   49.71254  0.6851   0.49326    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#ARIMA(2,1,0)
model3.da=Arima(train.diff, order=c(2,0,0),method="ML")
summary(model3.da) #AIC=927.92

## Series: train.diff 
## ARIMA(2,0,0) with non-zero mean 
## 
## Coefficients:
##           ar1      ar2     mean
##       -0.2124  -0.4482  36.4984
## s.e.   0.1179   0.1165  62.1634
## 
## sigma^2 = 625894:  log likelihood = -459.96
## AIC=927.92   AICc=928.69   BIC=936.1
## 
## Training set error measures:
##                     ME     RMSE      MAE      MPE     MAPE      MASE
## Training set -2.470455 770.0341 621.9286 58.38714 167.1116 0.5631317
##                     ACF1
## Training set -0.02837839

lmtest::coeftest(model3.da) #terdapat parameter yang tidak signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value Pr(>|z|)    
## ar1       -0.21240    0.11790 -1.8015 0.071630 .  
## ar2       -0.44821    0.11647 -3.8482 0.000119 ***
## intercept 36.49845   62.16337  0.5871 0.557111    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#model terpilih

#ARIMA(2,1,1)
model4.da=Arima(train.diff, order=c(2,0,1),method="ML")
summary(model4.da) #AIC=929.73

## Series: train.diff 
## ARIMA(2,0,1) with non-zero mean 
## 
## Coefficients:
##           ar1      ar2      ma1     mean
##       -0.1018  -0.4312  -0.1397  36.1618
## s.e.   0.2833   0.1294   0.3218  57.9669
## 
## sigma^2 = 635431:  log likelihood = -459.86
## AIC=929.73   AICc=930.91   BIC=939.95
## 
## Training set error measures:
##                     ME     RMSE      MAE      MPE    MAPE      MASE
## Training set -1.791302 768.6606 621.6785 61.63982 161.007 0.5629052
##                      ACF1
## Training set -0.002583302

lmtest::coeftest(model4.da) #terdapat parameter tidak signifikan

## 
## z test of coefficients:
## 
##           Estimate Std. Error z value  Pr(>|z|)    
## ar1       -0.10180    0.28328 -0.3593 0.7193352    
## ar2       -0.43117    0.12940 -3.3319 0.0008625 ***
## ma1       -0.13971    0.32183 -0.4341 0.6641995    
## intercept 36.16176   57.96694  0.6238 0.5327364    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#ARIMA(2,1,2)
model5.da=Arima(train.diff, order=c(2,0,2),method="ML")
summary(model5.da) #AIC=931.65

## Series: train.diff 
## ARIMA(2,0,2) with non-zero mean 
## 
## Coefficients:
##           ar1      ar2      ma1      ma2     mean
##       -0.1306  -0.3601  -0.1119  -0.0972  35.4339
## s.e.   0.3035   0.2986   0.3231   0.3600  54.9929
## 
## sigma^2 = 646606:  log likelihood = -459.82
## AIC=931.65   AICc=933.33   BIC=943.91
## 
## Training set error measures:
##                      ME     RMSE     MAE      MPE     MAPE     MASE
## Training set -0.1565411 768.0402 621.664 58.94574 166.1434 0.562892
##                      ACF1
## Training set -0.003933697

lmtest::coeftest(model5.da) #tidak ada parameter yang signifikan

## 
## z test of coefficients:
## 
##            Estimate Std. Error z value Pr(>|z|)
## ar1       -0.130621   0.303520 -0.4304   0.6669
## ar2       -0.360125   0.298649 -1.2058   0.2279
## ma1       -0.111903   0.323072 -0.3464   0.7291
## ma2       -0.097203   0.359968 -0.2700   0.7871
## intercept 35.433870  54.992949  0.6443   0.5194

Berdasarkan pendugaan parameter di atas, nilai AIC terkecil dimiliki oleh model ARIMA(0,1,2) namun ada parameter yang tidak signifikan. Namun, model ARIMA(2,1,0) memiliki nilai AIC ketiga terkecil dan parameter model ARIMA(2,1,0) juga seluruhnya lebih banyak signifikan daripada model lain jika intercept dihilangkan. Sehingga, model yang dipilih adalah model ARIMA(2,1,0).