diff --git a/quantitative_finance/risk_metrics.r b/quantitative_finance/risk_metrics.r new file mode 100644 index 00000000..d09b8fc0 --- /dev/null +++ b/quantitative_finance/risk_metrics.r @@ -0,0 +1,267 @@ +# Value at Risk (VaR) and Expected Shortfall (ES) Calculator +# Implements multiple VaR calculation methods and Expected Shortfall +# Features: Historical, Parametric, and Monte Carlo VaR/ES calculations + +library(R6) + +#' RiskMetrics Class +#' @description R6 class for calculating Value at Risk and Expected Shortfall +#' @details Implements multiple methods for VaR and ES calculation: +#' - Historical simulation (non-parametric) +#' - Parametric (variance-covariance) +#' - Monte Carlo simulation +#' Time complexity varies by method, documented per method +RiskMetrics <- R6Class( + "RiskMetrics", + + public = list( + #' @description Initialize risk calculator + #' @param returns Historical returns data + #' @param confidence_level Confidence level for VaR/ES (default 0.95) + #' @param time_horizon Time horizon in days (default 1) + initialize = function(returns = NULL, confidence_level = 0.95, time_horizon = 1) { + private$validate_parameters(confidence_level, time_horizon) + + self$returns <- returns + self$confidence_level <- confidence_level + self$time_horizon <- time_horizon + + if (!is.null(returns)) { + private$fit_distribution() + } + }, + + #' @description Calculate Historical VaR + #' @param portfolio_value Current portfolio value + #' @param method Calculation method ('historical', 'parametric', or 'monte_carlo') + #' @param n_simulations Number of simulations for Monte Carlo method + calculate_var = function(portfolio_value, method = "historical", n_simulations = 10000) { + if (is.null(self$returns)) { + stop("No returns data available. Please initialize with returns data.") + } + + method <- match.arg(method, c("historical", "parametric", "monte_carlo")) + + var_value <- switch(method, + "historical" = private$calculate_historical_var(portfolio_value), + "parametric" = private$calculate_parametric_var(portfolio_value), + "monte_carlo" = private$calculate_monte_carlo_var(portfolio_value, n_simulations) + ) + + # Scale VaR to time horizon + var_value * sqrt(self$time_horizon) + }, + + #' @description Calculate Expected Shortfall (Conditional VaR) + #' @param portfolio_value Current portfolio value + #' @param method Calculation method ('historical', 'parametric', or 'monte_carlo') + #' @param n_simulations Number of simulations for Monte Carlo method + calculate_es = function(portfolio_value, method = "historical", n_simulations = 10000) { + if (is.null(self$returns)) { + stop("No returns data available. Please initialize with returns data.") + } + + method <- match.arg(method, c("historical", "parametric", "monte_carlo")) + + es_value <- switch(method, + "historical" = private$calculate_historical_es(portfolio_value), + "parametric" = private$calculate_parametric_es(portfolio_value), + "monte_carlo" = private$calculate_monte_carlo_es(portfolio_value, n_simulations) + ) + + # Scale ES to time horizon + es_value * sqrt(self$time_horizon) + }, + + #' @description Generate risk report with multiple metrics + #' @param portfolio_value Current portfolio value + #' @param include_methods Which methods to include in report + generate_risk_report = function(portfolio_value, + include_methods = c("historical", "parametric", "monte_carlo")) { + results <- list() + + for (method in include_methods) { + results[[method]] <- list( + var = self$calculate_var(portfolio_value, method), + es = self$calculate_es(portfolio_value, method) + ) + } + + # Add distribution statistics + results$statistics <- list( + mean_return = mean(self$returns), + volatility = sd(self$returns), + skewness = private$calculate_skewness(), + kurtosis = private$calculate_kurtosis() + ) + + return(results) + }, + + #' @description Update returns data and recalculate distribution parameters + #' @param new_returns New returns data to use + update_returns = function(new_returns) { + self$returns <- new_returns + private$fit_distribution() + invisible(self) + }, + + # Public fields + returns = NULL, + confidence_level = NULL, + time_horizon = NULL + ), + + private = list( + # Distribution parameters + mean_return = NULL, + volatility = NULL, + + #' @description Fit distribution to returns data + fit_distribution = function() { + private$mean_return <- mean(self$returns) + private$volatility <- sd(self$returns) + }, + + #' @description Calculate Historical VaR + calculate_historical_var = function(portfolio_value) { + sorted_returns <- sort(self$returns) + index <- floor((1 - self$confidence_level) * length(sorted_returns)) + -sorted_returns[index] * portfolio_value + }, + + #' @description Calculate Parametric VaR + calculate_parametric_var = function(portfolio_value) { + z_score <- stats::qnorm(self$confidence_level) + portfolio_value * (z_score * private$volatility - private$mean_return) + }, + + #' @description Calculate Monte Carlo VaR + calculate_monte_carlo_var = function(portfolio_value, n_simulations) { + simulated_returns <- stats::rnorm( + n_simulations, + mean = private$mean_return, + sd = private$volatility + ) + sorted_returns <- sort(simulated_returns) + index <- floor((1 - self$confidence_level) * n_simulations) + -sorted_returns[index] * portfolio_value + }, + + #' @description Calculate Historical Expected Shortfall + calculate_historical_es = function(portfolio_value) { + sorted_returns <- sort(self$returns) + var_index <- floor((1 - self$confidence_level) * length(sorted_returns)) + tail_returns <- sorted_returns[1:var_index] + -mean(tail_returns) * portfolio_value + }, + + #' @description Calculate Parametric Expected Shortfall + calculate_parametric_es = function(portfolio_value) { + z_score <- stats::qnorm(self$confidence_level) + phi_z <- stats::dnorm(z_score) + lambda <- phi_z / (1 - self$confidence_level) + portfolio_value * (lambda * private$volatility - private$mean_return) + }, + + #' @description Calculate Monte Carlo Expected Shortfall + calculate_monte_carlo_es = function(portfolio_value, n_simulations) { + simulated_returns <- stats::rnorm( + n_simulations, + mean = private$mean_return, + sd = private$volatility + ) + sorted_returns <- sort(simulated_returns) + var_index <- floor((1 - self$confidence_level) * n_simulations) + tail_returns <- sorted_returns[1:var_index] + -mean(tail_returns) * portfolio_value + }, + + #' @description Calculate distribution skewness + calculate_skewness = function() { + r <- self$returns + n <- length(r) + m3 <- sum((r - mean(r))^3) / n + s3 <- sd(r)^3 + m3 / s3 + }, + + #' @description Calculate distribution kurtosis + calculate_kurtosis = function() { + r <- self$returns + n <- length(r) + m4 <- sum((r - mean(r))^4) / n + s4 <- sd(r)^4 + m4 / s4 - 3 # Excess kurtosis (normal = 0) + }, + + #' @description Validate input parameters + validate_parameters = function(confidence_level, time_horizon) { + if (confidence_level <= 0 || confidence_level >= 1) { + stop("Confidence level must be between 0 and 1") + } + if (time_horizon <= 0) { + stop("Time horizon must be positive") + } + } + ) +) + +# Demonstration +demonstrate_risk_metrics <- function() { + cat("=== Value at Risk and Expected Shortfall Demo ===\n\n") + + # Generate sample returns data + set.seed(42) + n_days <- 1000 + returns <- rnorm(n_days, mean = 0.0001, sd = 0.01) + + # Initialize calculator + risk_calc <- RiskMetrics$new( + returns = returns, + confidence_level = 0.95, + time_horizon = 1 + ) + + # Portfolio parameters + portfolio_value <- 1000000 # $1 million portfolio + + cat("Portfolio Parameters:\n") + cat(sprintf("Value: $%d\n", portfolio_value)) + cat(sprintf("Confidence Level: %.1f%%\n", risk_calc$confidence_level * 100)) + cat(sprintf("Time Horizon: %d day(s)\n\n", risk_calc$time_horizon)) + + # Calculate VaR using different methods + methods <- c("historical", "parametric", "monte_carlo") + + cat("Value at Risk (VaR) Results:\n") + for (method in methods) { + var_value <- risk_calc$calculate_var(portfolio_value, method) + cat(sprintf("%s VaR: $%.2f\n", tools::toTitleCase(method), var_value)) + } + cat("\n") + + cat("Expected Shortfall (ES) Results:\n") + for (method in methods) { + es_value <- risk_calc$calculate_es(portfolio_value, method) + cat(sprintf("%s ES: $%.2f\n", tools::toTitleCase(method), es_value)) + } + cat("\n") + + # Generate and display comprehensive risk report + cat("Comprehensive Risk Report:\n") + report <- risk_calc$generate_risk_report(portfolio_value) + + cat("\nDistribution Statistics:\n") + cat(sprintf("Mean Return: %.6f\n", report$statistics$mean_return)) + cat(sprintf("Volatility: %.6f\n", report$statistics$volatility)) + cat(sprintf("Skewness: %.6f\n", report$statistics$skewness)) + cat(sprintf("Excess Kurtosis: %.6f\n", report$statistics$kurtosis)) + + cat("\n=== Demo Complete ===\n") +} + +# Run demonstration if not in interactive mode +if (!interactive()) { + demonstrate_risk_metrics() +} \ No newline at end of file