1. Using the tools directly

There are some available package tools allow us to manage investment portfolio.

options(digits=4, scipen=100)
rm(list=ls(all=T))
Sys.setlocale('LC_ALL', 'C')
## [1] "C"
library(fPortfolio)
## Warning: package 'fPortfolio' was built under R version 4.0.5
## Warning: package 'timeDate' was built under R version 4.0.5
## Warning: package 'timeSeries' was built under R version 4.0.5
## Warning: package 'fBasics' was built under R version 4.0.5
## Warning: package 'fAssets' was built under R version 4.0.5
library(tseries)
## Warning: package 'tseries' was built under R version 4.0.5
library(quantmod)
library(modopt.matlab)
## Warning: package 'modopt.matlab' was built under R version 4.0.5
## Warning: package 'ROI.plugin.glpk' was built under R version 4.0.5

Select a group of stocks and set the observation period

DJ10 = c('MMM',  'AXP',  'AAPL', 'BA',  'CAT', 'CVX',
         'CSCO', 'KO',   'XOM', 'GE')
K = length(DJ10)
dd = as.Date(c('2015-10-01','2017-09-30'))

Download the historical stock prices and calculate the rate of return

mx =  sapply(DJ10, function(x) {
  G = getSymbols(x,from=dd[1],to=dd[2],auto.assign=F)
  as.numeric(Delt(G[,6])) })[-1,]
names(mx) = DJ10
ts = as.timeSeries(mx)

Set investment restrictions

cons = c('LongOnly')
spec = portfolioSpec()
setNFrontierPoints(spec) = 25
setSolver(spec) = "solveRquadprog"
frontier = portfolioFrontier(ts, spec, cons)

Draw the efficient frontier

par(cex=0.7)
tailoredFrontierPlot(frontier)

Calculate the best investment portfolio

par(cex=0.8)
weightsPlot(frontier,col=rainbow(length(DJ10)))




2. Manual optimization

Of course, we can also use R’s optimization (linear programming) suite to practice calculating the efficiency frontier and the best investment portfolio by ourselves.

Minimum Variance Portfolio, MVP

# Minimum Variance Portfolio, MVP
solution <- quadprog(
  H = cov(mx),      # objective function 
  f = rep(0, K),    # linear terms of objective
  A = NULL,         # Inequality Constrains (LHS)
  b = NULL,         # Inequality Constrains (RHS)
  Aeq = rep(1, K),  # Equality Constrains (LHS)  
  beq = 1,          # Equality Constrains (RHS)    
  lb = rep(0, K),   # lower bound
  ub = rep(1, K)    # upper bound  
  )
(portfolio = round(solution$x, 2))
##  [1] 0.20 0.14 0.08 0.00 0.00 0.00 0.00 0.46 0.10 0.03

Markowitz portfolio

Basically, the Markowitz portfolio is a minimally variable portfolio that limits average returns.

# define a optimization function
markowitz <- function(data, mu) {
  H = cov(data)
  f = rep(0, K)
  Aeq = rep(1, K)
  beq = 1
  A = -as.numeric(colMeans(data))
  b = -mu
  lb = rep(0, K)
  ub = rep(1, K)
  solution = quadprog(H, f, A, b, Aeq, beq, lb, ub)
  return(round(solution$x, 2))
  }

Practice calculating and draw the efficiency frontier of Markowitz

# calculate the frontier
P = 20    # number of points to draw
cm = colMeans(mx)
csd = apply(mx, 2, sd)
mu = c(); sd = c()
for(i in seq(mean(cm), max(cm)*0.99, length.out=P)) {
  portfolio = markowitz(mx, i)
  mu <- c(mu, mean(portfolio %*% t(mx)))
  sd <- c(sd,   sd(portfolio %*% t(mx))) }

# plot the frontier
plot(x = sd, y = mu, type='l', col='gold', lwd=3, 
     ylim=range(cm), xlim=c(min(sd),max(csd)),
     xlab="Standard Deviation", ylab="Return", 
     main="Efficient Frontier")
# points(x = sd, y = mu, col='pink', pch=19)
abline(h=seq(0,0.002,0.00025), v=seq(0, 0.02,0.001), col='lightgray', lty=3)
points(x=csd, y=cm, type='p', col='green', pch=19, cex=1.5)
text(x=csd, y=cm, labels=DJ10, col='darkgreen', pos=2, font=2)