I have computed the theoretical minimum variance portfolio using the 30 stocks in the Dow. The formula used is:
$$\underset{N\times 1}{\omega_{mvp}}=\frac{\lambda}{2}\cdot \Sigma^{-1}\iota=\frac{\Sigma^{-1}\iota}{\iota'\Sigma^{-1}\iota}$$
Where $\iota$ is a $N\times 1$ vector containing 1's.
For the data downloaded i get approximatly 0.0002712748, which i call Sigma_mvp in the script.
Then i generate one million different vectors containing 30 weights each for the assets. The weights can be negative, and i make sure they sum to one by dividing with the column sum.
My problem is: The smallest variance i manage to get with these random weights is 0.0004467729, so something must be wrong.
Any ideas? Hope the question is clear.
My code is provided below:
library(tidyquant)
library(tidyverse)
############ Getting Data for DOW ##############
tickers <- tq_index("DOW")
N <- length(tickers$symbol) # number of assets = 30
ones <- as.matrix(rep(1,N), ncol = 1) # one vector for later use
data <- tickers %>% tq_get(get = "stock.prices")
# calculate weekly returns
returns <- data %>%
group_by(symbol) %>%
tq_transmute(select = adjusted,
mutate_fun = to.weekly,
indexAt = "lastof") %>%
mutate(return = (log(adjusted) - log(lag(adjusted)))) %>%
na.omit()
# mean return vector
asset_returns <- returns %>% group_by(symbol) %>%
summarise(expected_return = mean(return)) %>%
select(expected_return) %>% as.matrix()
rownames(asset_returns) <- tickers$symbol %>% sort()
# create covariance matrix
Sigma <- returns %>%
select(-adjusted) %>%
pivot_wider(names_from = symbol, values_from = return) %>% # reorder data to a T x N matrix
na.omit() %>% # remove NA that got generated by "DOW"
select(-date) %>%
cov(use = "pairwise.complete.obs")
Sigma <- Sigma[rownames(asset_returns),rownames(asset_returns)] # reorder matrix to match asset_return vector sequence
############## Generating random portfolios ###################
# random weights
w_rdm <- matrix(runif(n = 1000000 * N, min = -3, max = 3), nrow = N)
w_rdm <- apply(w_rdm,2,function(x){x/sum(x)})
# Create points
eff_frontier_rdm <- matrix(0, nrow = 1000000, ncol = 2)
for(i in 1:ncol(w_rdm)){
eff_frontier_rdm[i,1] <- t(w_rdm[,i, drop = F]) %*% asset_returns
eff_frontier_rdm[i,2] <- t(w_rdm[,i, drop = F]) %*% Sigma %*% w_rdm[,i, drop = F]
}
colnames(eff_frontier_rdm) <- c("return", "variance")
eff_frontier_rdm <- eff_frontier_rdm %>% as_tibble()
# smallest variance achieved with random portoflios
min(eff_frontier_rdm$variance)
# Computing the minimum variance portfolio
lambda <- 2 / as.numeric((t(ones) %*% solve(Sigma) %*% ones))
w_mvp <- (solve(Sigma) %*% ones) * lambda/2
Sigma_mvp <- t(w_mvp) %*% Sigma %*% w_mvp
# theoretical min variance portoflio
Sigma_mvp
I have added a picture of the simulated portfolios with the code below, with the theoretically correct efficient frontier.
Based on the answers, i managed to create the following
Thanks for the different suggestions. Based on the above answers, I resolved my issue by the following procedure:
Given some covariance matrix Σ, and an expected return vector, I call asset_returns, I used the following steps:
Use the two mutual fund theorem to choose some random weight matrix, from a portfolio located on the frontier. Add some noise to each weight, by drawing from a standard normal. Normalize the created vector by dividing with the sum. Compute the portfolio variance and expected return with the created weight vector. After creating 1 million of such random weights, I manage to fill the area within the frontier. Creating more points would fill it all out.
The solution is based on the answers I got here on the post.
eff_frontier_rdm <- matrix(0, nrow = 1000000, ncol = 2)
for(i in 1:nrow(eff_frontier_rdm)){
c <- runif(1, min = -4, max = 4) # draw random number
w = c * w_mvp + (1-c) * port_2$w_eff # create weight
eps <- matrix(rnorm(N, mean = 0, sd = 0.1), ncol = 1)
w = w + eps
w = w / sum(w)
eff_frontier_rdm[i,1] <- t(w) %*% asset_returns
eff_frontier_rdm[i,2] <- t(w) %*% Sigma %*% w
}