pkgs <- c("forecast", "tseries", "FinTS", "rugarch",
"ggplot2", "dplyr", "lubridate", "nortest",
"lmtest", "zoo", "xts", "gridExtra")
for (p in pkgs) {
if (!requireNamespace(p, quietly = TRUE)) install.packages(p)
}
library(forecast)
library(tseries)
library(FinTS)
library(rugarch)
library(ggplot2)
library(dplyr)
library(lubridate)
library(nortest)
library(lmtest)
library(zoo)
library(xts)
library(gridExtra)df <- read.csv("C:\\Users\\aspir\\Documents\\Ekolan\\Data_Historis_GOTO (1).csv",
stringsAsFactors = FALSE)
df$Date <- as.Date(df$Date)
df <- df[order(df$Date), ]
close_ts <- df$Close
n <- length(close_ts)
cat("Jumlah observasi :", n, "\n")## Jumlah observasi : 481
## Periode : 2022-04-11 hingga 2024-04-05
## Statistik deskriptif harga Close:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 56.0 88.0 112.0 157.4 218.0 404.0
ggplot(df, aes(x = Date, y = Close)) +
geom_line(color = "#2c7bb6", linewidth = 0.7) +
labs(title = "Harga Penutupan Saham GOTO (GOTO.JK)",
subtitle = paste("Periode:", format(min(df$Date)), "–", format(max(df$Date))),
x = "Tanggal", y = "Harga (IDR)") +
theme_minimal(base_size = 13) +
theme(plot.title = element_text(face = "bold"))Uji Augmented Dickey-Fuller (ADF) digunakan untuk memeriksa apakah deret waktu bersifat stasioner. Hipotesis nol (H₀): data tidak stasioner (terdapat unit root).
##
## Augmented Dickey-Fuller Test
##
## data: close_ts
## Dickey-Fuller = -2.2508, Lag order = 7, p-value = 0.472
## alternative hypothesis: stationary
##
## -- Keputusan --
if (adf_level$p.value < 0.05) {
cat("H0 ditolak: Data STASIONER pada level (p =",
round(adf_level$p.value, 4), ")\n")
} else {
cat("H0 tidak ditolak: Data TIDAK STASIONER pada level (p =",
round(adf_level$p.value, 4), ") → perlu differencing\n")
}## H0 tidak ditolak: Data TIDAK STASIONER pada level (p = 0.472 ) → perlu differencing
##
## Augmented Dickey-Fuller Test
##
## data: close_diff
## Dickey-Fuller = -7.1166, Lag order = 7, p-value = 0.01
## alternative hypothesis: stationary
##
## -- Keputusan --
if (adf_diff$p.value < 0.05) {
cat("H0 ditolak: Data STASIONER setelah differencing 1 kali (p =",
round(adf_diff$p.value, 4), ") → d = 1\n")
} else {
cat("H0 tidak ditolak: Data masih TIDAK STASIONER (p =",
round(adf_diff$p.value, 4), ") → coba d = 2\n")
}## H0 ditolak: Data STASIONER setelah differencing 1 kali (p = 0.01 ) → d = 1
df_diff <- data.frame(
Date = df$Date[-1],
Diff = close_diff
)
ggplot(df_diff, aes(x = Date, y = Diff)) +
geom_line(color = "#d7191c", linewidth = 0.6) +
geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
labs(title = "Return/Differencing Pertama Harga Penutupan GOTO",
x = "Tanggal", y = "ΔClose (IDR)") +
theme_minimal(base_size = 13)ACF (Autocorrelation Function) dan PACF (Partial Autocorrelation Function) digunakan untuk mengidentifikasi orde AR (p) dan MA (q) yang sesuai pada model ARIMA.
par(mfrow = c(1, 2))
acf(close_ts, lag.max = 12, main = "ACF — Harga Close (Level)")
pacf(close_ts, lag.max = 12, main = "PACF — Harga Close (Level)")par(mfrow = c(1, 2))
acf(close_diff, lag.max = 12, main = "ACF — ΔClose (d=1)")
pacf(close_diff, lag.max = 12, main = "PACF — ΔClose (d=1)")Pembacaan: Berdasarkan pola ACF dan PACF pada data yang sudah di-differencing, kandidat model ARIMA diidentifikasi. Lag yang signifikan menunjukkan kandidat orde AR dan MA.
Beberapa kandidat model ARIMA diestimasi berdasarkan hasil correlogram. Parameter yang disertakan: dengan konstanta (c) dan tanpa konstanta. Model dengan MA lag ke-5 saja [5] membutuhkan parameter fixed.
y <- close_ts # data asli; differencing ditangani oleh order d=1
# ---- Model dengan konstanta (include.constant = TRUE) ----
fit_c110 <- Arima(y, order = c(1, 1, 0), include.constant = TRUE)
fit_c011 <- Arima(y, order = c(0, 1, 1), include.constant = TRUE)
fit_c111 <- Arima(y, order = c(1, 1, 1), include.constant = TRUE)
# ARIMA(0,1,[5]) dengan konstanta: MA hanya pada lag 5
fit_c015 <- Arima(y, order = c(0, 1, 5),
include.constant = TRUE,
fixed = c(0, 0, 0, 0, NA, NA))
# ARIMA(1,1,[5]) dengan konstanta: AR1 + MA hanya pada lag 5
fit_c115 <- Arima(y, order = c(1, 1, 5),
include.constant = TRUE,
fixed = c(NA, 0, 0, 0, 0, NA, NA))
# ---- Model tanpa konstanta ----
fit_110 <- Arima(y, order = c(1, 1, 0), include.constant = FALSE)
fit_011 <- Arima(y, order = c(0, 1, 1), include.constant = FALSE)
fit_111 <- Arima(y, order = c(1, 1, 1), include.constant = FALSE)
fit_015 <- Arima(y, order = c(0, 1, 5),
include.constant = FALSE,
fixed = c(0, 0, 0, 0, NA))
fit_115 <- Arima(y, order = c(1, 1, 5),
include.constant = FALSE,
fixed = c(NA, 0, 0, 0, 0, NA))model_list <- list(
"c ARIMA(1,1,0)" = fit_c110,
"c ARIMA(0,1,1)" = fit_c011,
"c ARIMA(1,1,1)" = fit_c111,
"c ARIMA(0,1,[5])" = fit_c015,
"c ARIMA(1,1,[5])" = fit_c115,
"ARIMA(1,1,0)" = fit_110,
"ARIMA(0,1,1)" = fit_011,
"ARIMA(1,1,1)" = fit_111,
"ARIMA(0,1,[5])" = fit_015,
"ARIMA(1,1,[5])" = fit_115
)
aic_df <- data.frame(
Model = names(model_list),
AIC = sapply(model_list, AIC),
BIC = sapply(model_list, BIC),
LogLik = sapply(model_list, logLik)
)
aic_df <- aic_df[order(aic_df$AIC), ]
rownames(aic_df) <- NULL
aic_df$AIC_per_obs <- aic_df$AIC / n
print(aic_df)## Model AIC BIC LogLik AIC_per_obs
## 1 c ARIMA(1,1,[5]) 3319.835 3336.530 -1655.918 6.901944
## 2 ARIMA(1,1,[5]) 3319.844 3332.366 -1656.922 6.901963
## 3 c ARIMA(1,1,0) 3322.908 3335.430 -1658.454 6.908333
## 4 c ARIMA(0,1,1) 3322.957 3335.478 -1658.478 6.908434
## 5 ARIMA(1,1,0) 3323.286 3331.634 -1659.643 6.909118
## 6 ARIMA(0,1,1) 3323.494 3331.842 -1659.747 6.909552
## 7 c ARIMA(1,1,1) 3324.862 3341.557 -1658.431 6.912395
## 8 ARIMA(1,1,1) 3325.283 3337.805 -1659.642 6.913271
## 9 c ARIMA(0,1,[5]) 3333.567 3346.088 -1663.783 6.930492
## 10 ARIMA(0,1,[5]) 3334.428 3342.776 -1665.214 6.932283
##
## ✔ Model terbaik berdasarkan AIC: c ARIMA(1,1,[5])
## AIC = 3319.835
ggplot(aic_df, aes(x = reorder(Model, AIC), y = AIC,
fill = (Model == best_model_name))) +
geom_col(show.legend = FALSE) +
scale_fill_manual(values = c("steelblue", "#e63946")) +
coord_flip() +
labs(title = "Perbandingan AIC Kandidat Model ARIMA",
subtitle = paste("Model terbaik:", best_model_name),
x = NULL, y = "AIC") +
theme_minimal(base_size = 12)# Gunakan model dengan AIC terkecil secara dinamis — tidak hardcode
best_fit <- model_list[[best_model_name]]
cat("=== Ringkasan Model Terbaik:", best_model_name, "===\n")## === Ringkasan Model Terbaik: c ARIMA(1,1,[5]) ===
## Series: y
## ARIMA(1,1,5) with drift
##
## Coefficients:
## ar1 ma1 ma2 ma3 ma4 ma5 drift
## 0.1797 0 0 0 0 0.0981 -0.6626
## s.e. 0.0449 0 0 0 0 0.0429 0.4650
##
## sigma^2 = 58.43: log likelihood = -1655.92
## AIC=3319.84 AICc=3319.92 BIC=3336.53
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.007967425 7.612315 4.740044 0.07831046 3.106406 0.9773287
## ACF1
## Training set 0.0004522709
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 0.179675 0.044917 4.0002 6.33e-05 ***
## ma5 0.098117 0.042902 2.2870 0.02219 *
## drift -0.662595 0.464956 -1.4251 0.15414
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Simpan orde ARMA untuk digunakan pada spesifikasi mean model ARCH
p_arma <- best_fit$arma[1] # AR order
q_arma <- best_fit$arma[2] # MA order
cat("\nOrde ARMA yang akan digunakan di mean model ARCH:",
"AR =", p_arma, "| MA =", q_arma, "\n")##
## Orde ARMA yang akan digunakan di mean model ARCH: AR = 1 | MA = 5
resid_arima <- residuals(best_fit)
# Shapiro-Wilk (untuk n ≤ 5000)
if (length(resid_arima) <= 5000) {
sw_test <- shapiro.test(resid_arima)
cat("Shapiro-Wilk Test:\n")
print(sw_test)
cat("Keputusan:", ifelse(sw_test$p.value < 0.05,
"H0 ditolak — Residual TIDAK berdistribusi normal\n",
"H0 tidak ditolak — Residual berdistribusi normal\n"))
}## Shapiro-Wilk Test:
##
## Shapiro-Wilk normality test
##
## data: resid_arima
## W = 0.86077, p-value < 2.2e-16
##
## Keputusan: H0 ditolak — Residual TIDAK berdistribusi normal
# Kolmogorov-Smirnov
ks_test <- ks.test(scale(resid_arima), "pnorm")
cat("\nKolmogorov-Smirnov Test:\n")##
## Kolmogorov-Smirnov Test:
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: scale(resid_arima)
## D = 0.11948, p-value = 2.17e-06
## alternative hypothesis: two-sided
##
## Jarque-Bera Test:
##
## Jarque Bera Test
##
## data: resid_arima
## X-squared = 1692, df = 2, p-value < 2.2e-16
cat("Keputusan:", ifelse(jb_test$p.value < 0.05,
"H0 ditolak — Residual TIDAK normal\n",
"H0 tidak ditolak — Residual normal\n"))## Keputusan: H0 ditolak — Residual TIDAK normal
par(mfrow = c(1, 2))
hist(resid_arima, breaks = 30, probability = TRUE,
main = "Histogram Residual", xlab = "Residual",
col = "#cce5ff", border = "white")
curve(dnorm(x, mean(resid_arima), sd(resid_arima)),
add = TRUE, col = "red", lwd = 2)
qqnorm(resid_arima, main = "Q-Q Plot Residual")
qqline(resid_arima, col = "red", lwd = 2)## === Uji Ljung-Box (White Noise) ===
for (lag in c(1, 6, 12, 24)) {
lb <- Box.test(resid_arima, lag = lag, type = "Ljung-Box", fitdf = 2)
cat(sprintf("Lag %-3d: Q = %.4f, p-value = %.4f → %s\n",
lag, lb$statistic, lb$p.value,
ifelse(lb$p.value > 0.05,
"Tidak ada autokorelasi (White Noise ✓)",
"Terdapat autokorelasi (bukan White Noise ✗)")))
}## Lag 1 : Q = 0.0001, p-value = NaN → NA
## Lag 6 : Q = 3.9172, p-value = 0.4173 → Tidak ada autokorelasi (White Noise ✓)
## Lag 12 : Q = 10.0066, p-value = 0.4399 → Tidak ada autokorelasi (White Noise ✓)
## Lag 24 : Q = 32.2665, p-value = 0.0730 → Tidak ada autokorelasi (White Noise ✓)
par(mfrow = c(1, 2))
acf(resid_arima, lag.max = 30, main = "ACF Residual")
pacf(resid_arima, lag.max = 30, main = "PACF Residual")Uji ARCH-LM (Lagrange Multiplier) digunakan untuk mendeteksi adanya heteroskedastisitas kondisional (efek ARCH) pada residual model ARIMA.
## === Uji ARCH-LM pada Residual ARIMA ===
for (lag in 1:12) {
arch_lm <- ArchTest(resid_arima, lags = lag)
cat(sprintf("Lag %-3d: Chi² = %.4f, p-value = %.4f → %s\n",
lag, arch_lm$statistic, arch_lm$p.value,
ifelse(arch_lm$p.value < 0.05,
"Terdapat efek ARCH ✗",
"Tidak ada efek ARCH ✓")))
}## Lag 1 : Chi² = 29.5696, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 2 : Chi² = 30.7980, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 3 : Chi² = 42.5202, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 4 : Chi² = 44.1316, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 5 : Chi² = 44.0574, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 6 : Chi² = 44.7274, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 7 : Chi² = 46.3282, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 8 : Chi² = 52.0703, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 9 : Chi² = 78.4336, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 10 : Chi² = 78.6379, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 11 : Chi² = 85.7399, p-value = 0.0000 → Terdapat efek ARCH ✗
## Lag 12 : Chi² = 87.5614, p-value = 0.0000 → Terdapat efek ARCH ✗
Karena terdapat efek ARCH pada residual, dilakukan pemodelan volatilitas menggunakan model ARCH. Pemilihan orde terbaik mempertimbangkan dua kriteria: (1) nilai AIC terkecil dan (2) signifikansi statistik semua parameter alpha. Model yang lolos keduanya diprioritaskan.
# Inisialisasi vektor penyimpan hasil
aic_arch <- rep(NA, 5)
sig_arch <- rep(FALSE, 5) # TRUE jika semua alpha signifikan (p < 0.05)
# Mean model mengikuti orde ARMA dari model ARIMA terbaik
cat(sprintf("Mean model yang digunakan: ARMA(%d, %d)\n\n", p_arma, q_arma))## Mean model yang digunakan: ARMA(1, 5)
for (q in 1:5) {
spec <- ugarchspec(
variance.model = list(model = "sGARCH", garchOrder = c(q, 0)),
mean.model = list(armaOrder = c(p_arma, q_arma),
include.mean = TRUE),
distribution.model = "norm"
)
tryCatch({
fit_arch_q <- ugarchfit(spec = spec, data = close_ts,
solver = "hybrid", out.sample = 0)
# AIC (per observasi, skala rugarch)
ic <- infocriteria(fit_arch_q)
aic_arch[q] <- ic["Akaike", ]
# Signifikansi parameter alpha1, ..., alphaq
coef_q <- coef(fit_arch_q)
se_q <- sqrt(diag(vcov(fit_arch_q)))
pval_q <- 2 * pnorm(abs(coef_q / se_q), lower.tail = FALSE)
alpha_names <- paste0("alpha", 1:q)
# Ambil hanya alpha yang benar-benar ada di koefisien
alpha_names <- alpha_names[alpha_names %in% names(pval_q)]
alpha_pval <- pval_q[alpha_names]
all_sig <- all(alpha_pval < 0.05, na.rm = TRUE)
sig_arch[q] <- all_sig
cat(sprintf(
"ARCH(%d): AIC = %.6f | Alpha signifikan: %s\n",
q, aic_arch[q],
ifelse(all_sig,
paste0("Semua signifikan ✓ (p: ",
paste(round(alpha_pval, 4), collapse = ", "), ")"),
paste0("Ada yang tidak signifikan ✗ (p: ",
paste(round(alpha_pval, 4), collapse = ", "), ")"))
))
}, error = function(e) {
cat(sprintf("ARCH(%d): Gagal diestimasi — %s\n", q, e$message))
aic_arch[q] <<- NA
sig_arch[q] <<- FALSE
})
}## ARCH(1): AIC = 6.656117 | Alpha signifikan: Semua signifikan ✓ (p: 0)
## ARCH(2): AIC = 6.650132 | Alpha signifikan: Semua signifikan ✓ (p: 0, NaN)
## ARCH(3): AIC = 6.491010 | Alpha signifikan: Ada yang tidak signifikan ✗ (p: 0, 0.2365, 0)
## ARCH(4): AIC = 6.492193 | Alpha signifikan: Ada yang tidak signifikan ✗ (p: 0, 0.2834, 0, 0.3218)
## ARCH(5): AIC = 6.458745 | Alpha signifikan: Ada yang tidak signifikan ✗ (p: 0, 0.2419, 0, 1, 0.0044)
# Tabel ringkasan perbandingan seluruh orde
arch_crit_df <- data.frame(
Orde = paste0("ARCH(", 1:5, ")"),
AIC = round(aic_arch, 6),
Semua_Alpha_Sig = ifelse(sig_arch, "Ya ✓", "Tidak ✗"),
Keterangan = ifelse(
is.na(aic_arch), "Gagal estimasi",
ifelse(sig_arch, "Memenuhi kedua kriteria", "Hanya lolos AIC")
)
)
print(arch_crit_df)## Orde AIC Semua_Alpha_Sig Keterangan
## 1 ARCH(1) 6.656117 Ya ✓ Memenuhi kedua kriteria
## 2 ARCH(2) 6.650132 Ya ✓ Memenuhi kedua kriteria
## 3 ARCH(3) 6.491010 Tidak ✗ Hanya lolos AIC
## 4 ARCH(4) 6.492193 Tidak ✗ Hanya lolos AIC
## 5 ARCH(5) 6.458745 Tidak ✗ Hanya lolos AIC
# Pilih: AIC terkecil di antara model yang semua alpha-nya signifikan
valid_idx <- which(!is.na(aic_arch) & sig_arch)
if (length(valid_idx) > 0) {
best_q <- valid_idx[which.min(aic_arch[valid_idx])]
cat(sprintf(
"\n✔ Model ARCH terbaik (AIC terkecil & semua alpha signifikan): ARCH(%d)\n",
best_q))
cat(" AIC =", round(aic_arch[best_q], 6), "\n")
cat(" Catatan: model ini memenuhi kriteria AIC sekaligus signifikansi parameter.\n")
} else {
# Fallback: pilih berdasarkan AIC saja dan beri peringatan
best_q <- which.min(aic_arch)
cat(sprintf(
"\n⚠ Tidak ada model dengan semua alpha signifikan.\n Dipilih berdasarkan AIC terkecil saja: ARCH(%d)\n",
best_q))
cat(" AIC =", round(aic_arch[best_q], 6), "\n")
}##
## ✔ Model ARCH terbaik (AIC terkecil & semua alpha signifikan): ARCH(2)
## AIC = 6.650132
## Catatan: model ini memenuhi kriteria AIC sekaligus signifikansi parameter.
# Warnai: merah = terpilih, biru = valid tapi bukan terbaik,
# abu = tidak memenuhi signifikansi
arch_df <- data.frame(
Orde = paste0("ARCH(", 1:5, ")"),
AIC = aic_arch,
Status = ifelse(
is.na(aic_arch), "Gagal",
ifelse(1:5 == best_q, "Terpilih",
ifelse(sig_arch, "Memenuhi", "Tidak Signifikan"))
)
)
ggplot(arch_df[!is.na(arch_df$AIC), ],
aes(x = Orde, y = AIC, fill = Status)) +
geom_col(show.legend = TRUE) +
scale_fill_manual(
values = c("Terpilih" = "#e63946",
"Memenuhi" = "#2a9d8f",
"Tidak Signifikan" = "steelblue",
"Gagal" = "gray70")) +
labs(title = "Perbandingan AIC Model ARCH(1)–ARCH(5)",
subtitle = paste0("Model terpilih: ARCH(", best_q,
") — AIC terkecil & semua alpha signifikan"),
x = NULL, y = "AIC (per obs)", fill = "Status") +
theme_minimal(base_size = 12)# Spesifikasi menggunakan orde terpilih (best_q) secara dinamis
spec_best <- ugarchspec(
variance.model = list(model = "sGARCH", garchOrder = c(best_q, 0)),
mean.model = list(armaOrder = c(p_arma, q_arma),
include.mean = TRUE),
distribution.model = "norm"
)
fit_best_arch <- ugarchfit(spec = spec_best, data = close_ts,
solver = "hybrid")
cat(sprintf("=== Model Terpilih: %s-ARCH(%d) ===\n\n",
best_model_name, best_q))## === Model Terpilih: c ARIMA(1,1,[5])-ARCH(2) ===
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(2,0)
## Mean Model : ARFIMA(1,0,5)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 373.590626 18.874356 19.79356 0.000000
## ar1 1.000000 0.000975 1025.23198 0.000000
## ma1 0.281670 0.041367 6.80901 0.000000
## ma2 0.181271 0.141182 1.28395 0.199159
## ma3 -0.085357 0.018641 -4.57888 0.000005
## ma4 0.011037 0.028132 0.39232 0.694820
## ma5 0.018698 0.080655 0.23183 0.816673
## omega 16.078988 10.761944 1.49406 0.135160
## alpha1 0.487211 0.055800 8.73143 0.000000
## alpha2 0.511788 0.402208 1.27245 0.203215
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 373.590626 248.908404 1.500916 0.13338
## ar1 1.000000 0.005933 168.558712 0.00000
## ma1 0.281670 0.432208 0.651699 0.51460
## ma2 0.181271 1.837221 0.098666 0.92140
## ma3 -0.085357 0.637731 -0.133845 0.89353
## ma4 0.011037 0.583793 0.018906 0.98492
## ma5 0.018698 1.116192 0.016752 0.98664
## omega 16.078988 134.895126 0.119196 0.90512
## alpha1 0.487211 0.905995 0.537764 0.59074
## alpha2 0.511788 5.054367 0.101257 0.91935
##
## LogLikelihood : -1589.357
##
## Information Criteria
## ------------------------------------
##
## Akaike 6.6501
## Bayes 6.7369
## Shibata 6.6493
## Hannan-Quinn 6.6843
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 3.377 6.611e-02
## Lag[2*(p+q)+(p+q)-1][17] 12.928 2.306e-09
## Lag[4*(p+q)+(p+q)-1][29] 18.022 1.477e-01
## d.o.f=6
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.01556 0.9007
## Lag[2*(p+q)+(p+q)-1][5] 2.83118 0.4387
## Lag[4*(p+q)+(p+q)-1][9] 4.47556 0.5102
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 2.703 0.500 2.000 0.1002
## ARCH Lag[5] 2.775 1.440 1.667 0.3241
## ARCH Lag[7] 3.274 2.315 1.543 0.4634
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 4.8799
## Individual Statistics:
## mu 0.005941
## ar1 0.121712
## ma1 0.506718
## ma2 0.059753
## ma3 0.243392
## ma4 0.309561
## ma5 0.197394
## omega 2.905915
## alpha1 1.206375
## alpha2 0.560335
##
## 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 1.2711 0.2043
## Negative Sign Bias 0.6080 0.5435
## Positive Sign Bias 0.9334 0.3511
## Joint Effect 1.8301 0.6084
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 72.22 3.914e-08
## 2 30 85.88 1.546e-07
## 3 40 105.36 5.132e-08
## 4 50 103.93 7.968e-06
##
##
## Elapsed time : 0.5622752
coef_arch <- coef(fit_best_arch)
se_arch <- sqrt(diag(vcov(fit_best_arch)))
z_arch <- coef_arch / se_arch
p_arch <- 2 * pnorm(abs(z_arch), lower.tail = FALSE)
sig_df <- data.frame(
Parameter = names(coef_arch),
Estimasi = round(coef_arch, 6),
Std_Error = round(se_arch, 6),
z_value = round(z_arch, 4),
p_value = round(p_arch, 4),
Signifikan = ifelse(p_arch < 0.05, "Ya ✓", "Tidak ✗")
)
cat(sprintf("=== Signifikansi Parameter ARCH(%d) ===\n\n", best_q))## === Signifikansi Parameter ARCH(2) ===
## Parameter Estimasi Std_Error z_value p_value Signifikan
## mu mu 373.590626 NaN NaN NaN <NA>
## ar1 ar1 1.000000 0.000975 1025.2320 0 Ya ✓
## ma1 ma1 0.281670 0.041367 6.8090 0 Ya ✓
## ma2 ma2 0.181271 NaN NaN NaN <NA>
## ma3 ma3 -0.085357 NaN NaN NaN <NA>
## ma4 ma4 0.011037 NaN NaN NaN <NA>
## ma5 ma5 0.018698 NaN NaN NaN <NA>
## omega omega 16.078988 NaN NaN NaN <NA>
## alpha1 alpha1 0.487211 0.055800 8.7314 0 Ya ✓
## alpha2 alpha2 0.511788 NaN NaN NaN <NA>
n_sig <- sum(p_arch < 0.05)
n_total <- length(p_arch)
alpha_idx <- grep("^alpha", names(p_arch))
n_alpha_sig <- sum(p_arch[alpha_idx] < 0.05)
cat(sprintf(
"\nRingkasan: %d dari %d parameter signifikan pada α = 5%%.\n",
n_sig, n_total))##
## Ringkasan: NA dari 10 parameter signifikan pada α = 5%.
## Parameter alpha (ARCH): NA dari 2 signifikan.
resid_arch <- residuals(fit_best_arch, standardize = TRUE)
cat(sprintf("=== Uji ARCH-LM pada Residual %s-ARCH(%d) ===\n\n",
best_model_name, best_q))## === Uji ARCH-LM pada Residual c ARIMA(1,1,[5])-ARCH(2) ===
for (lag in c(1, 5, 10, 15)) {
arch_lm2 <- ArchTest(resid_arch, lags = lag)
cat(sprintf("Lag %-3d: Chi² = %.4f, p-value = %.4f → %s\n",
lag, arch_lm2$statistic, arch_lm2$p.value,
ifelse(arch_lm2$p.value > 0.05,
"Tidak ada efek ARCH (model cukup ✓)",
"Masih terdapat efek ARCH ✗")))
}## Lag 1 : Chi² = 0.0154, p-value = 0.9011 → Tidak ada efek ARCH (model cukup ✓)
## Lag 5 : Chi² = 4.1271, p-value = 0.5313 → Tidak ada efek ARCH (model cukup ✓)
## Lag 10 : Chi² = 17.8188, p-value = 0.0581 → Tidak ada efek ARCH (model cukup ✓)
## Lag 15 : Chi² = 21.4241, p-value = 0.1238 → Tidak ada efek ARCH (model cukup ✓)
sigma_arch <- sigma(fit_best_arch)
vol_df <- data.frame(
Date = df$Date,
Close = close_ts,
Sigma = sigma_arch
)
p1 <- ggplot(vol_df, aes(x = Date, y = Close)) +
geom_line(color = "#2c7bb6") +
labs(title = "Harga Penutupan GOTO", x = NULL, y = "Harga (IDR)") +
theme_minimal(base_size = 11)
p2 <- ggplot(vol_df, aes(x = Date, y = Sigma)) +
geom_line(color = "#e63946") +
labs(title = paste0("Estimasi Volatilitas Kondisional — ARCH(", best_q, ")"),
x = "Tanggal", y = "σ_t (Std Dev Kondisional)") +
theme_minimal(base_size = 11)
grid.arrange(p1, p2, ncol = 1)Dilakukan peramalan harga saham GOTO menggunakan model gabungan terbaik untuk 30 hari ke depan.
h <- 30 # horizon peramalan (hari kerja)
fc_arima <- forecast(best_fit, h = h, level = c(80, 95))
autoplot(fc_arima) +
labs(title = paste0("Forecast Harga GOTO — ", best_model_name),
subtitle = paste("Horizon:", h, "hari ke depan"),
x = "Waktu", y = "Harga (IDR)") +
theme_minimal(base_size = 12)tanggal_fc <- seq(max(df$Date) + 1, by = "day", length.out = h * 2)
tanggal_fc <- tanggal_fc[!weekdays(tanggal_fc) %in% c("Saturday", "Sunday")][1:h]
fc_df <- data.frame(
Tanggal = tanggal_fc,
Forecast = round(as.numeric(fc_arima$mean), 2),
Lower_80 = round(as.numeric(fc_arima$lower[, 1]), 2),
Upper_80 = round(as.numeric(fc_arima$upper[, 1]), 2),
Lower_95 = round(as.numeric(fc_arima$lower[, 2]), 2),
Upper_95 = round(as.numeric(fc_arima$upper[, 2]), 2)
)
print(fc_df)## Tanggal Forecast Lower_80 Upper_80 Lower_95 Upper_95
## 1 2024-04-08 67.17 57.37 76.96 52.19 82.15
## 2 2024-04-09 66.74 51.59 81.89 43.57 89.91
## 3 2024-04-10 65.94 46.69 85.18 36.50 95.37
## 4 2024-04-11 65.44 42.80 88.09 30.81 100.08
## 5 2024-04-12 64.71 39.11 90.31 25.56 103.87
## 6 2024-04-15 64.04 35.37 92.71 20.20 107.88
## 7 2024-04-16 63.38 31.87 94.88 15.19 111.56
## 8 2024-04-17 62.71 28.59 96.84 10.52 114.90
## 9 2024-04-18 62.05 25.49 98.61 6.14 117.96
## 10 2024-04-19 61.39 22.55 100.23 1.99 120.79
## 11 2024-04-22 60.72 19.73 101.72 -1.97 123.42
## 12 2024-04-23 60.06 17.02 103.10 -5.76 125.89
## 13 2024-04-24 59.40 14.41 104.39 -9.41 128.21
## 14 2024-04-25 58.74 11.87 105.60 -12.94 130.41
## 15 2024-04-26 58.07 9.41 106.74 -16.35 132.50
## 16 2024-04-29 57.41 7.01 107.81 -19.67 134.49
## 17 2024-04-30 56.75 4.67 108.83 -22.90 136.40
## 18 2024-05-01 56.09 2.38 109.79 -26.05 138.22
## 19 2024-05-02 55.42 0.14 110.71 -29.13 139.97
## 20 2024-05-03 54.76 -2.06 111.58 -32.13 141.66
## 21 2024-05-06 54.10 -4.21 112.41 -35.08 143.28
## 22 2024-05-07 53.44 -6.33 113.20 -37.97 144.84
## 23 2024-05-08 52.77 -8.42 113.96 -40.81 146.35
## 24 2024-05-09 52.11 -10.47 114.69 -43.60 147.82
## 25 2024-05-10 51.45 -12.49 115.39 -46.34 149.23
## 26 2024-05-13 50.79 -14.48 116.05 -49.04 150.61
## 27 2024-05-14 50.12 -16.45 116.70 -51.69 151.94
## 28 2024-05-15 49.46 -18.39 117.31 -54.31 153.23
## 29 2024-05-16 48.80 -20.31 117.91 -56.89 154.49
## 30 2024-05-17 48.13 -22.21 118.48 -59.44 155.71
fc_arch <- ugarchforecast(fit_best_arch, n.ahead = h)
sigma_fc <- as.numeric(sigma(fc_arch))
mean_fc <- as.numeric(fitted(fc_arch))
fc_vol_df <- data.frame(
Hari = 1:h,
Tanggal = tanggal_fc,
Mean_FC = round(mean_fc, 4),
Sigma_FC = round(sigma_fc, 4),
Lower_95 = round(mean_fc - 1.96 * sigma_fc, 4),
Upper_95 = round(mean_fc + 1.96 * sigma_fc, 4)
)
print(fc_vol_df)## Hari Tanggal Mean_FC Sigma_FC Lower_95 Upper_95
## 1 1 2024-04-08 68.1168 4.1302 60.0216 76.2120
## 2 2 2024-04-09 67.8786 4.9762 58.1252 77.6320
## 3 3 2024-04-10 67.9273 6.0724 56.0254 79.8293
## 4 4 2024-04-11 67.9386 6.8351 54.5418 81.3353
## 5 5 2024-04-12 67.9226 7.5969 53.0327 82.8124
## 6 6 2024-04-15 67.9226 8.2527 51.7473 84.0978
## 7 7 2024-04-16 67.9226 8.8768 50.5240 85.3211
## 8 8 2024-04-17 67.9226 9.4513 49.3981 86.4471
## 9 9 2024-04-18 67.9226 9.9964 48.3297 87.5155
## 10 10 2024-04-19 67.9226 10.5110 47.3210 88.5241
## 11 11 2024-04-22 67.9226 11.0022 46.3583 89.4869
## 12 12 2024-04-23 67.9226 11.4716 45.4382 90.4069
## 13 13 2024-04-24 67.9226 11.9225 44.5545 91.2907
## 14 14 2024-04-25 67.9226 12.3566 43.7037 92.1414
## 15 15 2024-04-26 67.9226 12.7756 42.8823 92.9628
## 16 16 2024-04-29 67.9226 13.1811 42.0875 93.7576
## 17 17 2024-04-30 67.9226 13.5743 41.3170 94.5281
## 18 18 2024-05-01 67.9226 13.9561 40.5686 95.2765
## 19 19 2024-05-02 67.9225 14.3275 39.8407 96.0044
## 20 20 2024-05-03 67.9225 14.6893 39.1316 96.7135
## 21 21 2024-05-06 67.9225 15.0421 38.4400 97.4051
## 22 22 2024-05-07 67.9225 15.3866 37.7647 98.0803
## 23 23 2024-05-08 67.9225 15.7234 37.1047 98.7404
## 24 24 2024-05-09 67.9225 16.0529 36.4589 99.3862
## 25 25 2024-05-10 67.9225 16.3755 35.8265 100.0185
## 26 26 2024-05-13 67.9225 16.6917 35.2068 100.6383
## 27 27 2024-05-14 67.9225 17.0018 34.5989 101.2461
## 28 28 2024-05-15 67.9225 17.3062 34.0024 101.8427
## 29 29 2024-05-16 67.9225 17.6051 33.4165 102.4285
## 30 30 2024-05-17 67.9225 17.8988 32.8408 103.0042
hist_df <- data.frame(
Tanggal = df$Date,
Close = close_ts,
Type = "Aktual"
)
fc_plot_df <- data.frame(
Tanggal = tanggal_fc,
Mean = mean_fc,
Lower_95 = mean_fc - 1.96 * sigma_fc,
Upper_95 = mean_fc + 1.96 * sigma_fc
)
hist_recent <- tail(hist_df, 60)
ggplot() +
geom_line(data = hist_recent,
aes(x = Tanggal, y = Close),
color = "#2c7bb6", linewidth = 0.8) +
geom_ribbon(data = fc_plot_df,
aes(x = Tanggal, ymin = Lower_95, ymax = Upper_95),
fill = "#e63946", alpha = 0.2) +
geom_line(data = fc_plot_df,
aes(x = Tanggal, y = Mean),
color = "#e63946", linewidth = 0.9, linetype = "dashed") +
geom_vline(xintercept = as.numeric(max(df$Date)),
linetype = "dotted", color = "gray50") +
annotate("text",
x = max(df$Date) + 2,
y = max(hist_recent$Close),
label = "← Forecast", hjust = 0,
color = "#e63946", size = 3.5) +
labs(title = paste0("Forecast Harga Penutupan GOTO — ",
best_model_name, "-ARCH(", best_q, ")"),
subtitle = paste("Horizon:", h, "hari kerja | Pita: 95% CI"),
x = "Tanggal", y = "Harga (IDR)") +
theme_minimal(base_size = 12)## ============================================================
## RINGKASAN ANALISIS DERET WAKTU SAHAM GOTO
## ============================================================
## [1] Stasioneritas:
## • Data level tidak stasioner → differencing 1 kali (d = 1)
## • ADF setelah differencing: p = 0.01 → stasioner
## [2] Model ARIMA Terbaik: c ARIMA(1,1,[5])
## • Dipilih berdasarkan AIC terkecil di antara 10 kandidat model
## • AIC = 3319.835
## [3] Diagnostik Residual ARIMA:
## • Uji White Noise (Ljung-Box) : lihat output Bagian 5b
## • Uji Normalitas (JB/SW) : lihat output Bagian 5a
## • Efek ARCH terdeteksi : perlu pemodelan volatilitas
## [4] Pemilihan Model ARCH:
## • Kandidat: ARCH(1) hingga ARCH(5)
## • Kriteria: AIC terkecil + semua parameter alpha signifikan
valid_str <- if (length(valid_idx) > 0) {
paste0("ARCH(", valid_idx, ")", collapse = ", ")
} else "tidak ada"
cat(" • Orde yang memenuhi kedua kriteria:", valid_str, "\n")## • Orde yang memenuhi kedua kriteria: ARCH(1), ARCH(2)
## • Model terpilih: ARCH(2) — AIC = 6.650132
## [5] Diagnostik Residual ARCH:
## • Uji ARCH-LM pada residual standar ARCH(2): lihat Bagian 6
## • Jika semua lag p > 0.05 → model volatilitas sudah memadai
## [6] Forecast: 30 hari ke depan
## • Model gabungan: c ARIMA(1,1,[5])-ARCH(2)
## • Interval kepercayaan 95% disertakan
## ============================================================
Analisis dilakukan dengan R. Paket utama: forecast,
tseries, rugarch, FinTS.