Portfolio Optimization Comparison
Ethan Pinto
2023-11-28
Install necessary packages:
# install.packages("quantmod")
# install.packages("PerformanceAnalytics")
# install.packages("tseries")
# install.packages("R.matlab")
# install.packages("tidyquant")
# install.packages("tidyr")
# install.packages("dplyr")Load Libraries:
library(quantmod)## Loading required package: xts## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric## Loading required package: TTR## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(PerformanceAnalytics)##
## Attaching package: 'PerformanceAnalytics'## The following object is masked from 'package:graphics':
##
## legendlibrary(xts)
library(tseries)
library(quadprog)
library(ggplot2)
library(tidyquant)## Loading required package: lubridate##
## Attaching package: 'lubridate'## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, unionlibrary(tidyr)
library(dplyr)##
## ######################### 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: 'dplyr'## The following objects are masked from 'package:xts':
##
## first, last## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionThe first step is to pull in the historical stock data that we will use to create our portfolio and pass through the models:
stocks <- c("AAPL","MSFT", "GE", "LMT", "PFE", "NVO", "AMZN", "TSLA", "WFC", "BLK")
startDate <- as.Date("2018-01-01")
endDate <- as.Date("2022-10-31")
stockData <- list()
for (stock in stocks) {
stockData[[stock]] <- getSymbols(stock, src = "yahoo",
from = startDate,
to = endDate,
auto.assign = FALSE)
}
stockPrices <- do.call(merge, lapply(stockData@.Data, function(x) Ad(x)))Next, we will build the (equally weighted) portfolio:
prices <- do.call(merge, lapply(stockData, Ad))
returns <- Return.calculate(prices)
returns[is.na(returns)] <- 0
numStocks <- length(stockData)
weights <- rep(1/numStocks, numStocks)
portfolioReturns <- Return.portfolio(R = returns, weights = weights)
head(portfolioReturns)## portfolio.returns
## 2018-01-02 0.000000000
## 2018-01-03 0.006011433
## 2018-01-04 0.007244349
## 2018-01-05 0.007751192
## 2018-01-08 0.003792290
## 2018-01-09 0.003213404With this portfolio, we can now start testing our portfolio optimization methods.
The first method is Value-at-Risk and Conditional Value-at-Risk
confidenceLevel <- 0.95
VaR_Historical <- VaR(portfolioReturns, p = confidenceLevel, method = "historical")
CVaR_Historical <- ES(portfolioReturns, p = confidenceLevel, method = "historical")
cat("Historical VaR at", confidenceLevel*100, "% confidence level is:", VaR_Historical, "\n")## Historical VaR at 95 % confidence level is: -0.03115543cat("Historical CVaR at", confidenceLevel*100, "% confidence level is:", CVaR_Historical, "\n")## Historical CVaR at 95 % confidence level is: -0.04622044Next we will run a Monte Carlo simulation:
monteCarloOptimization <- function(returns, numPortfolios = 10000) {
numAssets <- ncol(returns)
results <- matrix(0, nrow = numPortfolios, ncol = 3 + numAssets) # Extend matrix to store weights
colnames(results) <- c("Return", "Risk", "Sharpe", paste("Weight", 1:numAssets, sep = "_"))
for (i in 1:numPortfolios) {
weights <- runif(numAssets)
weights <- weights / sum(weights)
portfolioReturn <- sum(apply(returns, 1, function(x) sum(x * weights)))
portfolioRisk <- sqrt(t(weights) %*% cov(returns) %*% weights)
sharpeRatio <- portfolioReturn / portfolioRisk
results[i,] <- c(portfolioReturn, portfolioRisk, sharpeRatio, weights)
}
return(results)
}
mcResults <- monteCarloOptimization(returns)
optimalPortfolio <- mcResults[which.max(mcResults[, "Sharpe"]),]
cat("Optimal Portfolio:\n")## Optimal Portfolio:cat(sprintf("Return: %f\n", optimalPortfolio["Return"]))## Return: 1.600277cat(sprintf("Risk: %f\n", optimalPortfolio["Risk"]))## Risk: 0.016983cat(sprintf("Sharpe Ratio: %f\n", optimalPortfolio["Sharpe"]))## Sharpe Ratio: 94.227205numAssets <- ncol(returns)
for (j in 1:numAssets) {
cat(sprintf("Weight of Asset %d: %f\n", j, optimalPortfolio[paste("Weight", j, sep = "_")]))
}## Weight of Asset 1: 0.089983
## Weight of Asset 2: 0.030653
## Weight of Asset 3: 0.003950
## Weight of Asset 4: 0.031873
## Weight of Asset 5: 0.206962
## Weight of Asset 6: 0.250239
## Weight of Asset 7: 0.085758
## Weight of Asset 8: 0.267519
## Weight of Asset 9: 0.020181
## Weight of Asset 10: 0.012883Next, we will do CAPM optimization. First, calculate the stock returns for further analysis:
stocks <- c("AAPL","MSFT", "GE", "LMT", "PFE", "NVO", "AMZN", "TSLA", "WFC", "BLK")
startDate <- as.Date("2018-01-01")
endDate <- as.Date("2022-10-31")
stockReturns <- list()
for (stock in stocks) {
tempReturns <- Delt(Ad(stockData[[stock]]))
tempReturns[1, ] <- 0 # Set the first day's return to 0
stockReturns[[stock]] <- tempReturns
}
mergedReturns <- do.call(merge, stockReturns)
head(mergedReturns)## Delt.1.arithmetic Delt.1.arithmetic.1 Delt.1.arithmetic.2
## 2018-01-02 0.000000000 0.0000000000 0.0000000000
## 2018-01-03 -0.000174141 0.0046538648 0.0094550831
## 2018-01-04 0.004645007 0.0088014679 0.0209365040
## 2018-01-05 0.011385554 0.0123979846 0.0005396223
## 2018-01-08 -0.003714225 0.0010206830 -0.0140236180
## 2018-01-09 -0.000114764 -0.0006794837 0.0153172270
## Delt.1.arithmetic.3 Delt.1.arithmetic.4 Delt.1.arithmetic.5
## 2018-01-02 0.000000000 0.000000000 0.0000000000
## 2018-01-03 0.008382108 0.007409472 0.0096013389
## 2018-01-04 0.017153680 0.002179333 -0.0034746850
## 2018-01-05 0.009059648 0.001902770 0.0045877822
## 2018-01-08 -0.004610340 -0.011123219 -0.0009132804
## 2018-01-09 0.007161003 -0.001097439 0.0031083693
## Delt.1.arithmetic.6 Delt.1.arithmetic.7 Delt.1.arithmetic.8
## 2018-01-02 0.000000000 0.000000000 0.000000000
## 2018-01-03 0.012775310 -0.010233113 0.007693562
## 2018-01-04 0.004476013 -0.008289978 0.012508093
## 2018-01-05 0.016162517 0.006229741 0.006738352
## 2018-01-08 0.014424679 0.062638220 -0.011314616
## 2018-01-09 0.004675710 -0.008085381 0.003545651
## Delt.1.arithmetic.9
## 2018-01-02 0.000000000
## 2018-01-03 0.010550846
## 2018-01-04 0.013201333
## 2018-01-05 0.008519984
## 2018-01-08 0.007458555
## 2018-01-09 0.008404273Next, we will need to pull in Market return data:
sp500Data <- getSymbols("^GSPC", src = "yahoo", from = startDate, to = endDate, auto.assign = FALSE)
if (is.null(sp500Data)) {
stop("Failed to retrieve S&P 500 data")
} else {
sp500Close <- Ad(sp500Data)
marketReturn <- Return.calculate(sp500Close)
marketReturn[1,] <- 0
head(marketReturn)
}## GSPC.Adjusted
## 2018-01-02 0.000000000
## 2018-01-03 0.006398819
## 2018-01-04 0.004028636
## 2018-01-05 0.007033767
## 2018-01-08 0.001662344
## 2018-01-09 0.001302932Next, we calculate the Beta using the CAPM method:
calculateBeta <- function(mergedReturns, marketReturn) {
combinedData <- cbind(coredata(mergedReturns), coredata(marketReturn))
covMat <- cov(combinedData)
beta <- covMat[1, 2] / var(coredata(marketReturn))
return(beta)
}
capmOptimization <- function(mergedReturns, marketReturn, riskFreeRate) {
numAssets <- ncol(mergedReturns)
betas <- numeric(numAssets)
for (i in 1:numAssets) {
betas[i] <- calculateBeta(mergedReturns[, i], marketReturn)
}
expectedReturns <- riskFreeRate + betas * (mean(marketReturn, na.rm = TRUE) - riskFreeRate)
covReturns <- cov(mergedReturns, use = "pairwise.complete.obs")
Dmat <- 2 * covReturns
dvec <- rep(0, numAssets)
Amat <- rbind(1, diag(numAssets))
bvec <- c(1, rep(0, numAssets))
portfolioOptim <- solve.QP(Dmat, dvec, t(Amat), bvec, meq = 1)
portfolioWeights <- portfolioOptim$solution
return(list(betas = betas, expectedReturns = expectedReturns, weights = portfolioWeights))
}
riskFreeRate <- 0.015
capmResults <- capmOptimization(mergedReturns, marketReturn, riskFreeRate)
format(capmResults, scientific = FALSE)## betas
## "1.2153067, 1.2030067, 1.1398328, 0.6943240, 0.6195502, 0.5809871, 1.0689141, 1.4330871, 1.1887334, 1.2212237"
## expectedReturns
## "-0.0027448360, -0.0025652424, -0.0016428319, 0.0048620942, 0.0059538748, 0.0065169399, -0.0006073399, -0.0059246723, -0.0023568364, -0.0028312301"
## weights
## "0.0000000000000000000000000, -0.0000000000000000001319786, 0.0174174226883696790790790, 0.2708489615816870865749877, 0.2682033891351394694346766, 0.2935596103868725603724954, 0.1279709420708866063343123, 0.0021416117236306591621597, 0.0198580624134140244774205, 0.0000000000000000026930492"print(capmResults)## $betas
## [1] 1.2153067 1.2030067 1.1398328 0.6943240 0.6195502 0.5809871 1.0689141
## [8] 1.4330871 1.1887334 1.2212237
##
## $expectedReturns
## [1] -0.0027448360 -0.0025652424 -0.0016428319 0.0048620942 0.0059538748
## [6] 0.0065169399 -0.0006073399 -0.0059246723 -0.0023568364 -0.0028312301
##
## $weights
## [1] 0.000000e+00 -1.319786e-19 1.741742e-02 2.708490e-01 2.682034e-01
## [6] 2.935596e-01 1.279709e-01 2.141612e-03 1.985806e-02 2.693049e-18Markowitz Optimization
markowitzOptimization <- function(mergedReturns, riskFreeRate) {
numAssets <- ncol(mergedReturns)
meanReturns <- colMeans(mergedReturns)
covReturns <- cov(mergedReturns)
Dmat <- 2 * covReturns
dvec <- rep(0, numAssets)
if (!is.matrix(Dmat) || nrow(Dmat) != ncol(Dmat) || nrow(Dmat) != length(dvec)) {
stop("Dmat or dvec dimension error")
}
Amat <- rbind(1, diag(numAssets))
bvec <- c(1, rep(0, numAssets))
if (ncol(Amat) != length(dvec) || nrow(Amat) != length(bvec)) {
stop("Amat or bvec dimension error")
}
sol <- solve.QP(Dmat, dvec, t(Amat), bvec, meq = 1)
weights <- sol$solution
return(weights)
}
riskFreeRate <- 0.015
optimalWeights <- markowitzOptimization(mergedReturns, riskFreeRate)
format(optimalWeights, scientific = FALSE)## [1] " 0.0000000000000000000000000" "-0.0000000000000000001319786"
## [3] " 0.0174174226883696790790790" " 0.2708489615816870865749877"
## [5] " 0.2682033891351394694346766" " 0.2935596103868725603724954"
## [7] " 0.1279709420708866063343123" " 0.0021416117236306591621597"
## [9] " 0.0198580624134140244774205" " 0.0000000000000000026930492"print(optimalWeights)## [1] 0.000000e+00 -1.319786e-19 1.741742e-02 2.708490e-01 2.682034e-01
## [6] 2.935596e-01 1.279709e-01 2.141612e-03 1.985806e-02 2.693049e-18(I also have a black litterman function i want to add here when I figure it out)
Mean Variance Frontier
meanVarianceFrontier <- function(mergedReturns, riskFreeRate, targetReturns) {
numAssets <- ncol(mergedReturns)
meanReturns <- colMeans(mergedReturns)
covReturns <- cov(mergedReturns)
Dmat <- 2 * covReturns
dvec <- rep(0, numAssets)
risks <- numeric(length(targetReturns))
weightsList <- vector("list", length(targetReturns))
for (i in seq_along(targetReturns)) {
Amat <- rbind(1, diag(numAssets), meanReturns)
bvec <- c(1, rep(0, numAssets), targetReturns[i])
sol <- solve.QP(Dmat, dvec, t(Amat), bvec, meq = 2)
risks[i] <- sqrt(t(sol$solution) %*% covReturns %*% sol$solution)
weightsList[[i]] <- sol$solution
}
efficientFrontier <- data.frame(Return = targetReturns, Risk = risks, Weights = weightsList)
return(efficientFrontier)
}
targetReturns <- seq(from = min(colMeans(mergedReturns)), to = max(colMeans(mergedReturns)), length.out = 10)
riskFreeRate <- 0.015
frontier <- meanVarianceFrontier(mergedReturns, riskFreeRate, targetReturns)
print(frontier)## Return Risk
## 1 -4.634187e-05 0.01225770
## 2 2.689091e-04 0.01225770
## 3 5.841600e-04 0.01225770
## 4 8.994110e-04 0.01280337
## 5 1.214662e-03 0.01523553
## 6 1.529913e-03 0.01916682
## 7 1.845164e-03 0.02391364
## 8 2.160415e-03 0.02922132
## 9 2.475666e-03 0.03499802
## 10 2.790917e-03 0.04107459
## Weights.c.0...1.31978639907573e.19..0.0174174226883697..0.270848961581687..
## 1 0.000000e+00
## 2 -1.319786e-19
## 3 1.741742e-02
## 4 2.708490e-01
## 5 2.682034e-01
## 6 2.935596e-01
## 7 1.279709e-01
## 8 2.141612e-03
## 9 1.985806e-02
## 10 2.693049e-18
## Weights.c.0...1.31978639907573e.19..0.0174174226883697..0.270848961581687...1
## 1 0.000000e+00
## 2 -1.319786e-19
## 3 1.741742e-02
## 4 2.708490e-01
## 5 2.682034e-01
## 6 2.935596e-01
## 7 1.279709e-01
## 8 2.141612e-03
## 9 1.985806e-02
## 10 2.693049e-18
## Weights.c.0...1.31978639907573e.19..0.0174174226883697..0.270848961581687...2
## 1 0.000000e+00
## 2 -1.319786e-19
## 3 1.741742e-02
## 4 2.708490e-01
## 5 2.682034e-01
## 6 2.935596e-01
## 7 1.279709e-01
## 8 2.141612e-03
## 9 1.985806e-02
## 10 2.693049e-18
## Weights.c..3.24544041989652e.18..0.0297718959886366...7.29493561413981e.18..
## 1 -3.245440e-18
## 2 2.977190e-02
## 3 -7.294936e-18
## 4 2.607006e-01
## 5 2.324713e-01
## 6 3.310897e-01
## 7 5.174021e-02
## 8 9.422625e-02
## 9 0.000000e+00
## 10 4.443706e-18
## Weights.c.4.69272179376424e.17..0.0785421256584486..4.80890527126556e.18..
## 1 4.692722e-17
## 2 7.854213e-02
## 3 4.808905e-18
## 4 1.950442e-01
## 5 1.371037e-01
## 6 3.662023e-01
## 7 -4.045212e-17
## 8 2.231077e-01
## 9 2.163452e-17
## 10 -2.486965e-17
## Weights.c.5.52105629765884e.17..0.0867261209731244..1.50096593021237e.18..
## 1 5.521056e-17
## 2 8.672612e-02
## 3 1.500966e-18
## 4 1.227461e-01
## 5 3.153111e-02
## 6 4.010511e-01
## 7 -3.047554e-17
## 8 3.579456e-01
## 9 -4.836212e-17
## 10 -1.631335e-17
## Weights.c.7.2138917269427e.17..0.0739225158645998...1.24835744175808e.18..
## 1 7.213892e-17
## 2 7.392252e-02
## 3 -1.248357e-18
## 4 1.489222e-02
## 5 0.000000e+00
## 6 4.086092e-01
## 7 -2.533174e-17
## 8 5.025761e-01
## 9 2.373321e-17
## 10 -1.197335e-17
## Weights.c.1.02077430195765e.16..0..7.26860060107578e.18..0...1.32502910183705e.17..
## 1 1.020774e-16
## 2 0.000000e+00
## 3 7.268601e-18
## 4 0.000000e+00
## 5 -1.325029e-17
## 6 3.279480e-01
## 7 -2.014579e-17
## 8 6.720520e-01
## 9 4.169394e-17
## 10 6.029592e-17
## Weights.c.1.25126038420057e.16..0..1.02951345012955e.17..0...6.07249133874413e.17..
## 1 1.251260e-16
## 2 0.000000e+00
## 3 1.029513e-17
## 4 0.000000e+00
## 5 -6.072491e-17
## 6 1.639740e-01
## 7 -2.696022e-17
## 8 8.360260e-01
## 9 5.400600e-17
## 10 1.451594e-16
## Weights.c.1.48174646644348e.16..0..1.33216684015151e.17..0...1.08199535756512e.16..
## 1 1.481746e-16
## 2 0.000000e+00
## 3 1.332167e-17
## 4 0.000000e+00
## 5 -1.081995e-16
## 6 4.163336e-16
## 7 -3.377464e-17
## 8 1.000000e+00
## 9 6.631805e-17
## 10 1.745117e-16Visualize the Mean Variance Frontier
ggplot(frontier, aes(x = Risk, y = Return)) +
geom_line(color = "blue") +
labs(title = "Mean-Variance Efficient Frontier",
x = "Portfolio Risk (Standard Deviation)",
y = "Expected Portfolio Return") +
theme_minimal()
Next, we pull in data from after our historical period to test each portfolio against: (Note, important to run this chunk again to overwrite the original (historical) data from the previous chunks)
stocks <- c("AAPL","MSFT", "GE", "LMT", "PFE", "NVO", "AMZN", "TSLA", "WFC", "BLK")
startDate <- as.Date("2023-01-01")
endDate <- as.Date("2023-10-31")
stockData <- list()
for (stock in stocks) {
stockData[[stock]] <- getSymbols(stock, src = "yahoo",
from = startDate,
to = endDate,
auto.assign = FALSE)
}
stockReturn <- list()
for (stock in stocks) {
tempReturns <- Delt(Ad(stockData[[stock]]))
tempReturns[1, ] <- 0
stockReturn[[stock]] <- tempReturns
}
stockReturns <- do.call(merge, stockReturn)
head(stockReturns)## Delt.1.arithmetic Delt.1.arithmetic.1 Delt.1.arithmetic.2
## 2023-01-03 0.000000000 0.000000000 0.000000000
## 2023-01-04 0.010314283 -0.043743227 0.058204185
## 2023-01-05 -0.010604730 -0.029637684 0.015527052
## 2023-01-06 0.036794126 0.011785399 0.009117706
## 2023-01-09 0.004088992 0.009736197 0.010147280
## 2023-01-10 0.004456326 0.007617200 0.035778142
## Delt.1.arithmetic.3 Delt.1.arithmetic.4 Delt.1.arithmetic.5
## 2023-01-03 0.000000000 0.000000000 0.000000000
## 2023-01-04 -0.002156975 -0.022044424 -0.001387312
## 2023-01-05 0.001196200 -0.009375697 -0.003802123
## 2023-01-06 -0.008028137 0.025372509 0.017541330
## 2023-01-09 -0.030111597 -0.049685684 -0.010170301
## 2023-01-10 0.007189739 -0.015912385 -0.026014578
## Delt.1.arithmetic.6 Delt.1.arithmetic.7 Delt.1.arithmetic.8
## 2023-01-03 0.000000000 0.000000000 0.0000000000
## 2023-01-04 -0.007923565 0.051248853 0.0205790589
## 2023-01-05 -0.023725589 -0.029039098 -0.0053927868
## 2023-01-06 0.035611152 0.024651090 0.0089579737
## 2023-01-09 0.014869874 0.059349011 -0.0095793148
## 2023-01-10 0.028731709 -0.007681374 -0.0007076527
## Delt.1.arithmetic.9
## 2023-01-03 0.000000000
## 2023-01-04 0.013847464
## 2023-01-05 -0.028618845
## 2023-01-06 0.052421341
## 2023-01-09 0.020311699
## 2023-01-10 0.005697328Equally Weighted Portfolio Returns (re-running this chunk using the new data)
prices <- do.call(merge, lapply(stockData, Ad))
returns <- Return.calculate(prices)
returns[is.na(returns)] <- 0
numStocks <- length(stockData)
weights <- rep(1/numStocks, numStocks)
portfolioReturns <- Return.portfolio(R = returns, weights = weights)
tail(portfolioReturns)## portfolio.returns
## 2023-10-23 0.0031114703
## 2023-10-24 0.0122825494
## 2023-10-25 -0.0098386475
## 2023-10-26 -0.0174554965
## 2023-10-27 0.0003512552
## 2023-10-30 0.0122524775Monte Carlo Optimization Returns
(get weights from function return)
mc_weights <- c(0.21859, 0.08281, 0.188986, 0.01652, 0.000796, 0.19554, 0.098302, 0.108662, 0.036544, 0.053249)
mc_portfolioReturns <- Return.portfolio(R = returns, weights = mc_weights)## Warning in Return.portfolio.geometric(R = R, weights = weights, wealth.index =
## wealth.index, : The weights for one or more periods do not sum up to 1:
## assuming a return of 0 for the residual weightstail(mc_portfolioReturns)## portfolio.returns
## 2023-10-23 0.0043549297
## 2023-10-24 0.0193693869
## 2023-10-25 -0.0136348172
## 2023-10-26 -0.0228415861
## 2023-10-27 0.0002376274
## 2023-10-30 0.0150845036CAPM Optimization Returns
capm_weights <- c(0.0000000000000000000000000, -0.0000000000000000001094735, 0.0174175751737773581995228, 0.2708490891210338369354815, 0.2682027762553574334880579, 0.2935597422412994927931607, 0.1279711440899859720143894, 0.0021415150238913334412860, 0.0198581580946546529253816, -0.0000000000000000160630695)
capm_portfolioReturns <- Return.portfolio(R = returns, weights = capm_weights)
tail(capm_portfolioReturns)## portfolio.returns
## 2023-10-23 0.007082206
## 2023-10-24 0.001191459
## 2023-10-25 -0.005836649
## 2023-10-26 -0.010727712
## 2023-10-27 -0.003255237
## 2023-10-30 0.023326329Markowitz Optimization
marko_weights <- c(0.0000000000000000000000000, -0.0000000000000000001094735, 0.0174175751737773581995228, 0.2708490891210338369354815, 0.2682027762553574334880579, 0.2935597422412994927931607,
0.1279711440899859720143894, 0.0021415150238913334412860, 0.0198581580946546529253816, -0.0000000000000000160630695)
marko_portfolioReturns <- Return.portfolio(R = returns, weights = marko_weights)
tail(marko_portfolioReturns)## portfolio.returns
## 2023-10-23 0.007082206
## 2023-10-24 0.001191459
## 2023-10-25 -0.005836649
## 2023-10-26 -0.010727712
## 2023-10-27 -0.003255237
## 2023-10-30 0.023326329Plotting the performance of each portfolio:
portfolioReturns_df <- as.data.frame(portfolioReturns)
portfolioReturns_tibble <- as_tibble(portfolioReturns_df, rownames = "Date")
mc_portfolioReturns_df <- as.data.frame(mc_portfolioReturns)
mc_portfolioReturns_tibble <- as_tibble(mc_portfolioReturns_df, rownames = "Date")
capm_portfolioReturns_df <- as.data.frame(capm_portfolioReturns)
capm_portfolioReturns_tibble <- as_tibble(capm_portfolioReturns_df, rownames = "Date")
marko_portfolioReturns_df <- as.data.frame(marko_portfolioReturns)
marko_portfolioReturns_tibble <- as_tibble(marko_portfolioReturns_df, rownames = "Date")
portfolioReturns_tibble <- portfolioReturns_tibble %>% mutate(Portfolio = "Portfolio 1") %>% rename(Returns = portfolio.returns)
mc_portfolioReturns_tibble <- mc_portfolioReturns_tibble %>% mutate(Portfolio = "MC Portfolio") %>% rename(Returns = portfolio.returns)
capm_portfolioReturns_tibble <- capm_portfolioReturns_tibble %>% mutate(Portfolio = "CAPM Portfolio") %>% rename(Returns = portfolio.returns)
marko_portfolioReturns_tibble <- marko_portfolioReturns_tibble %>% mutate(Portfolio = "Markowitz Portfolio") %>% rename(Returns = portfolio.returns)
combinedDf <- bind_rows(
portfolioReturns_tibble,
mc_portfolioReturns_tibble,
capm_portfolioReturns_tibble,
marko_portfolioReturns_tibble
)
# Define a color palette
color_palette <- c("Portfolio 1" = "blue",
"MC Portfolio" = "green",
"CAPM Portfolio" = "red",
"Markowitz Portfolio" = "purple")
ggplot(combinedDf, aes(x = as.Date(Date), y = Returns, color = Portfolio)) +
geom_line() +
labs(title = "Comparison of Portfolio Performance", x = "Date", y = "Returns") +
theme_minimal()
Calculate and graph the variance of each portfolio:
variance1 <- var(portfolioReturns, na.rm = TRUE)
variance2 <- var(mc_portfolioReturns, na.rm = TRUE)
variance3 <- var(capm_portfolioReturns, na.rm = TRUE)
variance4 <- var(marko_portfolioReturns, na.rm = TRUE)
variance_data <- data.frame(
Portfolio = c("Equally Weighted", "Monte Carlo", "CAPM", "Markowitz"),
Variance = c(variance1, variance2, variance3, variance4)
)
print(paste("Variance for Portfolio 1:", variance1))## [1] "Variance for Portfolio 1: 0.000128153727390074"print(paste("Variance for Portfolio 2:", variance2))## [1] "Variance for Portfolio 2: 0.000144182771129367"print(paste("Variance for Portfolio 3:", variance3))## [1] "Variance for Portfolio 3: 9.25578927399404e-05"print(paste("Variance for Portfolio 4:", variance4))## [1] "Variance for Portfolio 4: 9.25578927399404e-05"variance_data$Variance <- variance_data$Variance * 100
ggplot(variance_data, aes(x = Portfolio, y = Variance, fill = Portfolio)) +
geom_bar(stat = "identity") +
scale_y_continuous(labels = scales::percent_format(scale = 1)) +
labs(title = "Variance of Portfolios as Percentage",
x = "Portfolio",
y = "Variance (%)") +
theme_minimal() +
theme(legend.position = "none")
Calculate Beta of each portfolio for comparison (regression approach)
First, pull in the new market data:
sp500Data <- getSymbols("^GSPC", src = "yahoo", from = startDate, to = endDate, auto.assign = FALSE)
if (is.null(sp500Data)) {
stop("Failed to retrieve S&P 500 data")
} else {
sp500Close <- Ad(sp500Data)
marketReturn <- Return.calculate(sp500Close)
marketReturn[1,] <- 0
head(marketReturn)
}## GSPC.Adjusted
## 2023-01-03 0.0000000000
## 2023-01-04 0.0075389706
## 2023-01-05 -0.0116455289
## 2023-01-06 0.0228407810
## 2023-01-09 -0.0007676325
## 2023-01-10 0.0069782332calculate_beta1 <- function(portfolioReturns, marketReturn) {
model <- lm(portfolioReturns ~ marketReturn)
return(coef(model)[2])
}
calculate_beta2 <- function(mc_portfolioReturns, marketReturn) {
model <- lm(mc_portfolioReturns ~ marketReturn)
return(coef(model)[2])
}
calculate_beta3 <- function(capm_portfolioReturns, marketReturn) {
model <- lm(capm_portfolioReturns ~ marketReturn)
return(coef(model)[2])
}
calculate_beta4 <- function(marko_portfolioReturns, marketReturn) {
model <- lm(marko_portfolioReturns ~ marketReturn)
return(coef(model)[2])
}
beta1 <- calculate_beta1(portfolioReturns, marketReturn)
beta2 <- calculate_beta2(mc_portfolioReturns, marketReturn)
beta3 <- calculate_beta3(capm_portfolioReturns, marketReturn)
beta4 <- calculate_beta4(marko_portfolioReturns, marketReturn)
print(paste("Beta for Portfolio 1:", beta1))## [1] "Beta for Portfolio 1: 1.17609654011679"print(paste("Beta for Portfolio 2:", beta2))## [1] "Beta for Portfolio 2: 1.20002537294762"print(paste("Beta for Portfolio 3:", beta3))## [1] "Beta for Portfolio 3: 0.659289276710059"print(paste("Beta for Portfolio 4:", beta4))## [1] "Beta for Portfolio 4: 0.659289276710059"Calculate each portfolio’s Alpha value:
(commented out until I figure this piece out) calculate_alpha1 <- function(portfolioReturns, beta1, marketReturn, riskFreeRate) { averagePortfolioReturn <- mean(portfolioReturns) averageMarketReturn <- mean(marketReturn)
expectedReturn <- riskFreeRate + beta1 * (averageMarketReturn - riskFreeRate)
alpha <- averagePortfolioReturn - expectedReturn
return(alpha)}
calculate_alpha2 <- function(mc_portfolioReturns, beta2, marketReturn, riskFreeRate) { expectedReturn <- riskFreeRate + beta2 * (marketReturn - riskFreeRate) alpha <- mean(mc_portfolioReturns) - expectedReturn return(alpha) }
calculate_alpha3 <- function(capm_portfolioReturns, beta3, marketReturn, riskFreeRate) { expectedReturn <- riskFreeRate + beta3 * (marketReturn - riskFreeRate) alpha <- mean(capm_portfolioReturns) - expectedReturn return(alpha) }
calculate_alpha4 <- function(marko_portfolioReturns, beta4, marketReturn, riskFreeRate) { expectedReturn <- riskFreeRate + beta4 * (marketReturn - riskFreeRate) alpha <- mean(marko_portfolioReturns) - expectedReturn return(alpha) }
alpha1 <- calculate_alpha1(portfolioReturns, beta1, marketReturn, riskFreeRate) alpha2 <- calculate_alpha2(mc_portfolioReturns, beta2, marketReturn, riskFreeRate) alpha3 <- calculate_alpha3(capm_portfolioReturns, beta3, marketReturn, riskFreeRate) alpha4 <- calculate_alpha4(marko_portfolioReturns, beta4, marketReturn, riskFreeRate)
print(paste(“Alpha for Portfolio 1:”, alpha1)) print(paste(“Alpha for MC Portfolio:”, alpha2)) print(paste(“Alpha for CAPM Portfolio:”, alpha3)) print(paste(“Alpha for Markowitz Portfolio:”, alpha4))
Visualizing Alpha:
alpha_data <- data.frame( Portfolio = c(“Equally Weighted”, “Monte Carlo”, “CAPM”, “Markowitz”), Alpha = c(alpha1, alpha2, alpha3, alpha4) )
Print column names
print(names(alpha_data)) print(alpha1)
Print the entire data frame
print(alpha_data)
Basic ggplot call
ggplot(alpha_data, aes(x = Portfolio, y = Alpha)) + geom_bar(stat = “identity”)
Calculating the total return for each portfolio:
```r
investment <- 100000
totalReturn1 <- cumprod(1 + portfolioReturns) - 1
totalReturn2 <- cumprod(1 + mc_portfolioReturns) - 1
totalReturn3 <- cumprod(1 + capm_portfolioReturns) - 1
totalReturn4 <- cumprod(1 + marko_portfolioReturns) - 1
finalReturn1 <- investment * (1 + as.numeric(last(totalReturn1)))
finalReturn2 <- investment * (1 + as.numeric(last(totalReturn2)))
finalReturn3 <- investment * (1 + as.numeric(last(totalReturn3)))
finalReturn4 <- investment * (1 + as.numeric(last(totalReturn4)))
returns_data <- data.frame(
Portfolio = c("Equally Weighted", "Monte Carlo", "CAPM", "Markowitz"),
FinalValue = c(finalReturn1, finalReturn2, finalReturn3, finalReturn4)
)
print(finalReturn1)## [1] 126632.8print(finalReturn2)## [1] 146053.1print(finalReturn3)## [1] 109660.1print(finalReturn4)## [1] 109660.1ggplot(returns_data, aes(x = Portfolio, y = FinalValue, fill = Portfolio)) +
geom_bar(stat = "identity") +
labs(title = "Final Monetary Value of Portfolios",
x = "Portfolio",
y = "Final Value (USD)") +
theme_minimal() +
theme(legend.position = "none") +
scale_y_continuous(labels = scales::dollar_format())