# Install necessary packages if not already installed
packages <- c("quadprog", "tidyverse", "PerformanceAnalytics", "xts")
new_packages <- packages[!(packages %in% installed.packages()[,"Package"])]
if(length(new_packages)) install.packages(new_packages)
# Load the libraries
library(quadprog)
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.0.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(PerformanceAnalytics)
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
##
## ######################### Warning from 'xts' package ##########################
## # #
## # The dplyr lag() function breaks how base R's lag() function is supposed to #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or #
## # source() into this session won't work correctly. #
## # #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop #
## # dplyr from breaking base R's lag() function. #
## # #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning. #
## # #
## ###############################################################################
##
## Attaching package: 'xts'
##
## The following objects are masked from 'package:dplyr':
##
## first, last
##
##
## Attaching package: 'PerformanceAnalytics'
##
## The following object is masked from 'package:graphics':
##
## legend
library(xts)
set.seed(123) # For reproducibility
n_assets <- 4 # Number of assets
n_obs <- 100 # Number of observations
# Simulate asset returns with some correlation
returns <- matrix(rnorm(n_obs * n_assets, mean = 0.001, sd = 0.02),
nrow = n_obs, ncol = n_assets)
colnames(returns) <- c("Asset1", "Asset2", "Asset3", "Asset4")
# Convert to time series
dates <- Sys.Date() - seq(n_obs, 1)
returns_xts <- xts(returns, order.by = dates)
# Display first few rows
head(returns_xts)
## Asset1 Asset2 Asset3 Asset4
## 2024-12-22 -0.010209513 -0.0132081313 0.044976207 -0.013304844
## 2024-12-23 -0.003603550 0.0061376742 0.027248260 -0.014053779
## 2024-12-24 0.032174166 -0.0039338376 -0.004302901 -0.017770774
## 2024-12-25 0.002410168 -0.0059508520 0.011863881 -0.020050266
## 2024-12-26 0.003585755 -0.0180323713 -0.007286799 -0.007743191
## 2024-12-27 0.035301300 0.0000994455 -0.008524938 0.007623583
cov_matrix <- cov(returns_xts)
print(cov_matrix)
## Asset1 Asset2 Asset3 Asset4
## Asset1 3.332931e-04 -1.748843e-05 -4.480158e-05 -1.671857e-05
## Asset2 -1.748843e-05 3.740252e-04 1.123498e-05 1.761177e-05
## Asset3 -4.480158e-05 1.123498e-05 3.609080e-04 -1.770789e-05
## Asset4 -1.671857e-05 1.761177e-05 -1.770789e-05 4.316266e-04
risk_parity <- function(cov_matrix) {
n <- ncol(cov_matrix)
# Define the quadratic optimization problem
Dmat <- cov_matrix
dvec <- rep(0, n)
Amat <- cbind(rep(1, n), diag(n)) # Constraints: weights sum to 1, non-negative weights
bvec <- c(1, rep(0, n)) # Sum of weights = 1, and no short selling
# Solve quadratic programming problem
result <- solve.QP(Dmat, dvec, Amat, bvec, meq = 1)
return(result$solution) # Return optimal weights
}
rp_weights <- risk_parity(cov_matrix)
names(rp_weights) <- colnames(returns_xts)
print(rp_weights)
## Asset1 Asset2 Asset3 Asset4
## 0.3062252 0.2176182 0.2704764 0.2056803
# Compute asset volatilities and portfolio volatility
vols <- sqrt(diag(cov_matrix))
port_vol <- sqrt(t(rp_weights) %*% cov_matrix %*% rp_weights)
# Compute risk contributions
risk_contributions <- (rp_weights * vols) / port_vol
## Warning in (rp_weights * vols)/port_vol: Recycling array of length 1 in vector-array arithmetic is deprecated.
## Use c() or as.vector() instead.
print(risk_contributions) # Should be nearly equal for all assets
## Asset1 Asset2 Asset3 Asset4
## 0.6147521 0.4627982 0.5650327 0.4698862