PerformanceAnalytics 包

# 使用PerformanceAnalytics包计算ES

VaR_historical <- PerformanceAnalytics::VaR(returns, p = alpha , method = "historical")
VaR_gaussian <- PerformanceAnalytics::VaR(returns, p = alpha , method = "gaussian")
es_historical <- PerformanceAnalytics::ES(returns, p = alpha, method = "historical")
es_gaussian <- PerformanceAnalytics::ES(returns, p = alpha, method = "gaussian")

cat("PerformanceAnalytics包计算结果:\n")

## PerformanceAnalytics包计算结果:

cat("历史模拟法 VaR:", round(VaR_historical * 100, 2), "%\n")

## 历史模拟法 VaR: -2.58 %

cat("正态分布假设 VaR:", round(VaR_gaussian * 100, 2), "%\n")

## 正态分布假设 VaR: -2.64 %

cat("历史模拟法 ES:", round(es_historical * 100, 2), "%\n")

## 历史模拟法 ES: -3.91 %

cat("正态分布假设 ES:", round(es_gaussian * 100, 2), "%\n")

## 正态分布假设 ES: -3.32 %

risk_results[["performance"]] <- list(
  VaR_historical = VaR_historical,
  VaR_gaussian = VaR_gaussian,
  es_historical = es_historical,
  es_gaussian = es_gaussian
)

fGarch 包

# 使用fGarch包计算ES
library(fGarch)
# 拟合GARCH模型并预测
garch_model <- garchFit(~garch(1,1), data = returns, trace = FALSE)
forecast <- predict(garch_model, n.ahead = 1)
es_garch <- tail(fGarch::ES(garch_model), n = 1)
var_garch <- tail(fGarch::VaR(garch_model), n = 1)
cat("\nfGarch包计算结果:\n")

## 
## fGarch包计算结果:

cat("GARCH模型预测 VaR:", round(var_garch * 100, 2), "%\n")

## GARCH模型预测 VaR: 1.47 %

cat("GARCH模型预测 ES:", round(es_garch * 100, 2), "%\n")

## GARCH模型预测 ES: 1.86 %

risk_results[["fgarch"]] <- list(
  var_garch = var_garch,
  es_garch = es_garch
)

cvar包（针对连续分布的VaR/ES计算）

if (!require("cvar")) install.packages("cvar")

## Loading required package: cvar

## 
## Attaching package: 'cvar'

## The following objects are masked from 'package:PerformanceAnalytics':
## 
##     ES, VaR

library(cvar)

# 标准正态分布计算
gaussian_es <- list(
  qf = cvar::ES(qnorm, dist.type = "qf"),  # 使用分位数函数
  cdf = cvar::ES(pnorm, dist.type = "cdf")  # 使用累积分布函数
)

# t分布计算（df=4）
t_dist_es <- list(
  qf = cvar::ES(qt, dist.type = "qf", df = 4),
  cdf = cvar::ES(pt, dist.type = "cdf", df = 4)
)

risk_results[["cvar"]] <- list(
 gaussian_es = gaussian_es,
 t_dist_es = t_dist_es
)

5. rugarch 包

# 加载必要的包
library(rugarch)

## Loading required package: parallel

library(xts)  # 处理时间序列数据

# 定义GARCH模型（使用偏态t分布）
spec <- ugarchspec(
  variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
  mean.model = list(armaOrder = c(0, 0), include.mean = TRUE),
  distribution.model = "sstd"  # 偏态t分布
)

# 拟合模型
fit <- ugarchfit(spec, data = returns, solver = "hybrid")

# 使用固定参数重新过滤数据
spec2 <- spec
setfixed(spec2) <- as.list(coef(fit))
filt <- ugarchfilter(spec2, data = returns, n.old = 1000)

# 提取拟合值和波动率
mu <- fitted(filt)       # 条件均值
sigma_t <- sigma(filt)  # 条件波动率

# 获取分布参数
skew <- coef(fit)["skew"]
shape <- coef(fit)["shape"]

# 计算VaR (使用向量化方法)
p <- 0.05  # 置信水平(1-95%)
VaR <- mu + sigma_t * qdist("sstd", p = p, mu = 0, sigma = 1, skew = skew, shape = shape)

# 计算ES (使用数值积分)
es_integrand <- function(x) {
  qdist("sstd", p = x, mu = 0, sigma = 1, skew = skew, shape = shape)
}

# 向量化积分计算（对每个时间点单独计算ES）
ES <- numeric(length(returns))
for (i in 1:length(returns)) {
  int_result <- integrate(es_integrand, lower = 0, upper = p, subdivisions = 1000)
  ES[i] <- mu[i] + sigma_t[i] * int_result$value / p
}

# 确保所有向量长度一致（处理可能的长度差异）
min_len <- min(length(returns), length(VaR), length(ES))
actual <- as.numeric(tail(returns, min_len))
VaR <- tail(as.numeric(tail(VaR, min_len)),1)
ES <- tail(as.numeric(tail(ES, min_len)),1)

risk_results[["rugarch"]] <- list(
 VaR = VaR,
 ES = ES
)

4. gets包（基于ARX模型的条件风险度量）

# 加载必要的包
library(gets)

## 
## Attaching package: 'gets'

## The following objects are masked from 'package:cvar':
## 
##     ES, VaR

## The following objects are masked from 'package:fGarch':
## 
##     ES, VaR

## The following objects are masked from 'package:PerformanceAnalytics':
## 
##     ES, VaR

returns <- diff(log(Cl(IXIC))) 
returns <- as.vector(na.omit(returns))  # 强制转换为普通向量

# 设置置信水平
alpha <- 0.95

# 构建ARX模型（滞后1阶）
# 使用普通向量而非xts对象
arx_model <- arx(returns, mc=TRUE, ar=1)  # 直接使用向量格式

# 计算VaR和ES
risk_results[["gets"]] <- list(
  VaR = quantile(arx_model$residuals, probs=1-alpha),
  ES = mean(arx_model$residuals[arx_model$residuals <= quantile(arx_model$residuals, probs=1-alpha)])
)

# 打印结果
cat("ARX模型风险度量结果:\n")

## ARX模型风险度量结果:

cat("VaR (", alpha*100, "%):", round(risk_results[["gets"]]$VaR*100, 2), "%\n")

## VaR ( 95 %): -2.64 %

cat("ES (", alpha*100, "%):", round(risk_results[["gets"]]$ES*100, 2), "%\n")

## ES ( 95 %): -3.93 %

5. MSGARCH包（贝叶斯方法的风险预测）

if (!require("MSGARCH")) install.packages("MSGARCH")

## Loading required package: MSGARCH

library(MSGARCH)
spec_msgarch <- CreateSpec()
fit_msgarch <- FitML(spec_msgarch, data = returns)
risk_msgarch <- Risk(fit_msgarch, n.ahead = 1, level = alpha)
risk_results[["MSGARCH"]] <- list(
  VaR = as.numeric(risk_msgarch$VaR),
  ES = as.numeric(risk_msgarch$ES)
)

6. VaRES包（多种参数分布的风险度量）

if (!require("VaRES")) install.packages("VaRES")

## Loading required package: VaRES

library(VaRES)

returns <- diff(log(Cl(IXIC))) 
returns <- as.vector(na.omit(returns))  # 强制转换为普通向量

# 手动计算 VaR 和 ES（假设 df = 4）
df <- 4
alpha <- 0.05  # 95% 置信水平

# 计算 t 分布的 VaR
VaR_t <- -quantile(rt(n = 1e6, df = df), probs = alpha) * sd(returns)

# 计算 t 分布的 ES（解析解）
ES_t <- (dt(qt(alpha, df = df), df = df) / alpha) * 
        ((df + qt(alpha, df = df)^2) / (df - 2)) * 
        sd(returns)

risk_results[["VaRES"]] <- list(VaR = VaR_t, ES = ES_t)