1 Pendahuluan

1.1 Deskripsi Penelitian

Penelitian ini berfokus pada analisis risiko pasar dari portofolio saham subsektor perbankan di Indonesia yang terdiri dari lima aset utama:

Kode Saham Emiten
BBCA Bank Central Asia Tbk
BBNI Bank Negara Indonesia Tbk
BBRI Bank Rakyat Indonesia Tbk
BMRI Bank Mandiri Tbk
BRIS Bank Syariah Indonesia Tbk

Fokus utama analisis ini adalah mengukur potensi kerugian ekstrem menggunakan metrik Value at Risk (VaR), bukan untuk mencari komposisi portofolio yang optimal.

1.2 Latar Belakang & Motivasi

Pendekatan konvensional sering kali mengasumsikan bahwa return saham berdistribusi normal. Namun pada kenyataannya, data keuangan memiliki karakteristik:

  • Fat-tail (ekor gemuk): kejadian ekstrem terjadi lebih sering dari prediksi distribusi normal
  • Volatility clustering: periode volatilitas tinggi cenderung berkelompok
  • Leptokurtik: distribusi runcing dengan ekor tebal

Oleh karena itu, pendekatan Extreme Value Theory (EVT) digunakan untuk memodelkan perilaku data pada bagian ekor distribusi (tail behavior) secara lebih akurat, dengan dua metode utama:

  1. Block Maxima → distribusi Generalized Extreme Value (GEV)
  2. Peak Over Threshold (POT) → distribusi Generalized Pareto (GPD)

1.3 Alur Analisis

Data Harga Saham (Yahoo Finance)
         ↓
   Log Return Harian
         ↓
  Statistik Deskriptif & Uji Normalitas
         ↓
  Pembobotan Portofolio → Return Portofolio → Loss
         ↓
  ┌─────────────────────┬──────────────────────┐
  │   Block Maxima       │  Peak Over Threshold │
  │   (GEV Model)        │  (GPD Model)         │
  └────────┬────────────┴──────────┬───────────┘
           ↓                       ↓
     VaR GEV (95%)          VaR GPD (95%)
           ↓                       ↓
         Backtesting (Kupiec Test)
         ↓
   Kesimpulan & Rekomendasi

2 Load Packages & Persiapan Data

Langkah awal adalah memanggil semua library R yang diperlukan untuk manajemen data, visualisasi statistik, serta pemodelan nilai ekstrem. Data harga saham harian (adjusted closing price) diambil secara real-time dari Yahoo Finance mulai 1 Mei 2019 hingga 31 Mei 2025.

library(quantmod)
library(timeDate)
library(writexl)
library(moments)
library(evir)
library(dplyr)
library(lubridate)
library(xts)
library(ggplot2)
library(evmix)
library(tidyr)
library(scales)
library(eva)
library(extRemes)
library(gnFit)
library(knitr)
library(kableExtra)
library(tibble)
# Mengambil data saham dari Yahoo Finance
Saham <- c("BBCA.JK", "BBNI.JK", "BBRI.JK", "BMRI.JK", "BRIS.JK")
getSymbols(Saham, from = "2019-05-01", to = "2025-05-31", auto.assign = TRUE)
## [1] "BBCA.JK" "BBNI.JK" "BBRI.JK" "BMRI.JK" "BRIS.JK"
# Menggabungkan harga penutupan disesuaikan (Adjusted Close)
AdjCloses <- merge(
  Ad(BBCA.JK), Ad(BBNI.JK), Ad(BBRI.JK), Ad(BMRI.JK), Ad(BRIS.JK)
)
colnames(AdjCloses) <- c("BBCA", "BBNI", "BBRI", "BMRI", "BRIS")

# Menghapus data kosong pada hari libur bursa
AdjCloses <- na.omit(AdjCloses[isWeekday(index(AdjCloses)), ])

# Tampilkan beberapa baris pertama
head(AdjCloses) |> kable(caption = "Harga Penutupan Disesuaikan (6 Observasi Pertama)") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Harga Penutupan Disesuaikan (6 Observasi Pertama)
BBCA BBNI BBRI BMRI BRIS
4746.711 3426.591 2564.676 2434.270 512.6641
4693.053 3364.127 2564.676 2426.392 508.0456
4684.797 3301.663 2570.545 2410.636 503.4269
4639.394 3167.812 2482.512 2371.247 498.8083
4672.414 3194.583 2494.250 2410.636 498.8083
4705.436 3123.195 2476.643 2371.247 494.1897

3 Perhitungan Log Return

Harga saham ditransformasikan menjadi log return harian. Penggunaan log return lebih disukai karena:

  • Memiliki sifat aditif secara temporal
  • Menghasilkan distribusi data yang lebih stabil
  • Formula: \(r_t = \ln\left(\dfrac{P_t}{P_{t-1}}\right)\)
LogReturns <- na.omit(merge(
  dailyReturn(Ad(BBCA.JK), type = "log"),
  dailyReturn(Ad(BBNI.JK), type = "log"),
  dailyReturn(Ad(BBRI.JK), type = "log"),
  dailyReturn(Ad(BMRI.JK), type = "log"),
  dailyReturn(Ad(BRIS.JK), type = "log")
))
colnames(LogReturns) <- c("BBCA", "BBNI", "BBRI", "BMRI", "BRIS")
LogReturns <- na.omit(LogReturns)

head(LogReturns) |> kable(caption = "Log Return Harian (6 Observasi Pertama)") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Log Return Harian (6 Observasi Pertama)
BBCA BBNI BBRI BMRI BRIS
0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
-0.011368670 -0.018397212 0.000000000 -0.003241407 -0.009049767
-0.001760715 -0.018742244 0.002285839 -0.006514709 -0.009132535
-0.009738795 -0.041385094 -0.034846681 -0.016474908 -0.009216830
0.007092104 0.008415202 0.004716852 0.016474908 0.000000000
0.007042472 -0.022599873 -0.007083944 -0.016474908 -0.009302323

4 Analisis Deskriptif & Visualisasi Return

4.1 Statistik Deskriptif

Statistik deskriptif memberikan gambaran awal mengenai karakteristik pergerakan return masing-masing saham.

DeskriptifReturnSaham <- data.frame(
  Mean     = apply(LogReturns, 2, mean,     na.rm = TRUE),
  Min      = apply(LogReturns, 2, min,      na.rm = TRUE),
  Max      = apply(LogReturns, 2, max,      na.rm = TRUE),
  SD       = apply(LogReturns, 2, sd,       na.rm = TRUE),
  Variance = apply(LogReturns, 2, var,      na.rm = TRUE),
  Skewness = apply(LogReturns, 2, function(x) moments::skewness(x, na.rm = TRUE)),
  Kurtosis = apply(LogReturns, 2, function(x) moments::kurtosis(x, na.rm = TRUE))
)

DeskriptifReturnSaham <- DeskriptifReturnSaham %>%
  mutate(across(c(Mean, Min, Max, SD, Variance), ~ round(., 6))) %>%
  mutate(across(c(Skewness, Kurtosis), ~ round(., 3)))

DeskriptifReturnSaham |>
  kable(caption = "Statistik Deskriptif Log Return Harian") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = TRUE) |>
  column_spec(7, bold = TRUE, color = "white",
              background = ifelse(DeskriptifReturnSaham$Kurtosis > 3, "#e74c3c", "#27ae60"))
Statistik Deskriptif Log Return Harian
Mean Min Max SD Variance Skewness Kurtosis
BBCA 0.000430 -0.089153 0.159849 0.015871 0.000252 0.571 12.774
BBNI 0.000127 -0.124642 0.127927 0.021943 0.000481 0.151 7.115
BBRI 0.000306 -0.106733 0.186411 0.020970 0.000440 0.482 9.499
BMRI 0.000457 -0.139172 0.146721 0.021503 0.000462 -0.049 8.119
BRIS 0.001186 -0.155485 0.223144 0.035342 0.001249 2.017 14.258

Analisis: Nilai kurtosis yang jauh lebih besar dari 3 pada seluruh saham mengindikasikan karakteristik distribusi yang leptokurtik atau memiliki fat-tail. Hal ini mengonfirmasi bahwa pergerakan return saham perbankan memiliki peluang terjadinya nilai ekstrem yang lebih tinggi daripada distribusi normal standar.

4.2 Visualisasi Time Series Return

Datareturn <- data.frame(Tanggal = index(LogReturns), coredata(LogReturns))
Data_Long <- Datareturn %>%
  pivot_longer(cols = -Tanggal, names_to = "Saham", values_to = "Return")

ggplot(Data_Long, aes(x = Tanggal, y = Return, color = Saham)) +
  geom_line(linewidth = 0.4, alpha = 0.85) +
  facet_wrap(~ Saham, ncol = 2) +
  labs(
    title    = "Pergerakan Log Return Harian Saham Perbankan",
    subtitle = "Periode 1 Mei 2019 – 31 Mei 2025",
    x        = "Tahun",
    y        = "Log Return"
  ) +
  theme_minimal(base_size = 12) +
  theme(
    legend.position  = "none",
    strip.text       = element_text(face = "bold", size = 11),
    plot.title       = element_text(face = "bold", size = 14),
    plot.subtitle    = element_text(color = "gray50"),
    panel.grid.minor = element_blank()
  )

4.3 Distribusi Log Return

ggplot(Data_Long, aes(x = Return)) +
  geom_histogram(
    aes(y = after_stat(density), fill = Saham),
    color = "white", bins = 50, alpha = 0.75
  ) +
  geom_density(color = "midnightblue", linewidth = 0.9) +
  facet_wrap(~ Saham, scales = "free") +
  labs(
    title    = "Distribusi Log Return Harian",
    subtitle = "Histogram dengan Kurva Densitas",
    x        = "Log Return",
    y        = "Density"
  ) +
  theme_minimal(base_size = 12) +
  theme(
    legend.position  = "none",
    strip.text       = element_text(face = "bold", size = 11),
    plot.title       = element_text(face = "bold", size = 14),
    plot.subtitle    = element_text(color = "gray50")
  )

5 Uji Normalitas Formal

Untuk memperkuat indikasi dari nilai kurtosis, dilakukan pengujian formal menggunakan uji Shapiro-Wilk.

  • \(H_0\): Data berdistribusi normal
  • \(H_1\): Data tidak berdistribusi normal
  • Tolak \(H_0\) jika \(p\text{-value} < 0.05\)
UjiNormalitas <- data.frame(
  Saham       = colnames(LogReturns),
  W_Statistic = numeric(ncol(LogReturns)),
  P_Value     = numeric(ncol(LogReturns)),
  Kesimpulan  = character(ncol(LogReturns)),
  stringsAsFactors = FALSE
)

for (i in 1:ncol(LogReturns)) {
  data_saham <- as.numeric(LogReturns[, i])
  uji <- shapiro.test(data_saham)
  UjiNormalitas$W_Statistic[i] <- round(uji$statistic, 4)
  UjiNormalitas$P_Value[i]     <- uji$p.value
  UjiNormalitas$Kesimpulan[i]  <- ifelse(uji$p.value < 0.05,
                                         "Tolak H₀ (Tidak Normal)",
                                         "Gagal Tolak H₀ (Normal)")
}

UjiNormalitas |>
  kable(caption = "Hasil Uji Normalitas Shapiro-Wilk") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE) |>
  column_spec(4, bold = TRUE, color = "white",
              background = ifelse(UjiNormalitas$P_Value < 0.05, "#e74c3c", "#27ae60"))
Hasil Uji Normalitas Shapiro-Wilk
Saham W_Statistic P_Value Kesimpulan
BBCA 0.9260 0 Tolak H₀ (Tidak Normal)
BBNI 0.9469 0 Tolak H₀ (Tidak Normal)
BBRI 0.9399 0 Tolak H₀ (Tidak Normal)
BMRI 0.9426 0 Tolak H₀ (Tidak Normal)
BRIS 0.8228 0 Tolak H₀ (Tidak Normal)

Analisis: Nilai p-value yang sangat kecil (mendekati 0) pada semua saham menolak \(H_0\) pada tingkat signifikansi 5%. Hal ini menegaskan secara statistik bahwa data return tidak berdistribusi normal, sehingga pemodelan nilai ekstrem (EVT) sangat relevan dan valid untuk digunakan.

6 Pembobotan Portofolio & Transformasi Loss

6.1 Pembobotan Berbasis Return Historis

Komposisi bobot investasi pada masing-masing aset ditentukan secara proporsional berdasarkan rata-rata return historis.

meanreturn <- colMeans(LogReturns)
bobot      <- as.numeric(meanreturn / sum(meanreturn))
names(bobot) <- colnames(LogReturns)

# Tampilkan tabel bobot
data.frame(
  Saham = names(bobot),
  `Rata-rata Return` = round(meanreturn, 6),
  `Bobot (%)` = round(bobot * 100, 4)
) |>
  kable(caption = "Bobot Portofolio Berdasarkan Return Historis") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Bobot Portofolio Berdasarkan Return Historis
Saham Rata.rata.Return Bobot….
BBCA BBCA 0.000430 17.1562
BBNI BBNI 0.000127 5.0711
BBRI BBRI 0.000306 12.2204
BMRI BMRI 0.000457 18.2203
BRIS BRIS 0.001186 47.3321

6.2 Return & Transformasi Loss Portofolio

Return portofolio dihitung sebagai kombinasi linear berbobot, kemudian dikalikan \(-1\) untuk mendapatkan nilai loss (kerugian).

\[\text{Loss}_t = -1 \times R_{p,t} = -1 \times \sum_{i=1}^{5} w_i \cdot r_{i,t}\]

Mengapa transformasi ini wajib? EVT secara matematis dirancang untuk menganalisis ekor kanan distribusi (nilai positif besar). Dengan mengalikan return negatif dengan \(-1\), risiko kerugian yang awalnya berada di ekor kiri berpindah ke ekor kanan, sehingga siap dianalisis oleh algoritma EVT.

logreturns <- as.matrix(LogReturns)
Rp         <- as.numeric(logreturns %*% matrix(bobot, ncol = 1))

ReturnPortofolio <- data.frame(
  Tanggal = index(LogReturns),
  Return  = Rp
)
ReturnPortofolio$Loss <- -1 * ReturnPortofolio$Return

head(ReturnPortofolio) |>
  kable(caption = "Return dan Loss Portofolio (6 Observasi Pertama)") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Return dan Loss Portofolio (6 Observasi Pertama)
Tanggal Return Loss
2019-05-01 0.0000000 0.0000000
2019-05-02 -0.0077574 0.0077574
2019-05-03 -0.0064828 0.0064828
2019-05-06 -0.0153922 0.0153922
2019-05-07 0.0052217 -0.0052217
2019-05-08 -0.0082083 0.0082083

7 Identifikasi Nilai Ekstrem (EVT)

7.1 Block Maxima

Data dibagi menjadi sub-blok non-overlapping berukuran 15 observasi, kemudian nilai maksimum setiap blok diambil.

ReturnPortofolio$Tanggal <- as.Date(ReturnPortofolio$Tanggal)

BlockMaxima <- ReturnPortofolio %>%
  arrange(Tanggal) %>%
  mutate(BlockID = ceiling(row_number() / 15)) %>%
  group_by(BlockID) %>%
  slice_max(order_by = Loss, n = 1, with_ties = FALSE) %>%
  ungroup() %>%
  dplyr::select(Tanggal, Loss) %>%
  rename(Max_Loss = Loss)

cat(sprintf("Jumlah blok (Block Maxima): %d blok\n", nrow(BlockMaxima)))
## Jumlah blok (Block Maxima): 99 blok
# Visualisasi Block Maxima
ggplot(BlockMaxima, aes(x = Tanggal, y = Max_Loss)) +
  geom_col(fill = "#2980b9", alpha = 0.8, width = 10) +
  labs(
    title    = "Block Maxima Loss Portofolio",
    subtitle = "Nilai Kerugian Maksimum per Blok 15 Observasi",
    x        = "Tanggal",
    y        = "Maximum Loss"
  ) +
  theme_minimal(base_size = 12) +
  theme(plot.title = element_text(face = "bold"))

7.2 Peak Over Threshold (POT)

Menetapkan threshold \(u\) pada persentil ke-90 (mengambil 10% data kerugian terbesar), lalu mengambil seluruh observasi yang melampaui \(u\).

port_loss   <- ReturnPortofolio$Loss
loss_sorted <- sort(port_loss, decreasing = TRUE)
n_total     <- length(loss_sorted)
jumlah_ekstrem <- floor(0.10 * n_total)
threshold_u <- loss_sorted[jumlah_ekstrem]

cat(sprintf("Total observasi       : %d\n", n_total))
## Total observasi       : 1474
cat(sprintf("Jumlah data ekstrem   : %d (10%%)\n", jumlah_ekstrem))
## Jumlah data ekstrem   : 147 (10%)
cat(sprintf("Threshold (u)         : %.8f\n", threshold_u))
## Threshold (u)         : 0.01933157

7.3 Mean Residual Life Plot

Grafik MRL digunakan untuk memvalidasi secara visual apakah pemilihan threshold \(u\) berada pada area yang linier (stabil).

evmix::mrlplot(
  port_loss,
  main       = "Mean Residual Life Plot",
  xlab       = "Threshold u",
  ylab       = "Mean Excess",
  legend.loc = NULL
)
abline(v = threshold_u, col = "blue", lty = 2, lwd = 1.5)
legend("topright",
  legend   = c("Sample Mean Excess", "95% CI", "Threshold yang Dipilih"),
  col      = c("black", "black", "blue"),
  lty      = c(1, 2, 2),
  cex      = 0.7,
  bg       = "white",
  box.col  = "gray"
)

8 Pemodelan Distribusi GEV dan GPD

8.1 Pemodelan GEV (Block Maxima)

Menggunakan metode Maximum Likelihood Estimation (MLE) untuk mencocokkan data ke distribusi GEV.

gevfit <- extRemes::fevd(BlockMaxima$Max_Loss, type = "GEV", method = "MLE")
summary(gevfit)
## 
## extRemes::fevd(x = BlockMaxima$Max_Loss, type = "GEV", method = "MLE")
## 
## [1] "Estimation Method used: MLE"
## 
## 
##  Negative Log-Likelihood Value:  -272.4126 
## 
## 
##  Estimated parameters:
##   location      scale      shape 
## 0.02093909 0.01203623 0.15646214 
## 
##  Standard Error Estimates:
##    location       scale       shape 
## 0.001377596 0.001053116 0.083253022 
## 
##  Estimated parameter covariance matrix.
##               location         scale         shape
## location  1.897770e-06  7.821009e-07 -3.802284e-05
## scale     7.821009e-07  1.109053e-06 -1.233156e-05
## shape    -3.802284e-05 -1.233156e-05  6.931066e-03
## 
##  AIC = -538.8253 
## 
##  BIC = -531.0399

8.2 Pemodelan GPD (Peak Over Threshold)

gpdfit <- extRemes::fevd(
  ReturnPortofolio$Loss,
  threshold = threshold_u,
  type      = "GP",
  method    = "MLE"
)
summary(gpdfit)
## 
## extRemes::fevd(x = ReturnPortofolio$Loss, threshold = threshold_u, 
##     type = "GP", method = "MLE")
## 
## [1] "Estimation Method used: MLE"
## 
## 
##  Negative Log-Likelihood Value:  -450.8469 
## 
## 
##  Estimated parameters:
##       scale       shape 
##  0.01796278 -0.06833495 
## 
##  Standard Error Estimates:
##       scale       shape 
## 0.001840074 0.062215409 
## 
##  Estimated parameter covariance matrix.
##               scale         shape
## scale  3.385873e-06 -7.684277e-05
## shape -7.684277e-05  3.870757e-03
## 
##  AIC = -897.6937 
## 
##  BIC = -891.7265

8.3 Diagnostik Grafik GEV

par(mfrow = c(1, 2), mar = c(4, 4, 3, 1))
plot(gevfit, type = "density", main = "Density Plot – GEV")
plot(gevfit, type = "qq",      main = "Q-Q Plot – GEV")

8.4 Diagnostik Grafik GPD

par(mfrow = c(1, 2), mar = c(4, 4, 3, 1))
plot(gpdfit, type = "density", main = "Density Plot – GPD")
plot(gpdfit, type = "qq",      main = "Q-Q Plot – GPD")

9 Uji Kesesuaian Distribusi (Goodness of Fit)

Pengujian Anderson-Darling memberikan bobot lebih pada bagian ekor distribusi (tails), sehingga sangat sesuai untuk EVT.

  • \(H_0\): Data mengikuti distribusi yang dispesifikasikan
  • Diterima jika \(p\text{-value} > 0.05\)

9.1 Uji GEV

par_gev <- distill(gevfit)
vec_gev <- c(par_gev["location"], par_gev["scale"], par_gev["shape"])
test_gev_gnfit <- gnFit::gnfit(
  dat  = BlockMaxima$Max_Loss,
  dist = "gev",
  pr   = vec_gev
)

print(test_gev_gnfit)
## $Wpval
## [1] 0.04468385
## 
## $Apval
## [1] 0.06975282
## 
## $Cram
## [1] 0.1295
## 
## $Ander
## [1] 0.6946
## 
## attr(,"class")
## [1] "gnfit"

9.2 Uji GPD

par_gpd <- distill(gpdfit)
vec_gpd <- c(par_gpd["scale"], par_gpd["shape"])
adtestgpd <- gnFit::gnfit(
  dat       = ReturnPortofolio$Loss,
  dist      = "gpd",
  pr        = vec_gpd,
  threshold = threshold_u
)

print(adtestgpd)
## $Wpval
## [1] 0.2206736
## 
## $Apval
## [1] 0.2253517
## 
## $Cram
## [1] 0.0779
## 
## $Ander
## [1] 0.4864
## 
## attr(,"class")
## [1] "gnfit"

10 Estimasi Value at Risk (VaR)

Nilai VaR dihitung pada tingkat kepercayaan 95%, mengindikasikan batas kerugian maksimum portofolio dalam horizon waktu satu hari ke depan dengan peluang sebesar 5%.

10.1 Formula VaR-GEV

\[\text{VaR}_\text{GEV}(p) = \mu - \frac{\sigma}{\xi}\left[1 - (-\ln p)^{-\xi}\right]\]

10.2 Formula VaR-GPD

\[\text{VaR}_\text{GPD}(p) = u + \frac{\sigma}{\xi}\left[\left(\frac{n}{N_u}\cdot p\right)^{-\xi} - 1\right]\]

conf_level <- 0.95

# --- VaR GEV ---
p_gev     <- conf_level
mu_gev    <- par_gev["location"]
sigma_gev <- par_gev["scale"]
xi_gev    <- par_gev["shape"]
term_gev  <- (-log(p_gev))^(-xi_gev)
VaR_GEV   <- mu_gev - (sigma_gev / xi_gev) * (1 - term_gev)

# --- VaR GPD ---
p_gpd   <- 1 - conf_level
sigma   <- par_gpd["scale"]
xi      <- par_gpd["shape"]
u       <- threshold_u
Nu      <- sum(ReturnPortofolio$Loss > u)
n       <- length(ReturnPortofolio$Loss)
VaR_GPD <- u + (sigma / xi) * (((n / Nu) * p_gpd)^(-xi) - 1)

# Tampilkan hasil
tibble(
  Model       = c("VaR-GEV (Block Maxima)", "VaR-GPD (Peak Over Threshold)"),
  `Confidence Level` = c("95%", "95%"),
  `Nilai VaR`  = c(round(VaR_GEV, 6), round(VaR_GPD, 6)),
  Interpretasi = c(
    paste0("Max loss harian: ", round(VaR_GEV * 100, 4), "%"),
    paste0("Max loss harian: ", round(VaR_GPD * 100, 4), "%")
  )
) |>
  kable(caption = "Estimasi Value at Risk (VaR) pada Confidence Level 95%") |>
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE) |>
  row_spec(1:2, bold = TRUE)
Estimasi Value at Risk (VaR) pada Confidence Level 95%
Model Confidence Level Nilai VaR Interpretasi
VaR-GEV (Block Maxima) 95% 0.066447 Max loss harian: 6.6447%
VaR-GPD (Peak Over Threshold) 95% 0.031329 Max loss harian: 3.1329%

11 Backtesting — Validasi Model (Kupiec Test)

Kupiec Test (Proportion of Failures Test) menguji apakah frekuensi pelanggaran aktual konsisten dengan target teoritis kegagalan sebesar 5%.

  • \(H_0\): Tingkat pelanggaran aktual = target teoritis (5%)
  • Model valid jika \(p\text{-value} > 0.05\)

Statistik uji: \[LR_\text{POF} = -2\ln\left[\frac{p^{N_f}(1-p)^{T-N_f}}{\hat{p}^{N_f}(1-\hat{p})^{T-N_f}}\right] \sim \chi^2(1)\]

Actual_Loss <- ReturnPortofolio$Loss
T_Total     <- length(Actual_Loss)
p_expected  <- 1 - conf_level
N_Expected  <- ceiling(T_Total * p_expected)

# --- Backtesting GEV ---
Is_Violation_GEV <- ifelse(Actual_Loss > VaR_GEV, 1, 0)
N_Failures_GEV   <- sum(Is_Violation_GEV)
Failure_Rate_GEV <- N_Failures_GEV / T_Total
term1_gev <- (p_expected^N_Failures_GEV) * ((1 - p_expected)^(T_Total - N_Failures_GEV))
term2_gev <- (Failure_Rate_GEV^N_Failures_GEV) * ((1 - Failure_Rate_GEV)^(T_Total - N_Failures_GEV))
LR_POF_GEV <- -2 * log(term1_gev / term2_gev)
P_Value_GEV <- 1 - pchisq(LR_POF_GEV, df = 1)

# --- Backtesting GPD ---
Is_Violation_GPD <- ifelse(Actual_Loss > VaR_GPD, 1, 0)
N_Failures_GPD   <- sum(Is_Violation_GPD)
Failure_Rate_GPD <- N_Failures_GPD / T_Total
term1_gpd <- (p_expected^N_Failures_GPD) * ((1 - p_expected)^(T_Total - N_Failures_GPD))
term2_gpd <- (Failure_Rate_GPD^N_Failures_GPD) * ((1 - Failure_Rate_GPD)^(T_Total - N_Failures_GPD))
LR_POF_GPD <- -2 * log(term1_gpd / term2_gpd)
P_Value_GPD <- 1 - pchisq(LR_POF_GPD, df = 1)

# --- Tabel Komparasi ---
Tabel_Perbandingan <- tibble(
  Indikator = c(
    "Model VaR",
    "Total Hari Observasi",
    "Target Pelanggaran (5%)",
    "Jumlah Pelanggaran Aktual",
    "Tingkat Pelanggaran Aktual",
    "Nilai VaR Estimasi",
    "Kupiec Likelihood Ratio (LR)",
    "P-Value (Kupiec)",
    "Validitas Model"
  ),
  `VaR-GEV` = c(
    "Generalized Extreme Value",
    as.character(T_Total),
    paste0(N_Expected, " hari"),
    as.character(N_Failures_GEV),
    paste0(sprintf("%.4f", Failure_Rate_GEV * 100), "%"),
    sprintf("%.6f", VaR_GEV),
    sprintf("%.4f", LR_POF_GEV),
    sprintf("%.4f", P_Value_GEV),
    ifelse(P_Value_GEV > 0.05, "✅ Valid", "❌ Tidak Valid")
  ),
  `VaR-GPD` = c(
    "Generalized Pareto Distribution",
    as.character(T_Total),
    paste0(N_Expected, " hari"),
    as.character(N_Failures_GPD),
    paste0(sprintf("%.4f", Failure_Rate_GPD * 100), "%"),
    sprintf("%.6f", VaR_GPD),
    sprintf("%.4f", LR_POF_GPD),
    sprintf("%.4f", P_Value_GPD),
    ifelse(P_Value_GPD > 0.05, "✅ Valid", "❌ Tidak Valid")
  )
)

Tabel_Perbandingan |>
  kable(caption = "Perbandingan Validasi Backtesting Model VaR-EVT") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = TRUE) |>
  row_spec(9, bold = TRUE, background = "#eaf4fb")
Perbandingan Validasi Backtesting Model VaR-EVT
Indikator VaR-GEV VaR-GPD
Model VaR Generalized Extreme Value Generalized Pareto Distribution
Total Hari Observasi 1474 1474
Target Pelanggaran (5%) 74 hari 74 hari
Jumlah Pelanggaran Aktual 7 81
Tingkat Pelanggaran Aktual 0.4749% 5.4953%
Nilai VaR Estimasi 0.066447 0.031329
Kupiec Likelihood Ratio (LR) 103.5705 0.7384
P-Value (Kupiec) 0.0000 0.3902
Validitas Model ❌ Tidak Valid | Valid |

11.1 Visualisasi Pelanggaran VaR

df_bt <- data.frame(
  Tanggal   = ReturnPortofolio$Tanggal,
  Loss      = Actual_Loss,
  Viol_GEV  = Is_Violation_GEV,
  Viol_GPD  = Is_Violation_GPD
)

# Plot GEV
p1 <- ggplot(df_bt, aes(x = Tanggal, y = Loss)) +
  geom_line(color = "gray60", linewidth = 0.4) +
  geom_point(data = df_bt[df_bt$Viol_GEV == 1, ],
             aes(x = Tanggal, y = Loss), color = "#e74c3c", size = 1.5) +
  geom_hline(yintercept = VaR_GEV, color = "#2980b9", linetype = "dashed", linewidth = 0.9) +
  annotate("text", x = min(df_bt$Tanggal), y = VaR_GEV,
           label = paste0("VaR-GEV = ", round(VaR_GEV, 5)),
           hjust = 0, vjust = -0.5, color = "#2980b9", size = 3.2) +
  labs(title = "Backtesting VaR-GEV", subtitle = paste("Pelanggaran:", N_Failures_GEV, "hari"),
       x = NULL, y = "Loss") +
  theme_minimal(base_size = 11) +
  theme(plot.title = element_text(face = "bold"))

# Plot GPD
p2 <- ggplot(df_bt, aes(x = Tanggal, y = Loss)) +
  geom_line(color = "gray60", linewidth = 0.4) +
  geom_point(data = df_bt[df_bt$Viol_GPD == 1, ],
             aes(x = Tanggal, y = Loss), color = "#e67e22", size = 1.5) +
  geom_hline(yintercept = VaR_GPD, color = "#8e44ad", linetype = "dashed", linewidth = 0.9) +
  annotate("text", x = min(df_bt$Tanggal), y = VaR_GPD,
           label = paste0("VaR-GPD = ", round(VaR_GPD, 5)),
           hjust = 0, vjust = -0.5, color = "#8e44ad", size = 3.2) +
  labs(title = "Backtesting VaR-GPD", subtitle = paste("Pelanggaran:", N_Failures_GPD, "hari"),
       x = "Tanggal", y = "Loss") +
  theme_minimal(base_size = 11) +
  theme(plot.title = element_text(face = "bold"))

gridExtra::grid.arrange(p1, p2, ncol = 1)

12 Kesimpulan Akhir

Berdasarkan hasil lengkap analisis di atas:

1. Konfirmasi Non-Normalitas
Seluruh saham perbankan (BBCA, BBNI, BBRI, BMRI, BRIS) memiliki distribusi return yang tidak normal, dengan kurtosis tinggi (leptokurtik) dan terkonfirmasi oleh uji Shapiro-Wilk. Pendekatan EVT terbukti lebih tepat daripada model distribusi normal konvensional.

2. Estimasi VaR
- VaR-GEV (95%): 0.066447 — Dalam 100 hari perdagangan, terdapat kemungkinan 5 hari di mana kerugian portofolio melebihi nilai ini.
- VaR-GPD (95%): 0.031329 — Estimasi serupa dengan memanfaatkan seluruh data ekses di atas threshold.

3. Validasi Backtesting (Kupiec Test)
- Model VaR-GEV: tidak valid secara statistik (p-value = 0.0000)
- Model VaR-GPD: valid dan dapat diandalkan (p-value = 0.3902)

4. Rekomendasi
Pendekatan Peak Over Threshold (GPD) umumnya memberikan estimasi yang lebih dinamis karena memanfaatkan informasi data ekor secara penuh tanpa membuang fluktuasi harian dalam blok waktu tertentu seperti pada metode Block Maxima (GEV). Kedua model ini dapat digunakan sebagai acuan bagi investor dan manajer risiko dalam mengelola eksposur portofolio saham perbankan Indonesia.

Analisis ini dibuat untuk keperluan akademik.
Umar Sodiq | NIM. 22106010001
Mei 2026