lapply(c("dplyr","readxl","forecast","TTR","graphics","smooth","Mcomp","knitr","tseries",
"readr","TSA","ggplot2","MLmetrics","cowplot","gridExtra","gtable","grid","MASS",
"lmtest","GGally", "lubridate", "caret", "tseries", "tidyverse"), library, character.only = T)[[1]]
##
## 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
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## Loading required package: greybox
## Package "greybox", v1.0.6 loaded.
## This is package "smooth", v3.1.6
##
## Attaching package: 'smooth'
## The following object is masked from 'package:TTR':
##
## lags
## Registered S3 methods overwritten by 'TSA':
## method from
## fitted.Arima forecast
## plot.Arima forecast
##
## Attaching package: 'TSA'
## The following object is masked from 'package:readr':
##
## spec
## The following objects are masked from 'package:stats':
##
## acf, arima
## The following object is masked from 'package:utils':
##
## tar
##
## Attaching package: 'MLmetrics'
## The following objects are masked from 'package:greybox':
##
## MAE, MAPE, MSE
## The following object is masked from 'package:base':
##
## Recall
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
## Registered S3 method overwritten by 'GGally':
## method from
## +.gg ggplot2
## Warning: package 'lubridate' was built under R version 4.2.3
##
## Attaching package: 'lubridate'
## The following object is masked from 'package:cowplot':
##
## stamp
## The following object is masked from 'package:greybox':
##
## hm
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
## Loading required package: lattice
## Registered S3 method overwritten by 'lava':
## method from
## print.pcor greybox
##
## Attaching package: 'caret'
## The following objects are masked from 'package:MLmetrics':
##
## MAE, RMSE
## The following object is masked from 'package:greybox':
##
## MAE
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ tibble 3.1.8 ✔ stringr 1.5.0
## ✔ tidyr 1.2.1 ✔ forcats 0.5.2
## ✔ purrr 0.3.5
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ lubridate::as.difftime() masks base::as.difftime()
## ✖ gridExtra::combine() masks dplyr::combine()
## ✖ lubridate::date() masks base::date()
## ✖ dplyr::filter() masks stats::filter()
## ✖ lubridate::hm() masks greybox::hm()
## ✖ lubridate::intersect() masks base::intersect()
## ✖ dplyr::lag() masks stats::lag()
## ✖ purrr::lift() masks caret::lift()
## ✖ MASS::select() masks dplyr::select()
## ✖ lubridate::setdiff() masks base::setdiff()
## ✖ TSA::spec() masks readr::spec()
## ✖ tidyr::spread() masks greybox::spread()
## ✖ lubridate::stamp() masks cowplot::stamp()
## ✖ lubridate::union() masks base::union()
## [1] "dplyr" "stats" "graphics" "grDevices" "utils" "datasets"
## [7] "methods" "base"
#BELAWAN
# Set base path where the Excel files are located
belawan <- read_xlsx("C:/Users/User/Documents/Portofolio/Data Jumlah Penumpang Pelabuhan 2006-2023Belawan_2006_2023.xlsx")
belawan$`Total Keberangkatan` <- as.numeric(belawan$`Total Keberangkatan`)
belawan$`Total Kedatangan` <- as.numeric(belawan$`Total Kedatangan`)
str(belawan)
## tibble [213 × 5] (S3: tbl_df/tbl/data.frame)
## $ Tahun : num [1:213] 2006 2006 2006 2006 2006 ...
## $ Bulan : chr [1:213] "Januari" "Februari" "Maret" "April" ...
## $ Pelabuhan Utama : chr [1:213] "Belawan" "Belawan" "Belawan" "Belawan" ...
## $ Total Keberangkatan: num [1:213] 11610 5044 4483 4075 4404 ...
## $ Total Kedatangan : num [1:213] 7968 6707 4775 3469 3955 ...
# Memeriksa nilai NA, NaN, atau Inf
sum(is.na(belawan)) # Jumlah NA
## [1] 0
boxplot(belawan$`Total Keberangkatan` ~ belawan$Tahun)
data_belawan <- log(belawan$`Total Keberangkatan` + 1)
data.ts<-ts(data_belawan, frequency = 12, start = c(2006,1), end = c(2023,9))
# Box-Cox
#lambda <- BoxCox.lambda(data.ts)
# data.ts <- BoxCox(data.ts, lambda)
#Plot semua data
plot(data.ts,xlab ="Bulan", ylab = "Jumlah Penumpang", col="black", main = "Plot Deret Waktu Data")
points(data.ts)
seasonplot(data.ts,12,main="Jumlah Penumpang Pelabuhan Belawan", ylab="Air Passenger (miles)",year.labels = TRUE, col=rainbow(18))
fligner.test(data.ts ~ Tahun, data=belawan)
##
## Fligner-Killeen test of homogeneity of variances
##
## data: data.ts by Tahun
## Fligner-Killeen:med chi-squared = 26.324, df = 17, p-value = 0.06876
data.train <- ts(data.ts[1:170])
data.test <- ts(data.ts[171:213])
#Time Series Data
training.ts<-ts(data.train,frequency=12, start = c(2006,1), end = c(2020,2))
testing.ts<-ts(data.test,frequency=12, start = c(2020,3), end = c(2023,9))
#Plot Data
plot(training.ts, xlab ="Periode", ylab = "Jumlah Kedatangan", col="red", main = "Plot Data Training")
points(training.ts)
plot(testing.ts, xlab ="Periode", ylab = "Jumlah Kedatangan", col="red", main = "Plot Data Testing")
points(testing.ts)
## Plot Data Training dan testing Inflasi
ts.plot(data.ts, xlab = "Periode", ylab ="Jumlah Kedatangan",
main = "Plot Deret Waktu Data Jumlah kedatangan FKTP APBN")
lines(training.ts, col = "blue")
lines(testing.ts, col="Red")
legend("bottomleft",c("Data Training","Data Testing"),
lty=1, col=c("blue","red"), cex=0.8)
abline(v=2020.20, col=c("black"), lty=1, lwd=1)
acf(data.train, lag.max = 24, main = "Plot ACF")
pacf(data.train, lag.max = 24, main = "Plot PACF")
adf.test(data.train)
##
## Augmented Dickey-Fuller Test
##
## data: data.train
## Dickey-Fuller = -3.43, Lag order = 5, p-value = 0.05179
## alternative hypothesis: stationary
# DIFFERENCING
data.dif<-diff(data.train,differences = 1)
plot.ts(ts(data.dif,frequency=12, start = c(2015,1), end = c(2021,12)),lty=1,xlab = "Periode", ylab= "Data Inflasi Pembedaan 1", main="Plot Differencing Data Inflasi")
adf.test(data.dif)
## Warning in adf.test(data.dif): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: data.dif
## Dickey-Fuller = -10.141, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary
acf(data.dif, lag.max = 36, main = "Plot ACF Setelah Differencing 1 kali")
pacf(data.dif, lag.max = 48, main = "Plot PACF Setelah Differencing 1 kali")
eacf(data.dif)
## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x x o o o x o o o x o x o x
## 1 x x o o o x o o o x o x x x
## 2 x o o o o x o o o o o x x x
## 3 x o o o o x o o o o o x o x
## 4 x o o x x o o x o o o x x o
## 5 o x x x x o o x o x o x o x
## 6 x o x o x o x o o o o x o x
## 7 x x o o x x x x o o x o o o
ARIMA(0,1,2) ARIMA(0,1,3) ARIMA(1,1,2) ARIMA(1,1,3) ARIMA(2,1,1) ARIMA(2,1,2) ARIMA(2,1,3) ARIMA(3,1,1) ARIMA(3,1,2) ARIMA(3,1,3) ARIMA(5,1,0)
auto.arima(data.dif)
## Series: data.dif
## ARIMA(5,0,0) with zero mean
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5
## -0.6364 -0.7774 -0.6964 -0.5464 -0.4591
## s.e. 0.0684 0.0729 0.0778 0.0733 0.0684
##
## sigma^2 = 0.3776: log likelihood = -156.07
## AIC=324.14 AICc=324.66 BIC=342.92
model <- Arima(data.dif, order=c(0,0,2), method="ML")
model1 <- Arima(data.dif, order=c(0,0,3), method="ML")
model2 <- Arima(data.dif, order=c(1,0,2), method="ML")
model3 <- Arima(data.dif, order=c(1,0,3), method="ML")
model4 <- Arima(data.dif, order=c(2,0,1), method="ML")
model5 <- Arima(data.dif, order=c(2,0,2), method="ML")
model6 <- Arima(data.dif, order=c(2,0,3), method="ML")
model7 <- Arima(data.dif, order=c(3,0,1), method="ML")
model8 <- Arima(data.dif, order=c(3,0,2), method="ML")
model9 <- Arima(data.dif, order=c(3,0,3), method="ML")
model10 <- Arima(data.dif, order=c(5,0,0), method="ML")
summary(model)
## Series: data.dif
## ARIMA(0,0,2) with non-zero mean
##
## Coefficients:
## ma1 ma2 mean
## -0.6485 -0.3515 0.0023
## s.e. 0.0756 0.0731 0.0014
##
## sigma^2 = 0.4202: log likelihood = -167.35
## AIC=342.7 AICc=342.95 BIC=355.22
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.03184583 0.6424263 0.505148 Inf Inf 0.5316993 -0.02592382
coeftest(model)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.6485040 0.0756117 -8.5768 < 2.2e-16 ***
## ma2 -0.3514956 0.0730648 -4.8107 1.504e-06 ***
## intercept 0.0022704 0.0013539 1.6769 0.09356 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model1)
## Series: data.dif
## ARIMA(0,0,3) with non-zero mean
##
## Coefficients:
## ma1 ma2 ma3 mean
## -0.7169 -0.3799 0.0968 0.0023
## s.e. 0.0931 0.0732 0.0855 0.0012
##
## sigma^2 = 0.419: log likelihood = -166.73
## AIC=343.46 AICc=343.82 BIC=359.11
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.03396382 0.6395893 0.5025076 Inf Inf 0.5289201 0.02035327
coeftest(model1)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.7169437 0.0931369 -7.6977 1.385e-14 ***
## ma2 -0.3798852 0.0732192 -5.1883 2.122e-07 ***
## ma3 0.0968291 0.0855484 1.1319 0.25769
## intercept 0.0022860 0.0011878 1.9246 0.05428 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model2)
## Series: data.dif
## ARIMA(1,0,2) with non-zero mean
##
## Coefficients:
## ar1 ma1 ma2 mean
## -0.1112 -0.5626 -0.4374 0.0023
## s.e. 0.1575 0.1338 0.1324 0.0013
##
## sigma^2 = 0.4213: log likelihood = -167.1
## AIC=344.2 AICc=344.57 BIC=359.85
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.03255776 0.6413296 0.5042171 Inf Inf 0.5307195 -0.009005403
coeftest(model2)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.1112390 0.1574910 -0.7063 0.4799892
## ma1 -0.5625917 0.1338363 -4.2036 2.627e-05 ***
## ma2 -0.4373798 0.1323811 -3.3039 0.0009534 ***
## intercept 0.0022777 0.0012951 1.7587 0.0786328 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model3)
## Series: data.dif
## ARIMA(1,0,3) with non-zero mean
##
## Coefficients:
## ar1 ma1 ma2 ma3 mean
## 0.4612 -1.1859 -0.0678 0.3108 0.0047
## s.e. 0.1728 0.1586 0.1517 0.0663 0.0054
##
## sigma^2 = 0.4187: log likelihood = -164.71
## AIC=341.42 AICc=341.94 BIC=360.2
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.0120605 0.6374609 0.4870183 Inf Inf 0.5126167 0.03271562
coeftest(model3)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 0.4612348 0.1727609 2.6698 0.00759 **
## ma1 -1.1858570 0.1586157 -7.4763 7.645e-14 ***
## ma2 -0.0678148 0.1517438 -0.4469 0.65494
## ma3 0.3108217 0.0663499 4.6846 2.805e-06 ***
## intercept 0.0047366 0.0054354 0.8714 0.38351
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model4)
## Series: data.dif
## ARIMA(2,0,1) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 mean
## 0.2701 -0.3041 -0.9307 0.0037
## s.e. 0.0796 0.0776 0.0442 0.0037
##
## sigma^2 = 0.4108: log likelihood = -163.74
## AIC=337.49 AICc=337.86 BIC=353.14
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.01566578 0.6333281 0.48923 Inf Inf 0.5149447 -0.0349622
coeftest(model4)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 0.2701310 0.0796295 3.3923 0.000693 ***
## ar2 -0.3040776 0.0776356 -3.9167 8.976e-05 ***
## ma1 -0.9306814 0.0441885 -21.0616 < 2.2e-16 ***
## intercept 0.0037199 0.0036640 1.0153 0.309975
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model5)
## Series: data.dif
## ARIMA(2,0,2) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 ma2 mean
## 0.6370 -0.4282 -1.3327 0.4194 0.0045
## s.e. 0.1393 0.0784 0.1408 0.1385 0.0055
##
## sigma^2 = 0.4027: log likelihood = -161.47
## AIC=334.95 AICc=335.47 BIC=353.73
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.01157864 0.6251408 0.4734702 Inf Inf 0.4983565 -0.01152289
coeftest(model5)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 0.6370072 0.1392640 4.5741 4.783e-06 ***
## ar2 -0.4281611 0.0784342 -5.4589 4.792e-08 ***
## ma1 -1.3327314 0.1408451 -9.4624 < 2.2e-16 ***
## ma2 0.4194274 0.1384708 3.0290 0.002454 **
## intercept 0.0044537 0.0054683 0.8145 0.415381
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model6)
## Series: data.dif
## ARIMA(2,0,3) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 ma2 ma3 mean
## -0.9739 -0.5326 0.3228 -0.5697 -0.7531 0.0023
## s.e. 0.1451 0.1784 0.1514 0.1668 0.0969 0.0012
##
## sigma^2 = 0.4025: log likelihood = -162.44
## AIC=338.88 AICc=339.58 BIC=360.79
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.03287781 0.6230307 0.4909989 Inf Inf 0.5168065 0.002433312
coeftest(model6)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.9739431 0.1450689 -6.7137 1.898e-11 ***
## ar2 -0.5326489 0.1783520 -2.9865 0.0028219 **
## ma1 0.3227545 0.1514439 2.1312 0.0330742 *
## ma2 -0.5696944 0.1668352 -3.4147 0.0006385 ***
## ma3 -0.7530521 0.0969064 -7.7709 7.792e-15 ***
## intercept 0.0023372 0.0011966 1.9532 0.0507924 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model7)
## Series: data.dif
## ARIMA(3,0,1) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ma1 mean
## 0.1916 -0.2874 -0.1593 -0.8890 0.0042
## s.e. 0.0877 0.0760 0.0842 0.0514 0.0045
##
## sigma^2 = 0.4052: log likelihood = -162.02
## AIC=336.05 AICc=336.57 BIC=354.83
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.01333594 0.6270353 0.4772887 Inf Inf 0.5023756 -0.01562436
coeftest(model7)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 0.1915576 0.0876566 2.1853 0.0288655 *
## ar2 -0.2873985 0.0759666 -3.7832 0.0001548 ***
## ar3 -0.1593131 0.0841786 -1.8926 0.0584162 .
## ma1 -0.8889775 0.0514115 -17.2914 < 2.2e-16 ***
## intercept 0.0042030 0.0045205 0.9298 0.3524941
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model8)
## Series: data.dif
## ARIMA(3,0,2) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 mean
## -0.5889 -0.0732 -0.3350 -0.0738 -0.7772 0.0039
## s.e. 0.1299 0.1020 0.0769 0.1213 0.1173 0.0039
##
## sigma^2 = 0.4075: log likelihood = -162.05
## AIC=338.1 AICc=338.8 BIC=360.01
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.01475585 0.6269483 0.4825039 Inf Inf 0.507865 -0.02097627
coeftest(model8)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.5889447 0.1298596 -4.5352 5.754e-06 ***
## ar2 -0.0731575 0.1020368 -0.7170 0.4734
## ar3 -0.3349900 0.0769099 -4.3556 1.327e-05 ***
## ma1 -0.0738185 0.1212756 -0.6087 0.5427
## ma2 -0.7772287 0.1172670 -6.6279 3.406e-11 ***
## intercept 0.0038855 0.0039337 0.9878 0.3233
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model9)
## Series: data.dif
## ARIMA(3,0,3) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 ma3 mean
## -0.2959 0.1224 -0.4217 -0.4022 -0.7557 0.3314 0.0044
## s.e. 0.1646 0.1389 0.0816 0.1707 0.1362 0.1436 0.0054
##
## sigma^2 = 0.4036: log likelihood = -160.64
## AIC=337.29 AICc=338.19 BIC=362.33
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.01148596 0.6219606 0.4660784 Inf Inf 0.4905762 -0.001067397
coeftest(model9)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.2959278 0.1645786 -1.7981 0.07216 .
## ar2 0.1223986 0.1389416 0.8809 0.37835
## ar3 -0.4216525 0.0815765 -5.1688 2.356e-07 ***
## ma1 -0.4021696 0.1707087 -2.3559 0.01848 *
## ma2 -0.7557424 0.1362474 -5.5468 2.909e-08 ***
## ma3 0.3313898 0.1435848 2.3080 0.02100 *
## intercept 0.0044029 0.0054045 0.8147 0.41526
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model10)
## Series: data.dif
## ARIMA(5,0,0) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 mean
## -0.6369 -0.7780 -0.6969 -0.5467 -0.4594 0.0035
## s.e. 0.0684 0.0729 0.0778 0.0733 0.0683 0.0115
##
## sigma^2 = 0.3797: log likelihood = -156.02
## AIC=326.04 AICc=326.74 BIC=347.95
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.006282219 0.6051259 0.4519219 -Inf Inf 0.4756756 -0.03090491
coeftest(model10)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.6368975 0.0683826 -9.3137 < 2.2e-16 ***
## ar2 -0.7779711 0.0729434 -10.6654 < 2.2e-16 ***
## ar3 -0.6969115 0.0778069 -8.9569 < 2.2e-16 ***
## ar4 -0.5466870 0.0732836 -7.4599 8.660e-14 ***
## ar5 -0.4594229 0.0683392 -6.7227 1.784e-11 ***
## intercept 0.0035018 0.0114504 0.3058 0.7597
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Overfitting
model10b <- Arima(data.dif, order=c(6,0,0), method="ML")
model10c <- Arima(data.dif, order=c(5,0,1), method = "ML")
summary(model10b)
## Series: data.dif
## ARIMA(6,0,0) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 mean
## -0.6694 -0.8163 -0.7459 -0.6014 -0.5047 -0.0702 0.0035
## s.e. 0.0770 0.0841 0.0945 0.0947 0.0844 0.0772 0.0107
##
## sigma^2 = 0.3801: log likelihood = -155.61
## AIC=327.22 AICc=328.12 BIC=352.26
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.006370254 0.6035983 0.4528699 -Inf Inf 0.4766733 -0.00728627
coeftest(model10b)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.6693654 0.0770169 -8.6912 < 2.2e-16 ***
## ar2 -0.8163052 0.0841224 -9.7038 < 2.2e-16 ***
## ar3 -0.7458703 0.0945334 -7.8900 3.021e-15 ***
## ar4 -0.6014350 0.0947053 -6.3506 2.145e-10 ***
## ar5 -0.5047100 0.0843935 -5.9804 2.225e-09 ***
## ar6 -0.0701510 0.0772427 -0.9082 0.3638
## intercept 0.0034739 0.0106777 0.3253 0.7449
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model10c)
## Series: data.dif
## ARIMA(5,0,1) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ma1 mean
## -0.4788 -0.7003 -0.6044 -0.4783 -0.4081 -0.2019 0.0035
## s.e. 0.1521 0.0989 0.1125 0.0944 0.0888 0.1677 0.0102
##
## sigma^2 = 0.3792: log likelihood = -155.41
## AIC=326.83 AICc=327.73 BIC=351.87
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.006561576 0.6028646 0.4529101 -Inf Inf 0.4767157 0.003862243
coeftest(model10c)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.4787572 0.1520740 -3.1482 0.001643 **
## ar2 -0.7002507 0.0988572 -7.0835 1.406e-12 ***
## ar3 -0.6044196 0.1125189 -5.3717 7.799e-08 ***
## ar4 -0.4782802 0.0944448 -5.0641 4.103e-07 ***
## ar5 -0.4080629 0.0887664 -4.5970 4.285e-06 ***
## ma1 -0.2019196 0.1677205 -1.2039 0.228626
## intercept 0.0034989 0.0102289 0.3421 0.732302
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
tetap menggunakan model -> ARIMA(5,0,0)
#Diagnostik Sisaan
sisaan <- model10$residuals
#kenormalan
shapiro.test(sisaan)
##
## Shapiro-Wilk normality test
##
## data: sisaan
## W = 0.96676, p-value = 0.0004503
#sisaan saling bebas
Box.test(sisaan, type = "Ljung")
##
## Box-Ljung test
##
## data: sisaan
## X-squared = 0.1643, df = 1, p-value = 0.6852
# Forecasting
ramalan <- forecast(Arima(data.ts, order=c(5,1,0),method="ML",include.drift = TRUE),h = length(data.test))
data.ramalan <- ramalan$mean
data.ramalan.ts <- ts(data.ramalan, start = 2006, frequency = 12)
plot(ramalan,col="black",col.sub ="black",col.axis="black",
col.lab="black",col.main="black",lwd=2)
box(col="black",lwd=2)
hasilforecast<-matrix(data=c(data.ramalan[1:12]), nrow = 12, ncol = 1)
colnames(hasilforecast)<-c("Hasil Forecast")
head(hasilforecast)
## Hasil Forecast
## [1,] 8.168776
## [2,] 8.901080
## [3,] 8.986021
## [4,] 8.889284
## [5,] 8.912894
## [6,] 8.625313
hasilforecast <- matrix(data=exp(data.ramalan[1:15]), nrow=15, ncol=1)
colnames(hasilforecast) <- c("Hasil Forecast")
# Menampilkan baris pertama dari hasil forecast yang telah dikembalikan ke skala asli
hasilforecast
## Hasil Forecast
## [1,] 3529.023
## [2,] 7339.897
## [3,] 7990.602
## [4,] 7253.821
## [5,] 7427.125
## [6,] 5570.906
## [7,] 5054.168
## [8,] 5879.813
## [9,] 6602.481
## [10,] 6557.396
## [11,] 6511.253
## [12,] 6209.339
## [13,] 5923.222
## [14,] 5979.447
## [15,] 6207.540
Hasil forecast 15 bulan, sampai akhir 2024 -> awal 2025 dan selanjutnya plot forecast sudah lurus
# Contoh tanggal dalam format "dd/mm/yyyy"
start_date <- as.Date("01/10/2023", format="%d/%m/%Y")
# Membuat vektor tanggal untuk forecast dengan asumsi bahwa forecast dimulai pada bulan berikutnya
forecast_dates <- seq(from=start_date, length.out=15, by="month")
# Menambahkan vektor tanggal ke hasil forecast
hasilforecast <- matrix(data=exp(data.ramalan[1:15]), nrow=15, ncol=1)
colnames(hasilforecast) <- c("Hasil Forecast")
hasilforecast <- cbind(ForecastDate = forecast_dates, hasilforecast)
# Mengubah matrix menjadi data frame
hasilforecast_df <- as.data.frame(hasilforecast)
# Mengkonversi angka menjadi tanggal
hasilforecast_df$ForecastDate <- as.Date(hasilforecast_df$ForecastDate, origin = "1970-01-01")
# Menampilkan hasil forecast dengan tanggal yang benar
print(hasilforecast_df)
## ForecastDate Hasil Forecast
## 1 2023-10-01 3529.023
## 2 2023-11-01 7339.897
## 3 2023-12-01 7990.602
## 4 2024-01-01 7253.821
## 5 2024-02-01 7427.125
## 6 2024-03-01 5570.906
## 7 2024-04-01 5054.168
## 8 2024-05-01 5879.813
## 9 2024-06-01 6602.481
## 10 2024-07-01 6557.396
## 11 2024-08-01 6511.253
## 12 2024-09-01 6209.339
## 13 2024-10-01 5923.222
## 14 2024-11-01 5979.447
## 15 2024-12-01 6207.540
# Ubah data test menjadi dataframe
data_test_df <- data.frame(Actual = data.test)
# Ubah prediksi menjadi dataframe dan pastikan jumlah barisnya sama dengan data test
data_forecast_df <- data.frame(Forecast = head(data.ramalan, length(data.test)))
ramalan <- data.ramalan[1:43]
calculate_sMAPE <- function(actual, forecast) {
# Pastikan bahwa panjang vektor aktual dan forecast sama
if(length(actual) != length(forecast)) {
stop("The lengths of actual and forecast must be equal")
}
# Menghindari peringatan atau error ketika ada nilai nol
smape <- mean(2 * abs(forecast - actual) / (abs(actual) + abs(forecast) + .Machine$double.eps), na.rm = TRUE) * 100
return(smape)
}
# Validasi model
## Hitung error
error <- data.frame(data.test)-data.frame(data.ramalan[1:43])
## SSE (Sum Square Error)
SSE <- sum(error^2, na.rm = T)
## MSE (Mean Squared Error)
MSE<- sapply(error^2, mean, na.rm = T)
## RMSE (Root Mean Square Error)
RMSE1 <- sqrt(MSE)
# Menghitung rata-rata dari nilai aktual
mean_actual <- mean(data.test, na.rm = TRUE)
# Menghitung persentase RMSE
RMSE <- (RMSE1 / mean_actual) * 100
## MAD (Mean Absolute Deviation)
MAD <- sapply(abs(error), mean, na.rm = T)
## MAPE (Mean Absolute Percentage Error)
r.error <- (error/data.frame(data.test))*100 # Relative Error
#MAPE <- sapply(abs(r.error), mean, na.rm = T)
sMAPE <- calculate_sMAPE(data.test, ramalan)
akurasifarima <- data.frame(
"Ukuran Keakuratan" = c("SSE", "MSE", "sMAPE", "RMSE", "MAD"),
"Forecasting" = c(SSE, MSE, sMAPE, RMSE, MAD))
akurasifarima
## Ukuran.Keakuratan Forecasting
## 1 SSE 163.090388
## 2 MSE 3.792800
## 3 sMAPE 17.191577
## 4 RMSE 24.481428
## 5 MAD 1.192321
sMAPE 17% RMSE 24%