library(tseries)
## Warning: package 'tseries' was built under R version 4.4.3
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(timeSeries)
## Warning: package 'timeSeries' was built under R version 4.4.3
## Loading required package: timeDate
## Warning: package 'timeDate' was built under R version 4.4.3
## 
## Attaching package: 'timeSeries'
## The following objects are masked from 'package:graphics':
## 
##     lines, points
library(TSA)
## Warning: package 'TSA' was built under R version 4.4.3
## 
## Attaching package: 'TSA'
## The following objects are masked from 'package:timeDate':
## 
##     kurtosis, skewness
## The following objects are masked from 'package:stats':
## 
##     acf, arima
## The following object is masked from 'package:utils':
## 
##     tar
library(forecast)
## Warning: package 'forecast' was built under R version 4.4.3
## Registered S3 methods overwritten by 'forecast':
##   method       from
##   fitted.Arima TSA 
##   plot.Arima   TSA
set.seed(63)
# Panjang data
n <- 200
# Parameter ARIMA(p=1, d=1, q=1)
ar <- 0.5 # AR(1)
ma <- 0.2 # MA(1)
ts_arima <- arima.sim(model = list(order = c(1,1,1), ar = ar, ma = ma), n = n)
# Plot
ts.plot(ts_arima, main = "Simulasi Data ARIMA(1,1,1)")

acf(ts_arima)

pacf(ts_arima)

adf.test(ts_arima)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts_arima
## Dickey-Fuller = -2.9559, Lag order = 5, p-value = 0.1757
## alternative hypothesis: stationary

Melakukan differencing pada data ts_arima karena hasil ADF test menunjukkan bahwa data tidak stastioner.

diff1 <- diff(ts_arima)
acf(diff1)

pacf(diff1)

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

Karena hasil ADF test data differencing menunjukkan angka 0.01 < 0.05 maka data setelah differencing sudah stasioner dan siap untuk dilakukan pemodelan lanjut.

data.ts <- ts(diff1)
head(data.ts)
## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
## [1] -0.016718635  0.303778156  0.398897912 -0.006858533  0.165673724
## [6] -0.260238895
acf(data.ts)

pacf(data.ts)

eacf(data.ts)
## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x x o o o o o o o o o  o  o  o 
## 1 x o o o o o o o o o o  o  o  o 
## 2 x o x o o o o o o o o  o  o  o 
## 3 x x x o o o o o o o o  o  o  o 
## 4 x x x o o o o o o o o  o  o  o 
## 5 x o x x o o o o o o o  o  o  o 
## 6 x x x o o o o o o o o  o  o  o 
## 7 o x o o o o o o o o o  o  o  o
auto.arima(data.ts)
## Series: data.ts 
## ARIMA(1,0,1) with zero mean 
## 
## Coefficients:
##          ar1     ma1
##       0.4168  0.2705
## s.e.  0.1184  0.1347
## 
## sigma^2 = 0.9656:  log likelihood = -279.52
## AIC=565.04   AICc=565.16   BIC=574.93
auto.arima(ts_arima)
## Series: ts_arima 
## ARIMA(1,1,1) 
## 
## Coefficients:
##          ar1     ma1
##       0.4168  0.2705
## s.e.  0.1184  0.1347
## 
## sigma^2 = 0.9656:  log likelihood = -279.52
## AIC=565.04   AICc=565.16   BIC=574.93
arima(x = data.ts, order = c(0, 0, 0), method = "ML")
## 
## Call:
## arima(x = data.ts, order = c(0, 0, 0), method = "ML")
## 
## Coefficients:
##       intercept
##         -0.1131
## s.e.     0.0865
## 
## sigma^2 estimated as 1.496:  log likelihood = -324.04,  aic = 650.08
arima(x = data.ts, order = c(0, 0, 1), method = "ML")
## 
## Call:
## arima(x = data.ts, order = c(0, 0, 1), method = "ML")
## 
## Coefficients:
##          ma1  intercept
##       0.5968    -0.1114
## s.e.  0.0511     0.1128
## 
## sigma^2 estimated as 1.001:  log likelihood = -284.12,  aic = 572.24
arima(x = data.ts, order = c(1, 0, 0), method = "ML")
## 
## Call:
## arima(x = data.ts, order = c(1, 0, 0), method = "ML")
## 
## Coefficients:
##          ar1  intercept
##       0.5885    -0.1072
## s.e.  0.0568     0.1683
## 
## sigma^2 estimated as 0.9735:  log likelihood = -281.31,  aic = 566.62
arima(x = data.ts, order = c(2, 0, 0), method = "ML")
## 
## Call:
## arima(x = data.ts, order = c(2, 0, 0), method = "ML")
## 
## Coefficients:
##          ar1      ar2  intercept
##       0.6643  -0.1280    -0.1098
## s.e.  0.0699   0.0699     0.1486
## 
## sigma^2 estimated as 0.9573:  log likelihood = -279.65,  aic = 565.3
arima(x = data.ts, order = c(2, 0, 2), method = "ML")
## 
## Call:
## arima(x = data.ts, order = c(2, 0, 2), method = "ML")
## 
## Coefficients:
##          ar1     ar2     ma1     ma2  intercept
##       0.0119  0.2545  0.6714  0.0029    -0.1083
## s.e.  0.4269  0.1796  0.4253  0.1742     0.1560
## 
## sigma^2 estimated as 0.9452:  log likelihood = -278.41,  aic = 566.82

Berdasarkan hasil uji coba beberapa model, didapatkan kesimpulan bahwa model yang cocok untuk ‘data.ts’ adalah model ARIMA (2, 0, 0).

fit <- arima(data.ts, order = c(2, 0, 0), method = "ML")
Box.test(fit$residual,lag = 10, type = "Ljung-Box" )
## 
##  Box-Ljung test
## 
## data:  fit$residual
## X-squared = 6.6496, df = 10, p-value = 0.7581
ts.plot(fit$residuals, main="Residuals ARIMA(2,0,0)")
abline(h=0, col="red")

shapiro.test(fit$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  fit$residuals
## W = 0.99022, p-value = 0.1925

Setelah dilakukan uji Ljung-Box dan uji shapiro wilk, model ARIMA (2, 0, 0) cukup baik dalam menangkap pola data karena asumsi dasar residual terpenuh