# =========================
# LIBRARY
# =========================
invisible(lapply(c("readxl","dplyr","ggplot2","readr","moments","MASS","forecast","tseries","rugarch","FinTS","scales"), library, character.only = TRUE))
## Warning: package 'readxl' was built under R version 4.4.3
## Warning: package 'dplyr' was built under R version 4.4.3
##
## 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
## Warning: package 'ggplot2' was built under R version 4.4.3
## Warning: package 'readr' was built under R version 4.4.3
## Warning: package 'MASS' was built under R version 4.4.2
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
## Warning: package 'forecast' was built under R version 4.4.3
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## Warning: package 'tseries' was built under R version 4.4.3
## Warning: package 'rugarch' was built under R version 4.4.3
## Loading required package: parallel
## Warning: package 'FinTS' was built under R version 4.4.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.4.3
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Attaching package: 'FinTS'
## The following object is masked from 'package:forecast':
##
## Acf
## Warning: package 'scales' was built under R version 4.4.3
##
## Attaching package: 'scales'
## The following object is masked from 'package:readr':
##
## col_factor
# =========================
# 1. INPUT DATA, RETURN, DAN SINKRONISASI
# =========================
setwd("C:/Users/HP/Downloads/SKEMA B baru")
baca_data <- function(file) read_excel(file) %>% transmute(Date = as.Date(Tanggal), Close = as.numeric(Terakhir)) %>% filter(!is.na(Date), !is.na(Close)) %>% arrange(Date) %>% distinct(Date, .keep_all = TRUE)
ihsg <- baca_data("IHSG.xlsx"); sp500 <- baca_data("SP500.xlsx")
cat("Jumlah data IHSG :", nrow(ihsg), "\nJumlah data SP500 :", nrow(sp500), "\n")
## Jumlah data IHSG : 3763
## Jumlah data SP500 : 3890
ihsg_ret <- ihsg %>% mutate(Return_IHSG = log(Close / lag(Close))) %>% filter(!is.na(Return_IHSG)) %>% rename(Close_IHSG = Close)
sp500_ret <- sp500 %>% mutate(Return_SP500 = log(Close / lag(Close))) %>% filter(!is.na(Return_SP500)) %>% rename(Close_SP500 = Close)
idx_sp <- findInterval(as.numeric(ihsg_ret$Date) - 1e-9, as.numeric(sp500_ret$Date)); idx_sp[idx_sp == 0] <- NA_integer_
data_model <- ihsg_ret %>% mutate(Date_SP500_Exog = sp500_ret$Date[idx_sp], Close_SP500_Exog = sp500_ret$Close_SP500[idx_sp], Return_SP500_Exog = sp500_ret$Return_SP500[idx_sp], Lag_Days_SP500 = as.integer(Date - Date_SP500_Exog)) %>% filter(!is.na(Return_IHSG), !is.na(Return_SP500_Exog)) %>% dplyr::select(Date, Close_IHSG, Return_IHSG, Date_SP500_Exog, Close_SP500_Exog, Return_SP500_Exog, Lag_Days_SP500)
str(data_model); head(data_model); tail(data_model); print(table(data_model$Lag_Days_SP500))
## tibble [3,761 × 7] (S3: tbl_df/tbl/data.frame)
## $ Date : Date[1:3761], format: "2011-01-05" "2011-01-06" ...
## $ Close_IHSG : num [1:3761] 3784 3736 3631 3479 3455 ...
## $ Return_IHSG : num [1:3761] 0.00627 -0.01262 -0.02845 -0.04302 -0.00676 ...
## $ Date_SP500_Exog : Date[1:3761], format: "2011-01-04" "2011-01-05" ...
## $ Close_SP500_Exog : num [1:3761] 1270 1277 1274 1272 1270 ...
## $ Return_SP500_Exog: num [1:3761] -0.00134 0.00503 -0.0022 -0.00181 -0.00134 ...
## $ Lag_Days_SP500 : int [1:3761] 1 1 1 3 1 1 1 1 3 4 ...
## # A tibble: 6 × 7
## Date Close_IHSG Return_IHSG Date_SP500_Exog Close_SP500_Exog
## <date> <dbl> <dbl> <date> <dbl>
## 1 2011-01-05 3784. 0.00627 2011-01-04 1270.
## 2 2011-01-06 3736. -0.0126 2011-01-05 1277.
## 3 2011-01-07 3631. -0.0285 2011-01-06 1274.
## 4 2011-01-10 3479. -0.0430 2011-01-07 1272.
## 5 2011-01-11 3455. -0.00676 2011-01-10 1270.
## 6 2011-01-12 3555. 0.0284 2011-01-11 1274.
## # ℹ 2 more variables: Return_SP500_Exog <dbl>, Lag_Days_SP500 <int>
## # A tibble: 6 × 7
## Date Close_IHSG Return_IHSG Date_SP500_Exog Close_SP500_Exog
## <date> <dbl> <dbl> <date> <dbl>
## 1 2026-06-17 6221. -0.00549 2026-06-16 7511.
## 2 2026-06-18 6172. -0.00781 2026-06-17 7420.
## 3 2026-06-19 6177. 0.000777 2026-06-18 7501.
## 4 2026-06-22 6117. -0.00983 2026-06-18 7501.
## 5 2026-06-23 6101. -0.00251 2026-06-22 7473.
## 6 2026-06-24 5888. -0.0356 2026-06-23 7365.
## # ℹ 2 more variables: Return_SP500_Exog <dbl>, Lag_Days_SP500 <int>
##
## 1 2 3 4 5
## 2890 42 728 100 1
# =========================
# 2. FUNGSI TEMA DAN PLOT DASAR
# =========================
tema_ts <- function() theme_minimal(base_family = "serif") + theme(axis.text.x = element_text(angle = 45, hjust = 1, color = "black", size = 12), axis.text.y = element_text(color = "black", size = 12), axis.title = element_text(color = "black", size = 16), legend.title = element_text(color = "black", size = 12), legend.text = element_text(color = "black", size = 11), axis.line = element_line(color = "black", linewidth = 0.5), panel.grid.major = element_line(color = "grey60", linewidth = 0.4), panel.grid.minor = element_line(color = "grey75", linewidth = 0.4))
buat_periode <- function(b) {b <- as.Date(b); if(length(b) == 0) return(data.frame(xmin = as.Date(character()), xmax = as.Date(character()))); if(length(b) %% 2 != 0) stop("Jumlah tanggal batas_koreksi harus genap."); data.frame(xmin = b[c(TRUE,FALSE)], xmax = b[c(FALSE,TRUE)])}
plot_ts <- function(data, y_col, y_lab, batas_koreksi = as.Date(character()), warna = "blue", is_return = FALSE) {periode <- buat_periode(batas_koreksi); p <- ggplot(data, aes(Date, .data[[y_col]])); if(nrow(periode) > 0) p <- p + geom_rect(data = periode, aes(xmin = xmin, xmax = xmax, ymin = -Inf, ymax = Inf), inherit.aes = FALSE, fill = "red", alpha = 0.30); if(is_return) p <- p + geom_hline(yintercept = 0, color = "black", linewidth = 0.5); p <- p + geom_line(color = warna, linewidth = 0.8); if(length(batas_koreksi) > 0) p <- p + geom_vline(xintercept = batas_koreksi, linetype = "dashed", color = "red", linewidth = 0.8, alpha = 0.7); p + scale_x_date(date_breaks = "1 year", date_labels = "%Y", expand = expansion(mult = c(0.01,0.01))) + scale_y_continuous(labels = if(is_return) percent_format(accuracy = 1) else label_number(big.mark = ".", decimal.mark = ",")) + labs(x = "Date", y = y_lab) + tema_ts()}
batas_koreksi_ihsg <- as.Date(character()); batas_koreksi_sp500 <- as.Date(character())
plot.1 <- plot_ts(ihsg, "Close", "Closing Price IHSG", batas_koreksi_ihsg, "blue")
plot.2 <- plot_ts(sp500, "Close", "Closing Price S&P 500", batas_koreksi_sp500, "blue")
plot.3 <- plot_ts(data_model, "Return_IHSG", "Return IHSG", batas_koreksi_ihsg, "red", TRUE)
plot.4 <- plot_ts(data_model, "Return_SP500_Exog", "Synchronized Return S&P 500", batas_koreksi_sp500, "red", TRUE)
plot.1; plot.2; plot.3; plot.4
# =========================
# 3. DESKRIPTIF, SPLIT DATA, DAN ADF
# =========================
stat_deskriptif <- function(x) {x <- na.omit(x); data.frame(Observations = length(x), Mean = mean(x), Std_Dev = sd(x), Skewness = skewness(x), Kurtosis = kurtosis(x))}
deskriptif_return <- bind_rows(stat_deskriptif(data_model$Return_IHSG) %>% mutate(Variable = "Return IHSG"), stat_deskriptif(data_model$Return_SP500_Exog) %>% mutate(Variable = "Return S&P 500")) %>% dplyr::select(Variable, everything())
deskriptif_return
## Variable Observations Mean Std_Dev Skewness Kurtosis
## 1 Return IHSG 3761 0.0001192317 0.01077818 -0.6297416 11.45148
## 2 Return S&P 500 3761 0.0005127508 0.01067056 -0.7731835 17.31353
n <- nrow(data_model); n_train <- floor(0.90 * n); data_train <- data_model[1:n_train, ]; data_test <- data_model[(n_train + 1):n, ]
cat("Jumlah data latih =", nrow(data_train), "\nJumlah data uji =", nrow(data_test), "\n"); print(range(data_train$Date)); print(range(data_test$Date))
## Jumlah data latih = 3384
## Jumlah data uji = 377
## [1] "2011-01-05" "2024-11-08"
## [1] "2024-11-11" "2026-06-24"
tanggal_awal_uji <- min(data_test$Date)
plot.train.test <- ggplot(ihsg %>% mutate(Jenis_Data = ifelse(Date < tanggal_awal_uji, "Data Latih", "Data Uji")), aes(Date, Close)) + geom_line(aes(color = Jenis_Data), linewidth = 0.8) + geom_vline(xintercept = tanggal_awal_uji, linetype = "dashed", color = "black", linewidth = 0.8) + scale_color_manual(values = c("Data Latih" = "blue", "Data Uji" = "red")) + scale_x_date(date_breaks = "1 year", date_labels = "%Y") + scale_y_continuous(labels = label_number(big.mark = ".", decimal.mark = ",")) + labs(x = "Tahun", y = "Harga Penutupan IHSG", color = "Keterangan") + tema_ts()
plot.train.test
adf_rapi <- function(x, nama, lag_order = 15) {h <- suppressWarnings(adf.test(na.omit(x), alternative = "stationary", k = lag_order)); data.frame(Variable = nama, Test = "ADF", Test_Statistic = as.numeric(h$statistic), Lag_Order = lag_order, P_Value = h$p.value, Decision = ifelse(h$p.value < 0.05, "Stationary", "Not Stationary"))}
lag_adf <- 15
adf_train <- bind_rows(adf_rapi(data_train$Return_IHSG, "IHSG Return", lag_adf), adf_rapi(data_train$Return_SP500_Exog, "S&P 500 Return", lag_adf))
adf_train
## Variable Test Test_Statistic Lag_Order P_Value Decision
## 1 IHSG Return ADF -15.28783 15 0.01 Stationary
## 2 S&P 500 Return ADF -14.47395 15 0.01 Stationary
# =========================
# 4. PREWHITENING DAN CCF
# =========================
acf(data_train$Return_SP500_Exog, lag.max = 20, xlim = c(1,20), ylim = c(-0.15, 0.15))
pacf(data_train$Return_SP500_Exog, lag.max = 20)
spec_sp <- list(model_1 = list(order = c(4,0,2), fixed = c(NA,NA,0,NA,NA,NA,NA)), model_2 = list(order = c(4,0,0), fixed = c(NA,NA,0,NA,NA)), model_3 = list(order = c(0,0,2), fixed = c(NA,NA,NA)), model_4 = list(order = c(2,0,2), fixed = c(NA,NA,NA,NA,NA)), model_5 = list(order = c(1,0,2), fixed = c(NA,NA,NA,NA)), model_6 = list(order = c(2,0,1), fixed = c(NA,NA,NA,NA)), model_7 = list(order = c(1,0,1), fixed = c(NA,NA,NA)))
model_sp <- lapply(spec_sp, function(s) stats::arima(data_train$Return_SP500_Exog, order = s$order, fixed = s$fixed, method = "ML", include.mean = TRUE))
## Warning in stats::arima(data_train$Return_SP500_Exog, order = s$order, fixed =
## s$fixed, : some AR parameters were fixed: setting transform.pars = FALSE
## Warning in stats::arima(data_train$Return_SP500_Exog, order = s$order, fixed =
## s$fixed, : some AR parameters were fixed: setting transform.pars = FALSE
rekap_sp500 <- data.frame(Model = names(model_sp), AIC = sapply(model_sp, AIC), BIC = sapply(model_sp, BIC), Ljung_Box_pvalue = sapply(model_sp, function(m) Box.test(residuals(m), lag = 10, type = "Ljung-Box")$p.value))
rekap_sp500
## Model AIC BIC Ljung_Box_pvalue
## model_1 model_1 -21114.89 -21072.00 0.24014224
## model_2 model_2 -21118.45 -21087.81 0.21283411
## model_3 model_3 -21114.64 -21090.14 0.04416628
## model_4 model_4 -21113.77 -21077.01 0.11370100
## model_5 model_5 -21113.96 -21083.32 0.07863368
## model_6 model_6 -21113.04 -21082.41 0.06284383
## model_7 model_7 -21111.35 -21086.84 0.02316958
model_sp500_terbaik <- forecast::Arima(data_train$Return_SP500_Exog, order = c(4,0,0), fixed = c(NA,NA,0,NA,NA), include.mean = TRUE, method = "ML")
model_sp500_terbaik
## Series: data_train$Return_SP500_Exog
## ARIMA(4,0,0) with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 ar4 mean
## -0.0995 0.0662 0 -0.0400 5e-04
## s.e. 0.0171 0.0172 0 0.0171 2e-04
##
## sigma^2 = 0.0001139: log likelihood = 10564.22
## AIC=-21118.45 AICc=-21118.43 BIC=-21087.81
res_sp500 <- residuals(model_sp500_terbaik); Box.test(res_sp500, lag = 10, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: res_sp500
## X-squared = 13.198, df = 10, p-value = 0.2128
sp500_pw <- residuals(model_sp500_terbaik)
ihsg_pw <- residuals(forecast::Arima(data_train$Return_IHSG, model = model_sp500_terbaik))
prewhitened_data <- data.frame(Date = data_train$Date, SP500_Prewhitened = as.numeric(sp500_pw), IHSG_Prewhitened = as.numeric(ihsg_pw)) %>% filter(!is.na(SP500_Prewhitened), !is.na(IHSG_Prewhitened))
ccf_pw <- ccf(prewhitened_data$SP500_Prewhitened, prewhitened_data$IHSG_Prewhitened, lag.max = 20, main = "CCF Prewhitening: S&P 500 terhadap IHSG", ylab = "Cross Correlation")
ccf_table <- data.frame(Lag = as.numeric(ccf_pw$lag), CCF = as.numeric(ccf_pw$acf)) %>% mutate(Batas_Atas = 1.96 / sqrt(nrow(prewhitened_data)), Batas_Bawah = -1.96 / sqrt(nrow(prewhitened_data)), Signifikan = ifelse(abs(CCF) > Batas_Atas, "Signifikan", "Tidak Signifikan")) %>% arrange(Lag)
ccf_table
## Lag CCF Batas_Atas Batas_Bawah Signifikan
## 1 -20 8.188814e-03 0.03369309 -0.03369309 Tidak Signifikan
## 2 -19 4.799343e-03 0.03369309 -0.03369309 Tidak Signifikan
## 3 -18 -7.300100e-03 0.03369309 -0.03369309 Tidak Signifikan
## 4 -17 1.277615e-02 0.03369309 -0.03369309 Tidak Signifikan
## 5 -16 -9.298898e-03 0.03369309 -0.03369309 Tidak Signifikan
## 6 -15 1.331959e-02 0.03369309 -0.03369309 Tidak Signifikan
## 7 -14 1.706262e-02 0.03369309 -0.03369309 Tidak Signifikan
## 8 -13 -1.713465e-02 0.03369309 -0.03369309 Tidak Signifikan
## 9 -12 -4.867493e-03 0.03369309 -0.03369309 Tidak Signifikan
## 10 -11 1.807978e-02 0.03369309 -0.03369309 Tidak Signifikan
## 11 -10 4.375356e-03 0.03369309 -0.03369309 Tidak Signifikan
## 12 -9 -3.662920e-02 0.03369309 -0.03369309 Signifikan
## 13 -8 3.889919e-03 0.03369309 -0.03369309 Tidak Signifikan
## 14 -7 4.655443e-02 0.03369309 -0.03369309 Signifikan
## 15 -6 -7.939770e-03 0.03369309 -0.03369309 Tidak Signifikan
## 16 -5 -3.429580e-02 0.03369309 -0.03369309 Signifikan
## 17 -4 1.748013e-02 0.03369309 -0.03369309 Tidak Signifikan
## 18 -3 1.071124e-03 0.03369309 -0.03369309 Tidak Signifikan
## 19 -2 3.201359e-02 0.03369309 -0.03369309 Tidak Signifikan
## 20 -1 -2.242445e-02 0.03369309 -0.03369309 Tidak Signifikan
## 21 0 3.041923e-01 0.03369309 -0.03369309 Signifikan
## 22 1 2.174147e-01 0.03369309 -0.03369309 Signifikan
## 23 2 -1.419924e-02 0.03369309 -0.03369309 Tidak Signifikan
## 24 3 -3.631938e-03 0.03369309 -0.03369309 Tidak Signifikan
## 25 4 1.691475e-02 0.03369309 -0.03369309 Tidak Signifikan
## 26 5 -1.818233e-02 0.03369309 -0.03369309 Tidak Signifikan
## 27 6 3.087976e-02 0.03369309 -0.03369309 Tidak Signifikan
## 28 7 3.409496e-02 0.03369309 -0.03369309 Signifikan
## 29 8 3.938027e-02 0.03369309 -0.03369309 Signifikan
## 30 9 -5.862090e-03 0.03369309 -0.03369309 Tidak Signifikan
## 31 10 -5.537485e-03 0.03369309 -0.03369309 Tidak Signifikan
## 32 11 2.340867e-02 0.03369309 -0.03369309 Tidak Signifikan
## 33 12 -3.759728e-02 0.03369309 -0.03369309 Signifikan
## 34 13 -6.017329e-02 0.03369309 -0.03369309 Signifikan
## 35 14 -8.920971e-03 0.03369309 -0.03369309 Tidak Signifikan
## 36 15 -3.171234e-02 0.03369309 -0.03369309 Tidak Signifikan
## 37 16 3.208040e-02 0.03369309 -0.03369309 Tidak Signifikan
## 38 17 -1.543306e-03 0.03369309 -0.03369309 Tidak Signifikan
## 39 18 -9.242605e-03 0.03369309 -0.03369309 Tidak Signifikan
## 40 19 -7.901571e-05 0.03369309 -0.03369309 Tidak Signifikan
## 41 20 -2.141856e-02 0.03369309 -0.03369309 Tidak Signifikan
# =========================
# 5. DATA ARIMAX DAN IDENTIFIKASI ORDE
# =========================
data_train_arimax <- data_train %>% mutate(Xt = Return_SP500_Exog, Xt_1 = lag(Return_SP500_Exog, 1)) %>% filter(!is.na(Xt_1))
data_test_arimax <- data_test %>% mutate(Xt = Return_SP500_Exog, Xt_1 = lag(Return_SP500_Exog, 1)) %>% filter(!is.na(Xt_1))
y_train <- data_train_arimax$Return_IHSG; y_test <- data_test_arimax$Return_IHSG
xmat <- function(data, cols) {x <- as.matrix(data[, cols, drop = FALSE]); colnames(x) <- cols; x}
x_train_xt <- xmat(data_train_arimax, "Xt"); x_test_xt <- xmat(data_test_arimax, "Xt")
x_train_xt1 <- xmat(data_train_arimax, "Xt_1"); x_test_xt1 <- xmat(data_test_arimax, "Xt_1")
x_train_xt_xt1 <- xmat(data_train_arimax, c("Xt","Xt_1")); x_test_xt_xt1 <- xmat(data_test_arimax, c("Xt","Xt_1"))
acf(y_train, lag.max = 20, main = "ACF Return IHSG", xlab = "Lag", ylab = "ACF", xlim = c(1,20), ylim = c(-0.06,0.06))
pacf(y_train, lag.max = 20, main = "PACF Return IHSG", xlab = "Lag", ylab = "PACF")
# =========================
# 6. ESTIMASI KANDIDAT ARIMAX
# =========================
p_kandidat <- 0:5; q_kandidat <- 0:5; lag_ljung <- 15
skenario_xreg <- list(Xt = x_train_xt, Xt_1 = x_train_xt1, Xt_dan_Xt_1 = x_train_xt_xt1)
uji_signifikansi <- function(model) {est <- coef(model); vc <- model$var.coef; if(is.null(vc)) return(NULL); d <- diag(vc); if(any(is.na(d)) || any(d <= 0)) return(NULL); se <- sqrt(d); z <- est / se; p <- 2 * pnorm(abs(z), lower.tail = FALSE); data.frame(Parameter = names(est), Estimate = as.numeric(est), Std_Error = as.numeric(se), Z_value = as.numeric(z), p_value = as.numeric(p), row.names = NULL)}
ambil_pvalue <- function(obj) if(is.null(obj)) NA_real_ else obj$p.value
estimasi_kandidat_arimax <- function(y, xreg, nama_skenario, p_set = p_kandidat, q_set = q_kandidat, lag_lb = lag_ljung) {
kandidat <- expand.grid(p = p_set, q = q_set); hasil_model <- list(); rekap_model <- list(); nama_xreg <- colnames(xreg)
for(i in 1:nrow(kandidat)) {
p <- kandidat$p[i]; q <- kandidat$q[i]
fit <- tryCatch(forecast::Arima(y, order = c(p,0,q), xreg = xreg, include.mean = TRUE, method = "ML"), error = function(e) NULL)
if(!is.null(fit)) {
sig <- tryCatch(uji_signifikansi(fit), error = function(e) NULL)
lb15 <- tryCatch(Box.test(residuals(fit), lag = lag_lb, type = "Ljung-Box", fitdf = p + q), error = function(e) NULL)
param_arma <- if(!is.null(sig)) sig[grepl("^(ar|ma)", sig$Parameter, ignore.case = TRUE), , drop = FALSE] else NULL
arma_signif <- if(!is.null(param_arma) && nrow(param_arma) > 0) all(param_arma$p_value < 0.05) else TRUE
xreg_signif <- if(!is.null(sig)) all(sapply(nama_xreg, function(v) if(v %in% sig$Parameter) sig$p_value[sig$Parameter == v][1] < 0.05 else FALSE)) else FALSE
nama_model <- paste0("ARIMAX(", p, ",0,", q, ")"); hasil_model[[nama_model]] <- fit
rekap_model[[length(rekap_model) + 1]] <- data.frame(Skenario = nama_skenario, Model = nama_model, p = p, q = q, AIC = AIC(fit), BIC = BIC(fit), Ljung_Box_15 = ambil_pvalue(lb15), ARMA_Signif = arma_signif, Xreg_Signif = xreg_signif)
}
}
list(hasil_model = hasil_model, rekap_model = bind_rows(rekap_model) %>% arrange(AIC, BIC))
}
hasil_arimax <- setNames(lapply(names(skenario_xreg), function(nm) estimasi_kandidat_arimax(y_train, skenario_xreg[[nm]], nm)), names(skenario_xreg))
hasil_xt <- hasil_arimax$Xt; hasil_xt1 <- hasil_arimax$Xt_1; hasil_xt_xt1 <- hasil_arimax$Xt_dan_Xt_1
rekap_model_arimax <- bind_rows(lapply(hasil_arimax, `[[`, "rekap_model")) %>% arrange(AIC, BIC)
model_lolos_lb15_xreg <- rekap_model_arimax %>% filter(Ljung_Box_15 > 0.05, Xreg_Signif) %>% arrange(AIC, BIC)
rekap_model_arimax; model_lolos_lb15_xreg
## Skenario Model p q AIC BIC Ljung_Box_15 ARMA_Signif
## 1 Xt ARIMAX(4,0,4) 4 4 -21739.10 -21671.71 0.0345126434 TRUE
## 2 Xt ARIMAX(4,0,1) 4 1 -21738.10 -21689.09 0.2280520025 TRUE
## 3 Xt_dan_Xt_1 ARIMAX(4,0,4) 4 4 -21737.14 -21663.62 0.0351533998 TRUE
## 4 Xt ARIMAX(2,0,4) 2 4 -21737.12 -21681.98 0.2062515325 FALSE
## 5 Xt ARIMAX(1,0,4) 1 4 -21737.00 -21687.99 0.1855400912 TRUE
## 6 Xt ARIMAX(4,0,2) 4 2 -21736.54 -21681.40 0.1878644752 FALSE
## 7 Xt ARIMAX(5,0,1) 5 1 -21736.34 -21681.20 0.1751061994 FALSE
## 8 Xt ARIMAX(4,0,3) 4 3 -21736.33 -21675.07 0.2446771502 FALSE
## 9 Xt_dan_Xt_1 ARIMAX(4,0,1) 4 1 -21736.14 -21681.00 0.2267472682 TRUE
## 10 Xt ARIMAX(1,0,5) 1 5 -21735.63 -21680.49 0.1541749836 FALSE
## 11 Xt_dan_Xt_1 ARIMAX(1,0,4) 1 4 -21735.04 -21679.90 0.1838410786 TRUE
## 12 Xt_dan_Xt_1 ARIMAX(4,0,2) 4 2 -21735.01 -21673.74 0.2059149218 FALSE
## 13 Xt_dan_Xt_1 ARIMAX(2,0,4) 2 4 -21734.70 -21673.43 0.1933831459 FALSE
## 14 Xt ARIMAX(5,0,3) 5 3 -21734.45 -21667.06 0.1928742031 FALSE
## 15 Xt_dan_Xt_1 ARIMAX(4,0,3) 4 3 -21734.38 -21666.99 0.2437917670 FALSE
## 16 Xt_dan_Xt_1 ARIMAX(5,0,1) 5 1 -21734.37 -21673.10 0.1740272313 FALSE
## 17 Xt ARIMAX(5,0,2) 5 2 -21734.10 -21672.84 0.1145735947 TRUE
## 18 Xt_dan_Xt_1 ARIMAX(5,0,5) 5 5 -21733.79 -21648.02 0.0487756498 TRUE
## 19 Xt ARIMAX(3,0,5) 3 5 -21733.78 -21666.39 0.1262074844 TRUE
## 20 Xt_dan_Xt_1 ARIMAX(1,0,5) 1 5 -21733.66 -21672.39 0.1520709782 FALSE
## 21 Xt_dan_Xt_1 ARIMAX(5,0,4) 5 4 -21733.28 -21653.64 0.0262021712 TRUE
## 22 Xt ARIMAX(0,0,5) 0 5 -21733.20 -21684.19 0.0625592103 FALSE
## 23 Xt ARIMAX(2,0,5) 2 5 -21733.19 -21671.93 0.0928622289 TRUE
## 24 Xt ARIMAX(5,0,4) 5 4 -21732.73 -21659.21 0.1159741707 FALSE
## 25 Xt ARIMAX(4,0,5) 4 5 -21732.72 -21659.20 0.1273585357 FALSE
## 26 Xt_dan_Xt_1 ARIMAX(3,0,5) 3 5 -21732.54 -21659.03 0.1697897165 FALSE
## 27 Xt_dan_Xt_1 ARIMAX(5,0,3) 5 3 -21732.48 -21658.96 0.1779204681 FALSE
## 28 Xt ARIMAX(5,0,0) 5 0 -21732.29 -21683.28 0.0467953625 FALSE
## 29 Xt_dan_Xt_1 ARIMAX(5,0,2) 5 2 -21732.14 -21664.75 0.1131920010 TRUE
## 30 Xt ARIMAX(5,0,5) 5 5 -21731.87 -21652.23 0.1071547566 TRUE
## 31 Xt_dan_Xt_1 ARIMAX(2,0,5) 2 5 -21731.23 -21663.84 0.0918482622 TRUE
## 32 Xt_dan_Xt_1 ARIMAX(0,0,5) 0 5 -21731.21 -21676.07 0.0619300035 FALSE
## 33 Xt_dan_Xt_1 ARIMAX(4,0,5) 4 5 -21730.88 -21651.23 0.1239483060 TRUE
## 34 Xt ARIMAX(4,0,0) 4 0 -21730.78 -21687.90 0.0234060834 FALSE
## 35 Xt ARIMAX(0,0,4) 0 4 -21730.41 -21687.52 0.0213091690 FALSE
## 36 Xt_dan_Xt_1 ARIMAX(5,0,0) 5 0 -21730.31 -21675.18 0.0462351823 FALSE
## 37 Xt ARIMAX(2,0,2) 2 2 -21729.81 -21686.92 0.0244067877 FALSE
## 38 Xt ARIMAX(3,0,4) 3 4 -21729.38 -21668.11 0.0295648429 TRUE
## 39 Xt_dan_Xt_1 ARIMAX(4,0,0) 4 0 -21728.80 -21679.79 0.0231367115 FALSE
## 40 Xt_dan_Xt_1 ARIMAX(0,0,4) 0 4 -21728.42 -21679.40 0.0210843963 FALSE
## 41 Xt ARIMAX(3,0,2) 3 2 -21727.92 -21678.91 0.0179118673 FALSE
## 42 Xt ARIMAX(2,0,3) 2 3 -21727.87 -21678.86 0.0179000295 FALSE
## 43 Xt_dan_Xt_1 ARIMAX(2,0,2) 2 2 -21727.82 -21678.81 0.0249447844 FALSE
## 44 Xt_dan_Xt_1 ARIMAX(3,0,4) 3 4 -21727.41 -21660.01 0.0301476385 TRUE
## 45 Xt ARIMAX(3,0,3) 3 3 -21727.23 -21672.09 0.0129668598 TRUE
## 46 Xt ARIMAX(1,0,1) 1 1 -21727.09 -21696.46 0.0074192028 TRUE
## 47 Xt ARIMAX(2,0,1) 2 1 -21726.53 -21689.77 0.0047060585 FALSE
## 48 Xt ARIMAX(1,0,2) 1 2 -21726.49 -21689.73 0.0046991381 FALSE
## 49 Xt_dan_Xt_1 ARIMAX(3,0,2) 3 2 -21725.97 -21670.83 0.0176683907 FALSE
## 50 Xt_dan_Xt_1 ARIMAX(2,0,3) 2 3 -21725.89 -21670.75 0.0171858261 FALSE
## 51 Xt_dan_Xt_1 ARIMAX(1,0,1) 1 1 -21725.17 -21688.41 0.0073318464 TRUE
## 52 Xt ARIMAX(3,0,1) 3 1 -21724.89 -21682.00 0.0031149017 FALSE
## 53 Xt ARIMAX(1,0,3) 1 3 -21724.72 -21681.84 0.0029933189 FALSE
## 54 Xt_dan_Xt_1 ARIMAX(2,0,1) 2 1 -21724.55 -21681.67 0.0047686961 FALSE
## 55 Xt_dan_Xt_1 ARIMAX(1,0,2) 1 2 -21724.50 -21681.62 0.0047449456 FALSE
## 56 Xt ARIMAX(0,0,3) 0 3 -21724.08 -21687.32 0.0020855269 FALSE
## 57 Xt ARIMAX(0,0,2) 0 2 -21723.92 -21693.29 0.0018517563 FALSE
## 58 Xt ARIMAX(3,0,0) 3 0 -21723.42 -21686.67 0.0016211699 FALSE
## 59 Xt ARIMAX(2,0,0) 2 0 -21723.07 -21692.43 0.0014098080 FALSE
## 60 Xt_dan_Xt_1 ARIMAX(3,0,1) 3 1 -21722.91 -21673.90 0.0030901451 FALSE
## 61 Xt_dan_Xt_1 ARIMAX(1,0,3) 1 3 -21722.76 -21673.75 0.0030145690 FALSE
## 62 Xt_dan_Xt_1 ARIMAX(0,0,3) 0 3 -21722.09 -21679.20 0.0020601232 FALSE
## 63 Xt_dan_Xt_1 ARIMAX(0,0,2) 0 2 -21721.93 -21685.17 0.0018358946 FALSE
## 64 Xt_dan_Xt_1 ARIMAX(3,0,0) 3 0 -21721.44 -21678.56 0.0015981226 FALSE
## 65 Xt_dan_Xt_1 ARIMAX(2,0,0) 2 0 -21721.08 -21684.32 0.0013896020 FALSE
## 66 Xt ARIMAX(0,0,0) 0 0 -21721.01 -21702.63 0.0007689282 TRUE
## 67 Xt ARIMAX(0,0,1) 0 1 -21720.48 -21695.98 0.0007995737 FALSE
## 68 Xt ARIMAX(1,0,0) 1 0 -21720.37 -21695.86 0.0007659966 FALSE
## 69 Xt_dan_Xt_1 ARIMAX(0,0,0) 0 0 -21719.47 -21694.96 0.0008719408 TRUE
## 70 Xt_dan_Xt_1 ARIMAX(0,0,1) 0 1 -21718.84 -21688.20 0.0008703596 FALSE
## 71 Xt_dan_Xt_1 ARIMAX(1,0,0) 1 0 -21718.73 -21688.10 0.0008387847 FALSE
## 72 Xt_dan_Xt_1 ARIMAX(3,0,3) 3 3 -21717.80 -21656.54 0.0006836181 TRUE
## 73 Xt_1 ARIMAX(5,0,4) 5 4 -21461.67 -21388.15 0.0652177115 TRUE
## 74 Xt_1 ARIMAX(4,0,5) 4 5 -21461.57 -21388.05 0.0626271775 TRUE
## 75 Xt_1 ARIMAX(5,0,5) 5 5 -21454.31 -21374.66 0.0054020628 TRUE
## 76 Xt_1 ARIMAX(4,0,4) 4 4 -21448.67 -21381.28 0.0135978796 TRUE
## 77 Xt_1 ARIMAX(4,0,1) 4 1 -21442.98 -21393.97 0.0157962350 FALSE
## 78 Xt_1 ARIMAX(4,0,2) 4 2 -21442.49 -21387.35 0.0164745643 FALSE
## 79 Xt_1 ARIMAX(1,0,4) 1 4 -21442.36 -21393.35 0.0141034949 FALSE
## 80 Xt_1 ARIMAX(5,0,1) 5 1 -21441.20 -21386.06 0.0100724383 FALSE
## 81 Xt_1 ARIMAX(1,0,5) 1 5 -21440.58 -21385.44 0.0088933910 FALSE
## 82 Xt_1 ARIMAX(3,0,5) 3 5 -21440.36 -21372.97 0.0124205560 TRUE
## 83 Xt_1 ARIMAX(0,0,4) 0 4 -21440.23 -21397.34 0.0055856769 FALSE
## 84 Xt_1 ARIMAX(5,0,3) 5 3 -21440.08 -21372.69 0.0163556968 TRUE
## 85 Xt_1 ARIMAX(4,0,0) 4 0 -21440.07 -21397.19 0.0052159155 FALSE
## 86 Xt_1 ARIMAX(2,0,4) 2 4 -21439.98 -21384.84 0.0072422192 FALSE
## 87 Xt_1 ARIMAX(2,0,2) 2 2 -21439.90 -21397.02 0.0053252732 TRUE
## 88 Xt_1 ARIMAX(3,0,3) 3 3 -21439.86 -21384.73 0.0077600709 FALSE
## 89 Xt_1 ARIMAX(0,0,5) 0 5 -21439.68 -21390.66 0.0049921264 FALSE
## 90 Xt_1 ARIMAX(5,0,0) 5 0 -21439.27 -21390.26 0.0040792145 FALSE
## 91 Xt_1 ARIMAX(5,0,2) 5 2 -21439.01 -21377.75 0.0051673967 TRUE
## 92 Xt_1 ARIMAX(0,0,1) 0 1 -21438.83 -21414.32 0.0025329047 TRUE
## 93 Xt_1 ARIMAX(2,0,1) 2 1 -21438.49 -21401.73 0.0029681986 TRUE
## 94 Xt_1 ARIMAX(1,0,0) 1 0 -21438.47 -21413.97 0.0023066958 TRUE
## 95 Xt_1 ARIMAX(2,0,3) 2 3 -21438.46 -21389.45 0.0034526696 TRUE
## 96 Xt_1 ARIMAX(2,0,5) 2 5 -21438.40 -21377.13 0.0045818110 FALSE
## 97 Xt_1 ARIMAX(1,0,2) 1 2 -21438.05 -21401.29 0.0025575695 TRUE
## 98 Xt_1 ARIMAX(2,0,0) 2 0 -21437.87 -21407.24 0.0019665572 FALSE
## 99 Xt_1 ARIMAX(3,0,2) 3 2 -21437.77 -21388.76 0.0028475995 TRUE
## 100 Xt_1 ARIMAX(0,0,2) 0 2 -21437.67 -21407.04 0.0018497282 FALSE
## 101 Xt_1 ARIMAX(3,0,0) 3 0 -21437.62 -21400.86 0.0022117559 FALSE
## 102 Xt_1 ARIMAX(3,0,1) 3 1 -21437.34 -21394.45 0.0023735849 FALSE
## 103 Xt_1 ARIMAX(4,0,3) 4 3 -21437.30 -21376.04 0.0030512322 TRUE
## 104 Xt_1 ARIMAX(1,0,1) 1 1 -21437.28 -21406.65 0.0016706575 FALSE
## 105 Xt_1 ARIMAX(0,0,3) 0 3 -21437.11 -21400.35 0.0018670548 FALSE
## 106 Xt_1 ARIMAX(1,0,3) 1 3 -21436.89 -21394.01 0.0020280132 FALSE
## 107 Xt_1 ARIMAX(3,0,4) 3 4 -21436.88 -21375.62 0.0026235701 TRUE
## 108 Xt_1 ARIMAX(0,0,0) 0 0 -21431.73 -21413.35 0.0003135110 TRUE
## Xreg_Signif
## 1 FALSE
## 2 TRUE
## 3 FALSE
## 4 TRUE
## 5 TRUE
## 6 TRUE
## 7 TRUE
## 8 TRUE
## 9 FALSE
## 10 TRUE
## 11 FALSE
## 12 FALSE
## 13 FALSE
## 14 TRUE
## 15 FALSE
## 16 FALSE
## 17 FALSE
## 18 FALSE
## 19 FALSE
## 20 FALSE
## 21 FALSE
## 22 TRUE
## 23 FALSE
## 24 TRUE
## 25 TRUE
## 26 FALSE
## 27 FALSE
## 28 TRUE
## 29 FALSE
## 30 FALSE
## 31 FALSE
## 32 FALSE
## 33 FALSE
## 34 TRUE
## 35 TRUE
## 36 FALSE
## 37 TRUE
## 38 FALSE
## 39 FALSE
## 40 FALSE
## 41 TRUE
## 42 TRUE
## 43 FALSE
## 44 FALSE
## 45 TRUE
## 46 TRUE
## 47 TRUE
## 48 TRUE
## 49 FALSE
## 50 FALSE
## 51 FALSE
## 52 TRUE
## 53 TRUE
## 54 FALSE
## 55 FALSE
## 56 TRUE
## 57 TRUE
## 58 TRUE
## 59 TRUE
## 60 FALSE
## 61 FALSE
## 62 FALSE
## 63 FALSE
## 64 FALSE
## 65 FALSE
## 66 TRUE
## 67 TRUE
## 68 TRUE
## 69 FALSE
## 70 FALSE
## 71 FALSE
## 72 FALSE
## 73 TRUE
## 74 TRUE
## 75 TRUE
## 76 FALSE
## 77 TRUE
## 78 TRUE
## 79 TRUE
## 80 TRUE
## 81 TRUE
## 82 FALSE
## 83 TRUE
## 84 FALSE
## 85 TRUE
## 86 TRUE
## 87 TRUE
## 88 TRUE
## 89 TRUE
## 90 TRUE
## 91 FALSE
## 92 TRUE
## 93 TRUE
## 94 TRUE
## 95 FALSE
## 96 TRUE
## 97 TRUE
## 98 TRUE
## 99 FALSE
## 100 TRUE
## 101 TRUE
## 102 TRUE
## 103 FALSE
## 104 TRUE
## 105 TRUE
## 106 TRUE
## 107 FALSE
## 108 TRUE
## Skenario Model p q AIC BIC Ljung_Box_15 ARMA_Signif
## 1 Xt ARIMAX(4,0,1) 4 1 -21738.10 -21689.09 0.22805200 TRUE
## 2 Xt ARIMAX(2,0,4) 2 4 -21737.12 -21681.98 0.20625153 FALSE
## 3 Xt ARIMAX(1,0,4) 1 4 -21737.00 -21687.99 0.18554009 TRUE
## 4 Xt ARIMAX(4,0,2) 4 2 -21736.54 -21681.40 0.18786448 FALSE
## 5 Xt ARIMAX(5,0,1) 5 1 -21736.34 -21681.20 0.17510620 FALSE
## 6 Xt ARIMAX(4,0,3) 4 3 -21736.33 -21675.07 0.24467715 FALSE
## 7 Xt ARIMAX(1,0,5) 1 5 -21735.63 -21680.49 0.15417498 FALSE
## 8 Xt ARIMAX(5,0,3) 5 3 -21734.45 -21667.06 0.19287420 FALSE
## 9 Xt ARIMAX(0,0,5) 0 5 -21733.20 -21684.19 0.06255921 FALSE
## 10 Xt ARIMAX(5,0,4) 5 4 -21732.73 -21659.21 0.11597417 FALSE
## 11 Xt ARIMAX(4,0,5) 4 5 -21732.72 -21659.20 0.12735854 FALSE
## 12 Xt_1 ARIMAX(5,0,4) 5 4 -21461.67 -21388.15 0.06521771 TRUE
## 13 Xt_1 ARIMAX(4,0,5) 4 5 -21461.57 -21388.05 0.06262718 TRUE
## Xreg_Signif
## 1 TRUE
## 2 TRUE
## 3 TRUE
## 4 TRUE
## 5 TRUE
## 6 TRUE
## 7 TRUE
## 8 TRUE
## 9 TRUE
## 10 TRUE
## 11 TRUE
## 12 TRUE
## 13 TRUE
# =========================
# 7. MODEL ARIMAX TERBAIK DAN DIAGNOSTIK
# =========================
model_arimax_terbaik <- hasil_xt$hasil_model[["ARIMAX(4,0,1)"]]
x_train_terbaik <- x_train_xt; x_test_terbaik <- x_test_xt
cat("Skenario terbaik: Xt\nModel ARIMAX terbaik: ARIMAX(4,0,1)\nAIC:", AIC(model_arimax_terbaik), "\nBIC:", BIC(model_arimax_terbaik), "\n")
## Skenario terbaik: Xt
## Model ARIMAX terbaik: ARIMAX(4,0,1)
## AIC: -21738.1
## BIC: -21689.09
summary(model_arimax_terbaik); uji_signifikansi(model_arimax_terbaik)
## Series: y
## Regression with ARIMA(4,0,1) errors
##
## Coefficients:
## ar1 ar2 ar3 ar4 ma1 intercept Xt
## -0.6728 -0.056 -0.0533 -0.0742 0.6510 1e-04 0.2836
## s.e. 0.1017 0.021 0.0212 0.0177 0.1007 1e-04 0.0160
##
## sigma^2 = 9.456e-05: log likelihood = 10877.05
## AIC=-21738.1 AICc=-21738.06 BIC=-21689.09
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -4.204684e-06 0.009714347 0.006784578 NaN Inf 0.6697186
## ACF1
## Training set 0.0007546374
## Parameter Estimate Std_Error Z_value p_value
## 1 ar1 -6.728417e-01 0.1017060436 -6.6155528 3.701659e-11
## 2 ar2 -5.597164e-02 0.0210315382 -2.6613193 7.783510e-03
## 3 ar3 -5.328292e-02 0.0211512969 -2.5191324 1.176444e-02
## 4 ar4 -7.422414e-02 0.0177204471 -4.1886153 2.806617e-05
## 5 ma1 6.509507e-01 0.1007468543 6.4612510 1.038409e-10
## 6 intercept 5.915819e-05 0.0001497484 0.3950507 6.928055e-01
## 7 Xt 2.835807e-01 0.0159566029 17.7719942 1.164754e-70
res_arimax <- na.omit(as.numeric(residuals(model_arimax_terbaik))); res_arimax2 <- res_arimax^2; fitdf_arimax <- 4 + 1
ljung_res_15 <- Box.test(res_arimax, lag = 15, type = "Ljung-Box", fitdf = fitdf_arimax)
jb_res <- tseries::jarque.bera.test(res_arimax)
par(mfrow = c(1,2)); acf(res_arimax2, lag.max = 30, main = "ACF Residual Kuadrat ARIMAX"); pacf(res_arimax2, lag.max = 30, main = "PACF Residual Kuadrat ARIMAX"); par(mfrow = c(1,1))
arch_lm_15 <- FinTS::ArchTest(res_arimax, lags = 15)
rekap_diagnostik_residual <- data.frame(Uji = c("Ljung-Box Residual lag 15","Jarque-Bera Normalitas Residual","ARCH-LM lag 15"), Statistik = c(as.numeric(ljung_res_15$statistic), as.numeric(jb_res$statistic), as.numeric(arch_lm_15$statistic)), P_Value = c(ljung_res_15$p.value, jb_res$p.value, arch_lm_15$p.value), Keputusan = c(ifelse(ljung_res_15$p.value > 0.05, "Tidak terdapat autokorelasi residual", "Terdapat autokorelasi residual"), ifelse(jb_res$p.value > 0.05, "Residual berdistribusi normal", "Residual tidak berdistribusi normal"), ifelse(arch_lm_15$p.value > 0.05, "Tidak terdapat efek ARCH/GARCH", "Terdapat efek ARCH/GARCH")))
rekap_diagnostik_residual
## Uji Statistik P_Value
## 1 Ljung-Box Residual lag 15 12.92224 2.280520e-01
## 2 Jarque-Bera Normalitas Residual 8973.72566 0.000000e+00
## 3 ARCH-LM lag 15 664.95295 5.434552e-132
## Keputusan
## 1 Tidak terdapat autokorelasi residual
## 2 Residual tidak berdistribusi normal
## 3 Terdapat efek ARCH/GARCH
# =========================
# 8. ARIMAX-GARCH(1,1)
# =========================
p_arimax <- 4; q_arimax <- 1; P_garch <- 1; Q_garch <- 1; lag_ljung <- 15
y_garch <- as.numeric(y_train); xreg_garch <- as.matrix(x_train_terbaik); stopifnot(length(y_garch) == nrow(xreg_garch))
spec_arimax_garch_11 <- rugarch::ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(P_garch,Q_garch)), mean.model = list(armaOrder = c(p_arimax,q_arimax), include.mean = TRUE, external.regressors = xreg_garch), distribution.model = "norm")
model_arimax_garch_terbaik <- rugarch::ugarchfit(spec = spec_arimax_garch_11, data = y_garch, solver = "hybrid")
model_arimax_garch_terbaik
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(4,0,1)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.000187 0.000111 1.6785 0.093254
## ar1 -0.675781 0.222168 -3.0418 0.002352
## ar2 -0.105452 0.024085 -4.3783 0.000012
## ar3 -0.094599 0.028050 -3.3725 0.000745
## ar4 -0.056395 0.018587 -3.0341 0.002412
## ma1 0.634412 0.222062 2.8569 0.004278
## mxreg1 0.252301 0.015536 16.2402 0.000000
## omega 0.000003 0.000002 1.8617 0.062644
## alpha1 0.103854 0.009733 10.6705 0.000000
## beta1 0.860659 0.005853 147.0387 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.000187 0.000234 0.79956 0.423966
## ar1 -0.675781 0.205084 -3.29515 0.000984
## ar2 -0.105452 0.026596 -3.96499 0.000073
## ar3 -0.094599 0.028801 -3.28457 0.001021
## ar4 -0.056395 0.022259 -2.53363 0.011289
## ma1 0.634412 0.206089 3.07834 0.002082
## mxreg1 0.252301 0.025042 10.07506 0.000000
## omega 0.000003 0.000014 0.22201 0.824307
## alpha1 0.103854 0.034720 2.99117 0.002779
## beta1 0.860659 0.101088 8.51398 0.000000
##
## LogLikelihood : 11309.17
##
## Information Criteria
## ------------------------------------
##
## Akaike -6.6800
## Bayes -6.6619
## Shibata -6.6800
## Hannan-Quinn -6.6735
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 1.365 0.2426
## Lag[2*(p+q)+(p+q)-1][14] 5.113 1.0000
## Lag[4*(p+q)+(p+q)-1][24] 9.476 0.8805
## d.o.f=5
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 3.124 0.07715
## Lag[2*(p+q)+(p+q)-1][5] 3.546 0.31628
## Lag[4*(p+q)+(p+q)-1][9] 4.854 0.45110
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.03916 0.500 2.000 0.8431
## ARCH Lag[5] 1.05128 1.440 1.667 0.7178
## ARCH Lag[7] 1.88540 2.315 1.543 0.7416
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 7.9528
## Individual Statistics:
## mu 0.27224
## ar1 0.22155
## ar2 0.13800
## ar3 0.07936
## ar4 0.76100
## ma1 0.21916
## mxreg1 3.20693
## omega 0.08314
## alpha1 0.15000
## beta1 0.09303
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 2.29 2.54 3.05
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 0.01242 9.901e-01
## Negative Sign Bias 3.44084 5.870e-04 ***
## Positive Sign Bias 1.25723 2.088e-01
## Joint Effect 23.72328 2.853e-05 ***
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 107.5 2.270e-14
## 2 30 129.0 1.402e-14
## 3 40 143.6 6.952e-14
## 4 50 146.9 1.046e-11
##
##
## Elapsed time : 0.8878551
coef_garch <- model_arimax_garch_terbaik@fit$matcoef
tabel_koef_garch <- data.frame(Parameter = rownames(coef_garch), Estimate = coef_garch[,1], Std_Error = coef_garch[,2], t_value = coef_garch[,3], p_value = coef_garch[,4], Keputusan = ifelse(coef_garch[,4] < 0.05, "Signifikan", "Tidak Signifikan"), row.names = NULL)
tabel_koef_garch
## Parameter Estimate Std_Error t_value p_value
## 1 mu 1.868400e-04 1.113153e-04 1.678475 9.325444e-02
## 2 ar1 -6.757810e-01 2.221679e-01 -3.041758 2.352007e-03
## 3 ar2 -1.054520e-01 2.408493e-02 -4.378339 1.195875e-05
## 4 ar3 -9.459902e-02 2.805009e-02 -3.372504 7.448805e-04
## 5 ar4 -5.639537e-02 1.858688e-02 -3.034150 2.412142e-03
## 6 ma1 6.344124e-01 2.220624e-01 2.856910 4.277872e-03
## 7 mxreg1 2.523009e-01 1.553559e-02 16.240184 0.000000e+00
## 8 omega 3.194582e-06 1.715942e-06 1.861708 6.264429e-02
## 9 alpha1 1.038535e-01 9.732738e-03 10.670534 0.000000e+00
## 10 beta1 8.606585e-01 5.853280e-03 147.038663 0.000000e+00
## Keputusan
## 1 Tidak Signifikan
## 2 Signifikan
## 3 Signifikan
## 4 Signifikan
## 5 Signifikan
## 6 Signifikan
## 7 Signifikan
## 8 Tidak Signifikan
## 9 Signifikan
## 10 Signifikan
std_resid <- na.omit(as.numeric(rugarch::residuals(model_arimax_garch_terbaik, standardize = TRUE))); std_resid2 <- std_resid^2
uji_lb_std <- Box.test(std_resid, lag = lag_ljung, type = "Ljung-Box", fitdf = p_arimax + q_arimax)
uji_lb_std2 <- Box.test(std_resid2, lag = lag_ljung, type = "Ljung-Box")
uji_arch_std <- FinTS::ArchTest(std_resid, lags = lag_ljung)
uji_jb_std <- tseries::jarque.bera.test(std_resid)
ringkasan_diagnostik_garch <- data.frame(Uji = c("Ljung-Box standardized residual lag 15","Ljung-Box standardized residual kuadrat lag 15","ARCH-LM standardized residual lag 15","Jarque-Bera standardized residual"), Statistik = c(as.numeric(uji_lb_std$statistic), as.numeric(uji_lb_std2$statistic), as.numeric(uji_arch_std$statistic), as.numeric(uji_jb_std$statistic)), P_Value = c(uji_lb_std$p.value, uji_lb_std2$p.value, uji_arch_std$p.value, uji_jb_std$p.value), Keputusan = c(ifelse(uji_lb_std$p.value > 0.05, "Tidak ada autokorelasi standardized residual", "Masih ada autokorelasi standardized residual"), ifelse(uji_lb_std2$p.value > 0.05, "Tidak ada autokorelasi standardized residual kuadrat", "Masih ada autokorelasi standardized residual kuadrat"), ifelse(uji_arch_std$p.value > 0.05, "Tidak ada efek ARCH/GARCH sisa", "Masih ada efek ARCH/GARCH sisa"), ifelse(uji_jb_std$p.value > 0.05, "Standardized residual berdistribusi normal", "Standardized residual tidak berdistribusi normal")))
ringkasan_diagnostik_garch
## Uji Statistik P_Value
## 1 Ljung-Box standardized residual lag 15 12.300282 0.2654623
## 2 Ljung-Box standardized residual kuadrat lag 15 9.942555 0.8233377
## 3 ARCH-LM standardized residual lag 15 8.557149 0.8995122
## 4 Jarque-Bera standardized residual 1765.585351 0.0000000
## Keputusan
## 1 Tidak ada autokorelasi standardized residual
## 2 Tidak ada autokorelasi standardized residual kuadrat
## 3 Tidak ada efek ARCH/GARCH sisa
## 4 Standardized residual tidak berdistribusi normal
# =========================
# 9. PREDIKSI DATA UJI: ROLLING ONE-STEP-AHEAD
# =========================
n_test <- nrow(x_test_terbaik); pred_mean <- pred_sigma <- rep(NA_real_, n_test)
for(i in 1:n_test) {
y_hist <- c(as.numeric(data_train_arimax$Return_IHSG), if(i == 1) NULL else as.numeric(data_test_arimax$Return_IHSG[1:(i-1)]))
x_hist <- if(i == 1) as.matrix(x_train_terbaik) else rbind(as.matrix(x_train_terbaik), as.matrix(x_test_terbaik[1:(i-1), , drop = FALSE]))
x_next <- as.matrix(x_test_terbaik[i, , drop = FALSE]); colnames(x_next) <- colnames(x_train_terbaik)
spec_i <- rugarch::ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(P_garch,Q_garch)), mean.model = list(armaOrder = c(p_arimax,q_arimax), include.mean = TRUE, external.regressors = x_hist), distribution.model = "norm")
fit_i <- tryCatch(rugarch::ugarchfit(spec = spec_i, data = y_hist, solver = "hybrid"), error = function(e) NULL)
if(!is.null(fit_i) && fit_i@fit$convergence == 0) {fc_i <- tryCatch(rugarch::ugarchforecast(fit_i, n.ahead = 1, external.forecasts = list(mregfor = x_next)), error = function(e) NULL); if(!is.null(fc_i)) {pred_mean[i] <- as.numeric(fitted(fc_i)); pred_sigma[i] <- as.numeric(sigma(fc_i))}}
if(i %% 25 == 0 || i == n_test) cat("Prediksi ke-", i, " dari ", n_test, " selesai\n", sep = "")
}
## Prediksi ke-25 dari 376 selesai
## Prediksi ke-50 dari 376 selesai
## Prediksi ke-75 dari 376 selesai
## Prediksi ke-100 dari 376 selesai
## Prediksi ke-125 dari 376 selesai
## Prediksi ke-150 dari 376 selesai
## Prediksi ke-175 dari 376 selesai
## Prediksi ke-200 dari 376 selesai
## Prediksi ke-225 dari 376 selesai
## Prediksi ke-250 dari 376 selesai
## Prediksi ke-275 dari 376 selesai
## Prediksi ke-300 dari 376 selesai
## Prediksi ke-325 dari 376 selesai
## Prediksi ke-350 dari 376 selesai
## Prediksi ke-375 dari 376 selesai
## Prediksi ke-376 dari 376 selesai
hasil_uji <- data.frame(Date = data_test_arimax$Date[1:n_test], Return_Aktual = as.numeric(data_test_arimax$Return_IHSG[1:n_test]), Return_Prediksi = pred_mean, Volatilitas_Prediksi = pred_sigma) %>% filter(!is.na(Return_Prediksi), !is.na(Volatilitas_Prediksi)) %>% mutate(Volatilitas_Aktual = abs(Return_Aktual))
head(hasil_uji); tail(hasil_uji)
## Date Return_Aktual Return_Prediksi Volatilitas_Prediksi
## 1 2024-11-12 0.007612909 3.031242e-03 0.011194613
## 2 2024-11-13 -0.001820834 5.294784e-04 0.010630117
## 3 2024-11-14 -0.012960108 -6.911561e-05 0.010044875
## 4 2024-11-15 -0.007415263 -1.364952e-03 0.010355328
## 5 2024-11-18 -0.003774608 -1.787614e-03 0.009962973
## 6 2024-11-19 0.008573680 2.285253e-03 0.009437662
## Volatilitas_Aktual
## 1 0.007612909
## 2 0.001820834
## 3 0.012960108
## 4 0.007415263
## 5 0.003774608
## 6 0.008573680
## Date Return_Aktual Return_Prediksi Volatilitas_Prediksi
## 371 2026-06-17 -0.0054874770 -0.0034345328 0.02930924
## 372 2026-06-18 -0.0078108503 -0.0065935345 0.02752599
## 373 2026-06-19 0.0007773607 0.0020848062 0.02583810
## 374 2026-06-22 -0.0098342806 0.0026283842 0.02429807
## 375 2026-06-23 -0.0025143203 0.0006675728 0.02315296
## 376 2026-06-24 -0.0356413512 -0.0026774994 0.02177303
## Volatilitas_Aktual
## 371 0.0054874770
## 372 0.0078108503
## 373 0.0007773607
## 374 0.0098342806
## 375 0.0025143203
## 376 0.0356413512
plot_dua_garis <- function(data, y1, y2, lab1, lab2, ylab) {d <- bind_rows(data %>% transmute(Date, Nilai = .data[[y1]], Keterangan = lab1), data %>% transmute(Date, Nilai = .data[[y2]], Keterangan = lab2)); ggplot(d, aes(Date, Nilai, color = Keterangan)) + geom_line(linewidth = 0.7) + scale_color_manual(values = setNames(c("blue","red"), c(lab1,lab2))) + labs(x = "Date", y = ylab, color = "Keterangan") + tema_ts()}
plot.return.uji <- plot_dua_garis(hasil_uji, "Return_Aktual", "Return_Prediksi", "Actual Return", "Prediction Return", "Return")
plot.vol.uji <- plot_dua_garis(hasil_uji, "Volatilitas_Aktual", "Volatilitas_Prediksi", "Actual Volatility", "Prediction Volatility", "Volatility")
plot.return.uji; plot.vol.uji
ringkasan_prediksi <- data.frame(Jumlah_Data_Uji = nrow(hasil_uji), Rata_rata_Return_Aktual = mean(hasil_uji$Return_Aktual), Rata_rata_Return_Prediksi = mean(hasil_uji$Return_Prediksi), Rata_rata_Volatilitas_Aktual = mean(hasil_uji$Volatilitas_Aktual), Rata_rata_Volatilitas_Prediksi = mean(hasil_uji$Volatilitas_Prediksi), Minimum_Volatilitas_Prediksi = min(hasil_uji$Volatilitas_Prediksi), Maksimum_Volatilitas_Prediksi = max(hasil_uji$Volatilitas_Prediksi))
ringkasan_prediksi
## Jumlah_Data_Uji Rata_rata_Return_Aktual Rata_rata_Return_Prediksi
## 1 376 -0.0005595846 0.0005095166
## Rata_rata_Volatilitas_Aktual Rata_rata_Volatilitas_Prediksi
## 1 0.01033608 0.01232308
## Minimum_Volatilitas_Prediksi Maksimum_Volatilitas_Prediksi
## 1 0.006374188 0.03178697
evaluasi_akurasi <- function(actual, pred) {idx <- is.finite(actual) & is.finite(pred); actual <- as.numeric(actual[idx]); pred <- as.numeric(pred[idx]); mae <- mean(abs(actual - pred)); rmse <- sqrt(mean((actual - pred)^2)); jangkauan <- max(actual) - min(actual); data.frame(n = length(actual), MAE = mae, RMSE = rmse, Jangkauan = jangkauan, Rel_MAE_Persen = ifelse(jangkauan == 0, NA, mae / jangkauan * 100), Rel_RMSE_Persen = ifelse(jangkauan == 0, NA, rmse / jangkauan * 100))}
ringkasan_evaluasi <- bind_rows(evaluasi_akurasi(y_train, fitted(model_arimax_garch_terbaik)) %>% mutate(Data = "Latih", Komponen = "Return"), evaluasi_akurasi(abs(y_train), sigma(model_arimax_garch_terbaik)) %>% mutate(Data = "Latih", Komponen = "Volatilitas"), evaluasi_akurasi(hasil_uji$Return_Aktual, hasil_uji$Return_Prediksi) %>% mutate(Data = "Uji", Komponen = "Return"), evaluasi_akurasi(hasil_uji$Volatilitas_Aktual, hasil_uji$Volatilitas_Prediksi) %>% mutate(Data = "Uji", Komponen = "Volatilitas")) %>% dplyr::select(Data, Komponen, everything())
ringkasan_evaluasi
## Data Komponen n MAE RMSE Jangkauan Rel_MAE_Persen
## 1 Latih Return 3383 0.006755425 0.009733083 0.19003853 3.554766
## 2 Latih Volatilitas 3383 0.005189624 0.006956661 0.09704184 5.347821
## 3 Uji Return 376 0.010067145 0.014727597 0.15530952 6.481988
## 4 Uji Volatilitas 376 0.007546261 0.010718750 0.08229083 9.170233
## Rel_RMSE_Persen
## 1 5.121636
## 2 7.168723
## 3 9.482740
## 4 13.025448
# =========================
# 10. VALUE AT RISK: NORMAL, FHS, STUDENT-t
# =========================
std_resid <- na.omit(as.numeric(rugarch::residuals(model_arimax_garch_terbaik, standardize = TRUE)))
ringkasan_std_resid <- data.frame(Mean = mean(std_resid), Std_Dev = sd(std_resid), Skewness = moments::skewness(std_resid), Kurtosis = moments::kurtosis(std_resid))
uji_jb_std_var <- tseries::jarque.bera.test(std_resid)
ringkasan_std_resid; uji_jb_std_var
## Mean Std_Dev Skewness Kurtosis
## 1 -0.01885228 1.004539 -0.7338683 6.220455
##
## Jarque Bera Test
##
## data: std_resid
## X-squared = 1765.6, df = 2, p-value < 2.2e-16
alpha_95 <- 0.05; alpha_99 <- 0.01
q_norm_95 <- qnorm(alpha_95); q_norm_99 <- qnorm(alpha_99)
q_fhs_95 <- as.numeric(quantile(std_resid, alpha_95, na.rm = TRUE)); q_fhs_99 <- as.numeric(quantile(std_resid, alpha_99, na.rm = TRUE))
fit_t <- tryCatch(MASS::fitdistr(std_resid, densfun = "t"), error = function(e) NULL)
## Warning in dt((x - m)/s, df, log = TRUE): NaNs produced
if(!is.null(fit_t)) {est_t <- fit_t$estimate; loc_t <- if("m" %in% names(est_t)) as.numeric(est_t["m"]) else 0; scale_t <- if("s" %in% names(est_t)) abs(as.numeric(est_t["s"])) else 1; df_t <- if("df" %in% names(est_t)) as.numeric(est_t["df"]) else NA_real_; q_t_95 <- loc_t + scale_t * qt(alpha_95, df_t); q_t_99 <- loc_t + scale_t * qt(alpha_99, df_t)} else {loc_t <- scale_t <- df_t <- q_t_95 <- q_t_99 <- NA_real_}
rekap_kuantil_var <- data.frame(Metode = c("Normal 95%","Normal 99%","FHS 95%","FHS 99%","Student-t 95%","Student-t 99%"), Alpha = c(alpha_95,alpha_99,alpha_95,alpha_99,alpha_95,alpha_99), Kuantil_Standardized_Residual = c(q_norm_95,q_norm_99,q_fhs_95,q_fhs_99,q_t_95,q_t_99))
param_student_t <- data.frame(Location = loc_t, Scale = scale_t, df = df_t)
rekap_kuantil_var; param_student_t
## Metode Alpha Kuantil_Standardized_Residual
## 1 Normal 95% 0.05 -1.644854
## 2 Normal 99% 0.01 -2.326348
## 3 FHS 95% 0.05 -1.687754
## 4 FHS 99% 0.01 -2.947449
## 5 Student-t 95% 0.05 -1.556590
## 6 Student-t 99% 0.01 -2.558209
## Location Scale df
## 1 0.02491875 0.8047289 5.645247
hasil_var <- hasil_uji %>% mutate(Loss_Aktual = -Return_Aktual, VaR_Normal_95 = -(Return_Prediksi + Volatilitas_Prediksi * q_norm_95), VaR_Normal_99 = -(Return_Prediksi + Volatilitas_Prediksi * q_norm_99), VaR_FHS_95 = -(Return_Prediksi + Volatilitas_Prediksi * q_fhs_95), VaR_FHS_99 = -(Return_Prediksi + Volatilitas_Prediksi * q_fhs_99), VaR_t_95 = -(Return_Prediksi + Volatilitas_Prediksi * q_t_95), VaR_t_99 = -(Return_Prediksi + Volatilitas_Prediksi * q_t_99))
head(hasil_var); tail(hasil_var)
## Date Return_Aktual Return_Prediksi Volatilitas_Prediksi
## 1 2024-11-12 0.007612909 3.031242e-03 0.011194613
## 2 2024-11-13 -0.001820834 5.294784e-04 0.010630117
## 3 2024-11-14 -0.012960108 -6.911561e-05 0.010044875
## 4 2024-11-15 -0.007415263 -1.364952e-03 0.010355328
## 5 2024-11-18 -0.003774608 -1.787614e-03 0.009962973
## 6 2024-11-19 0.008573680 2.285253e-03 0.009437662
## Volatilitas_Aktual Loss_Aktual VaR_Normal_95 VaR_Normal_99 VaR_FHS_95
## 1 0.007612909 -0.007612909 0.01538226 0.02301132 0.01586251
## 2 0.001820834 0.001820834 0.01695551 0.02419987 0.01741154
## 3 0.012960108 0.012960108 0.01659147 0.02343699 0.01702239
## 4 0.007415263 0.007415263 0.01839795 0.02545505 0.01884220
## 5 0.003774608 0.003774608 0.01817525 0.02496495 0.01860266
## 6 0.008573680 -0.008573680 0.01323832 0.01967003 0.01364320
## VaR_FHS_99 VaR_t_95 VaR_t_99
## 1 0.02996431 0.01439419 0.02560692
## 2 0.03080225 0.01601726 0.02666459
## 3 0.02967587 0.01570487 0.02576601
## 4 0.03188675 0.01748396 0.02785605
## 5 0.03115297 0.01729588 0.02727499
## 6 0.02553177 0.01240532 0.02185826
## Date Return_Aktual Return_Prediksi Volatilitas_Prediksi
## 371 2026-06-17 -0.0054874770 -0.0034345328 0.02930924
## 372 2026-06-18 -0.0078108503 -0.0065935345 0.02752599
## 373 2026-06-19 0.0007773607 0.0020848062 0.02583810
## 374 2026-06-22 -0.0098342806 0.0026283842 0.02429807
## 375 2026-06-23 -0.0025143203 0.0006675728 0.02315296
## 376 2026-06-24 -0.0356413512 -0.0026774994 0.02177303
## Volatilitas_Aktual Loss_Aktual VaR_Normal_95 VaR_Normal_99 VaR_FHS_95
## 371 0.0054874770 0.0054874770 0.05164394 0.07161802 0.05290131
## 372 0.0078108503 0.0078108503 0.05186976 0.07062857 0.05305063
## 373 0.0007773607 -0.0007773607 0.04041509 0.05802361 0.04152355
## 374 0.0098342806 0.0098342806 0.03733838 0.05389737 0.03838077
## 375 0.0025143203 0.0025143203 0.03741565 0.05319426 0.03840892
## 376 0.0356413512 0.0356413512 0.03849095 0.05332915 0.03942502
## VaR_FHS_99 VaR_t_95 VaR_t_99
## 371 0.08982202 0.04905702 0.07841371
## 372 0.08772499 0.04944023 0.07701078
## 373 0.07407168 0.03813454 0.06401447
## 374 0.06898892 0.03519375 0.05953115
## 375 0.06757458 0.03537210 0.05856254
## 376 0.06685240 0.03656919 0.05837747
# =========================
# 11. PLOT VaR DENGAN KETERANGAN DAN EXCEEDANCE
# =========================
tema_var <- function() theme_minimal(base_family = "serif") + theme(plot.title = element_text(size = 16, face = "bold", color = "black"), plot.subtitle = element_text(size = 11, color = "black"), axis.text.x = element_text(angle = 45, hjust = 1, size = 11, color = "black"), axis.text.y = element_text(size = 11, color = "black"), axis.title = element_text(size = 14, color = "black"), legend.position = "right", legend.title = element_text(size = 12, face = "bold", color = "black"), legend.text = element_text(size = 11, color = "black"), panel.grid.minor = element_blank(), panel.grid.major = element_line(color = "grey80", linewidth = 0.4), axis.line = element_line(color = "black", linewidth = 0.4))
warna_plot <- c("Actual Loss" = "black", "VaR Normal" = "dodgerblue2", "VaR FHS" = "yellow3", "VaR Student-t" = "green2", "Exceedance" = "darkred")
plot_var_individu <- function(data, var_col, metode, tingkat) {
data_plot <- data %>% mutate(VaR = .data[[var_col]], Exceedance = Loss_Aktual > VaR)
ggplot(data_plot, aes(Date)) + geom_line(aes(y = Loss_Aktual, color = "Actual Loss"), linewidth = 0.5) + geom_line(aes(y = VaR, color = metode), linewidth = 0.7) + geom_point(data = data_plot %>% filter(Exceedance), aes(y = Loss_Aktual, color = "Exceedance"), size = 2.3) + scale_color_manual(name = "Keterangan", values = warna_plot, breaks = c("Actual Loss", metode, "Exceedance")) + scale_x_date(date_breaks = "3 months", date_labels = "%Y-%m") + labs(title = paste0(metode, " ", tingkat), subtitle = "Titik exceedance menunjukkan saat actual loss melebihi nilai VaR", x = "Date", y = "Loss / VaR") + guides(color = guide_legend(title = "Keterangan", override.aes = list(linewidth = c(0.8,0.8,0), shape = c(NA,NA,16), size = c(NA,NA,3)))) + tema_var()
}
plot_var_gabungan <- function(data, tingkat = c("95","99")) {
tingkat <- match.arg(tingkat)
data_long <- if(tingkat == "95") bind_rows(transmute(data, Date, Loss_Aktual, Metode = "VaR Normal", VaR = VaR_Normal_95), transmute(data, Date, Loss_Aktual, Metode = "VaR FHS", VaR = VaR_FHS_95), transmute(data, Date, Loss_Aktual, Metode = "VaR Student-t", VaR = VaR_t_95)) else bind_rows(transmute(data, Date, Loss_Aktual, Metode = "VaR Normal", VaR = VaR_Normal_99), transmute(data, Date, Loss_Aktual, Metode = "VaR FHS", VaR = VaR_FHS_99), transmute(data, Date, Loss_Aktual, Metode = "VaR Student-t", VaR = VaR_t_99))
data_long <- data_long %>% mutate(Exceedance = Loss_Aktual > VaR)
ggplot() + geom_line(data = distinct(data_long, Date, Loss_Aktual), aes(Date, Loss_Aktual, color = "Actual Loss"), linewidth = 0.4) + geom_line(data = data_long, aes(Date, VaR, color = Metode), linewidth = 0.7) + geom_point(data = data_long %>% filter(Exceedance), aes(Date, Loss_Aktual, color = "Exceedance"), size = 2.3) + scale_color_manual(name = "Keterangan", values = warna_plot, breaks = c("Actual Loss","VaR Normal","VaR FHS","VaR Student-t","Exceedance")) + scale_x_date(date_breaks = "3 months", date_labels = "%Y-%m") + labs(title = paste0("Comparison of VaR Estimates at ", tingkat, "% Confidence Level"), subtitle = "Titik exceedance menunjukkan saat actual loss melebihi nilai VaR", x = "Date", y = "Loss / VaR") + guides(color = guide_legend(title = "Keterangan", override.aes = list(linewidth = c(0.8,0.8,0.8,0.8,0), shape = c(NA,NA,NA,NA,16), size = c(NA,NA,NA,NA,3)))) + tema_var()
}
plot_normal_95 <- plot_var_individu(hasil_var, "VaR_Normal_95", "VaR Normal", "95%")
plot_fhs_95 <- plot_var_individu(hasil_var, "VaR_FHS_95", "VaR FHS", "95%")
plot_t_95 <- plot_var_individu(hasil_var, "VaR_t_95", "VaR Student-t", "95%")
plot_normal_99 <- plot_var_individu(hasil_var, "VaR_Normal_99", "VaR Normal", "99%")
plot_fhs_99 <- plot_var_individu(hasil_var, "VaR_FHS_99", "VaR FHS", "99%")
plot_t_99 <- plot_var_individu(hasil_var, "VaR_t_99", "VaR Student-t", "99%")
plot_gabungan_95 <- plot_var_gabungan(hasil_var, "95")
plot_gabungan_99 <- plot_var_gabungan(hasil_var, "99")
plot_normal_95; plot_fhs_95; plot_t_95; plot_normal_99; plot_fhs_99; plot_t_99; plot_gabungan_95; plot_gabungan_99
# =========================
# 12. BACKTESTING VaR: KUPIEC UNCONDITIONAL COVERAGE TEST
# =========================
kupiec_test <- function(loss, var, alpha) {
idx <- is.finite(loss) & is.finite(var); loss <- as.numeric(loss[idx]); var <- as.numeric(var[idx]); violation <- ifelse(loss > var, 1, 0)
n <- length(violation); x <- sum(violation); p_hat <- x / n; eps <- 1e-10; p_hat_adj <- min(max(p_hat, eps), 1 - eps); alpha_adj <- min(max(alpha, eps), 1 - eps)
logL_H0 <- (n - x) * log(1 - alpha_adj) + x * log(alpha_adj); logL_H1 <- (n - x) * log(1 - p_hat_adj) + x * log(p_hat_adj); LR_Kupiec <- -2 * (logL_H0 - logL_H1); p_value <- 1 - pchisq(LR_Kupiec, df = 1)
data.frame(Jumlah_Data = n, Pelanggaran = x, Pelanggaran_Harapan = n * alpha, Rasio_Pelanggaran = p_hat, LR_Kupiec = LR_Kupiec, P_Value_Kupiec = p_value, Keputusan_Kupiec = ifelse(p_value > 0.05, "Lolos Kupiec", "Tidak lolos Kupiec"))
}
buat_backtest_kupiec <- function(data, nama_metode, kolom_var, alpha) data.frame(Metode = nama_metode, Tingkat = paste0((1 - alpha) * 100, "%"), kupiec_test(data$Loss_Aktual, data[[kolom_var]], alpha))
bt_spec <- data.frame(Metode = c("Normal","FHS","Student-t","Normal","FHS","Student-t"), Kolom = c("VaR_Normal_95","VaR_FHS_95","VaR_t_95","VaR_Normal_99","VaR_FHS_99","VaR_t_99"), Alpha = c(alpha_95,alpha_95,alpha_95,alpha_99,alpha_99,alpha_99))
rekap_backtest_kupiec <- bind_rows(lapply(1:nrow(bt_spec), function(i) buat_backtest_kupiec(hasil_var, bt_spec$Metode[i], bt_spec$Kolom[i], bt_spec$Alpha[i]))) %>% arrange(Tingkat, desc(P_Value_Kupiec))
rekap_backtest_kupiec
## Metode Tingkat Jumlah_Data Pelanggaran Pelanggaran_Harapan
## 1 FHS 95% 376 27 18.80
## 2 Normal 95% 376 29 18.80
## 3 Student-t 95% 376 32 18.80
## 4 FHS 99% 376 5 3.76
## 5 Student-t 99% 376 9 3.76
## 6 Normal 99% 376 15 3.76
## Rasio_Pelanggaran LR_Kupiec P_Value_Kupiec Keputusan_Kupiec
## 1 0.07180851 3.3366189 6.775370e-02 Lolos Kupiec
## 2 0.07712766 5.0335378 2.486107e-02 Tidak lolos Kupiec
## 3 0.08510638 8.1341743 4.343863e-03 Tidak lolos Kupiec
## 4 0.01329787 0.3743248 5.406563e-01 Lolos Kupiec
## 5 0.02393617 5.3046129 2.126903e-02 Tidak lolos Kupiec
## 6 0.03989362 19.3718042 1.075838e-05 Tidak lolos Kupiec