Downloading data

First lets load our packages

Next lets select a few stocks to build our portfolios.

We will choose the following 4 stocks 430090.BJ(同辉信息、设为2) 430300.BJ(辰光医疗、设为5) 430476.BJ(海能技术、设为8) 430418.BJ(苏轴股份、设为6)

Read the price data.

##选取投资组合对象
#430090.BJ(同辉信息、设为2) 430300.BJ(辰光医疗、设为5) 430476.BJ(海能技术、设为8) 430418.BJ(苏轴股份、设为6) 
##读取数据
price <-read.csv("C:/Users/86138/Desktop/投资组合管理/data.csv")
#修改数据格式
price[,1] = as.POSIXct(price[,1],format="%Y/%m/%d")
price = xts(price[,-1], order.by=price[,1])
##绘制走势图
names(price)<-c("2","5","8","6")
dygraph(price)

Next we will calculate the daily returns for these stocks. We will use the logarithmic returns.

pr_daily <- to.daily(price, indexAt = "lastof", OHLC = FALSE) 
#计算回报
rt_daily <- na.omit(Return.calculate(pr_daily, method = "log"))
dygraph(rt_daily)

Next let’s calculate the mean daily returns for each asset. Calculate the covariance matrix for all these stocks. We will annualize it by multiplying by 252.

mean_ret <- colMeans(rt_daily)
print(round(mean_ret, 4))
##       2       5       8       6 
## -0.0043 -0.0013 -0.0013  0.0010
cov_mat <- cov(rt_daily) * 252
print(round(cov_mat,4))
##        2      5      8      6
## 2 0.3079 0.0272 0.0473 0.0270
## 5 0.0272 0.0950 0.0141 0.0177
## 8 0.0473 0.0141 0.0920 0.0131
## 6 0.0270 0.0177 0.0131 0.1250

To calculate the portfolio returns and risk (standard deviation) we will us need

# Calculate the random weights
tick <- c("2","5","8","6")
wts <- runif(n = length(tick))
wts <- wts/sum(wts)
print(wts)
## [1] 0.29028400 0.07830176 0.36596221 0.26545203
# Calculate the portfolio returns
port_returns <- (sum(wts * mean_ret) + 1)^252 - 1
port_risk <- sqrt(t(wts) %*% (cov_mat %*% wts))
sharpe_ratio <- port_returns/port_risk

We have everything we need to perform our optimization. All we need now is to run this code on 10000 random portfolios. For that we will use a for loop.

Before we do that, we need to create empty vectors and matrix for storing our values.

# Creating a matrix to store the weights
num_port <- 10000

# Creating a matrix to store the weights

all_wts <- matrix(nrow = num_port,
                  ncol = length(tick))

# Creating an empty vector to store
# Portfolio returns

port_returns <- vector('numeric', length = num_port)

# Creating an empty vector to store
# Portfolio Standard deviation

port_risk <- vector('numeric', length = num_port)

# Creating an empty vector to store
# Portfolio Sharpe Ratio

sharpe_ratio <- vector('numeric', length = num_port)

Next lets run the for loop 5000 times.

for (i in seq_along(port_returns)) {
  
  wts <- runif(length(tick))
  wts <- wts/sum(wts)
  
  # Storing weight in the matrix
  all_wts[i,] <- wts
  
  # Portfolio returns
  
  port_ret <- sum(wts * mean_ret)
  port_ret <- ((port_ret + 1)^252) - 1
  
  # Storing Portfolio Returns values
  port_returns[i] <- port_ret
  
  
  # Creating and storing portfolio risk
  port_sd <- sqrt(t(wts) %*% (cov_mat  %*% wts))
  port_risk[i] <- port_sd
  
  # Creating and storing Portfolio Sharpe Ratios
  # Assuming 0% Risk free rate
  
  sr <- port_ret/port_sd
  sharpe_ratio[i] <- sr
}

All the heavy lifting has been done and now we can create a data table to store all the values together.

# Storing the values in the table
portfolio_values <- tibble(Return = port_returns,
                           Risk = port_risk,
                           SharpeRatio = sharpe_ratio)


library(tibble)
# Converting matrix to a tibble and changing column names
all_wts <- as_tibble(all_wts)
## Warning: The `x` argument of `as_tibble.matrix()` must have unique column names if
## `.name_repair` is omitted as of tibble 2.0.0.
## ℹ Using compatibility `.name_repair`.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
colnames(all_wts) <- colnames(pr_daily)
# Combing all the values together
portfolio_values <-  as_tibble(cbind(all_wts, portfolio_values))

We have the weights in each asset with the risk and returns along with the Sharpe ratio of each portfolio.

Next lets look at the portfolios that matter the most.

min_var <- portfolio_values[which.min(portfolio_values$Risk),]
max_sr <- portfolio_values[which.max(portfolio_values$SharpeRatio),]

Plot the weights of each portfolio. First with the minimum variance portfolio.

p <- min_var %>%
  gather(1:4, key = Asset,
         value = Weights) %>%
  mutate(Asset = as.factor(Asset)) %>%
  ggplot(aes(x = fct_reorder(Asset,Weights), y = Weights, fill = Asset)) +
  geom_bar(stat = 'identity') +
  theme_minimal() +
  labs(x = 'Assets', y = 'Weights', title = "Minimum Variance Portfolio Weights") +
  scale_y_continuous(labels = scales::percent) 


ggplotly(p)

Next lets look at the tangency portfolio or the the portfolio with the highest sharpe ratio.

p <- max_sr %>%
  gather(1:4, key = Asset,
         value = Weights) %>%
  mutate(Asset = as.factor(Asset)) %>%
  ggplot(aes(x = fct_reorder(Asset,Weights), y = Weights, fill = Asset)) +
  geom_bar(stat = 'identity') +
  theme_minimal() +
  labs(x = 'Assets', y = 'Weights', title = "Tangency Portfolio Weights") +
  scale_y_continuous(labels = scales::percent) 


ggplotly(p)

Finally lets plot all the random portfolios and visualize the efficient frontier.

p <- portfolio_values %>%
  ggplot(aes(x = Risk, y = Return, color = SharpeRatio)) +
  geom_point() +
  theme_classic() +
  scale_y_continuous(labels = scales::percent) +
  scale_x_continuous(labels = scales::percent) +
  labs(x = 'Annualized Risk',
       y = 'Annualized Returns',
       title = "Portfolio Optimization & Efficient Frontier") +
  geom_point(aes(x = Risk,
                 y = Return), data = min_var, color = 'red') +
  geom_point(aes(x = Risk,
                 y = Return), data = max_sr, color = 'red') #+
  #annotate('text', x = 0.20, y = 0.42, label = "Tangency Portfolio") +
  #annotate('text', x = 0.18, y = 0.01, label = "Minimum variance portfolio") +
  #annotate(geom = 'segment', x = 0.14, xend = 0.135,  y = 0.01, 
  #         yend = 0.06, color = 'red', arrow = arrow(type = "open")) +
  #annotate(geom = 'segment', x = 0.22, xend = 0.2275,  y = 0.405, 
  #         yend = 0.365, color = 'red', arrow = arrow(type = "open"))
  

ggplotly(p)