Portfolio Optimization in R: error

Question

I don't have much experience in R but I am trying to run a constrained portfolio optimization with R based on fPortfolio package. Using the function portfolioFrontier I get an error message:

Error in `colnames<-`(`*tmp*`, value = names(getMu(Data))) : 
attempt to set 'colnames' on an object with less than two dimensions

My code:

#specs
spec <- portfolioSpec()
setNFrontierPoints(spec)<-200
setRiskFreeRate(spec)<-0.001366179
constraints<-c("Long only", "maxW[1]=.05","maxW[2]=.25")
setTargetReturn(Spec) = mean(colMeans(as.timeSeries(portfolioXTS)))
#analysis
portfolioConstraints(as.timeSeries(portfolioXTS), spec, constraints)
portfolioData(as.timeSeries(portfolioXTS), spec = portfolioSpec())
efficientPortfolio(as.timeSeries(portfolioXTS), spec, constraints)

frontier <- portfolioFrontier(as.timeSeries(portfolioXTS),spec, 
constraints,include.mvl = TRUE, title="Traditional Asset Mix")

tangencyPortfolio(as.timeSeries(portfolioXTS), spec, constraints)
minvariancePortfolio(as.timeSeries(portfolioXTS), spec, constraints)
print(frontier)

I found on the forum that someone had a similar problem but there was no solution. I know that the number of rows must be greater than the number of columns but this requirement is satisfied.

I'd appreciate your help!

Data Sample:

dput

> dput(head(portfolioXTS))
structure(c(0.00156666666666667, 0.00158333333333333, 0.00169166666666667, 
0.0016, 0.00158025, 0.00155, 0.00673371680501522, 0.00900662025078258, 
-0.0240575754236355, 0.0248764429466077, 0.00560202723890391, 
0.0141166753754614, -0.0301780036878394, -0.00848891285361475, 
0.0444382350834758, -0.0335723121416733, 0.00238564176285805, 
-0.0603502320765891, 0.032911309501205, 0.0157910651749094, 0.0593749401684385, 
0.0061188298291901, -0.0166768949719553, -0.0753692400590273), .indexTZ = "UTC", .indexCLASS = "Date", tclass = "Date", tzone = c(TZ = "UTC"), class = c("xts", 
"zoo"), index = structure(c(1009843200, 1012521600, 1014940800, 
1017532800, 1020124800, 1022803200), tzone = "UTC", tclass = "Date"), .Dim = c(6L, 
4L), .Dimnames = list(NULL, c("Cash", "Bonds", "Equities_dev", 
"Equities_em")))

Answer 1

How about this option?

# Economist at Large
# Modern Portfolio Theory
# Use solve.QP to solve for efficient frontier
# Last Edited 5/3/13

# This file uses the solve.QP function in the quadprog package to solve for the
# efficient frontier.
# Since the efficient frontier is a parabolic function, we can find the solution
# that minimizes portfolio variance and then vary the risk premium to find
# points along the efficient frontier. Then simply find the portfolio with the
# largest Sharpe ratio (expected return / sd) to identify the most
# efficient portfolio

library(stockPortfolio) # Base package for retrieving returns
library(ggplot2) # Used to graph efficient frontier

## Warning: package 'ggplot2' was built under R version 3.3.3

library(reshape2) # Used to melt the data

## Warning: package 'reshape2' was built under R version 3.3.3

library(quadprog) #Needed for solve.QP

# Create the portfolio using ETFs, incl. hypothetical non-efficient allocation
stocks <- c(
  "VTSMX" = .0,
  "SPY" = .20,
  "EFA" = .10,
  "IWM" = .10,
  "VWO" = .30,
  "LQD" = .20,
  "HYG" = .10)

# Retrieve returns, from earliest start date possible (where all stocks have
# data) through most recent date
returns <- getReturns(names(stocks[-1]), freq="week") #Currently, drop index

#### Efficient Frontier function ####
eff.frontier <- function (returns, short="no", max.allocation=NULL,
                          risk.premium.up=.5, risk.increment=.005){
  # return argument should be a m x n matrix with one column per security
  # short argument is whether short-selling is allowed; default is no (short
  # selling prohibited)max.allocation is the maximum % allowed for any one
  # security (reduces concentration) risk.premium.up is the upper limit of the
  # risk premium modeled (see for loop below) and risk.increment is the
  # increment (by) value used in the for loop

  covariance <- cov(returns)
  print(covariance)
  n <- ncol(covariance)

  # Create initial Amat and bvec assuming only equality constraint
  # (short-selling is allowed, no allocation constraints)
  Amat <- matrix (1, nrow=n)
  bvec <- 1
  meq <- 1

  # Then modify the Amat and bvec if short-selling is prohibited
  if(short=="no"){
    Amat <- cbind(1, diag(n))
    bvec <- c(bvec, rep(0, n))
  }

  # And modify Amat and bvec if a max allocation (concentration) is specified
  if(!is.null(max.allocation)){
    if(max.allocation > 1 | max.allocation <0){
      stop("max.allocation must be greater than 0 and less than 1")
    }
    if(max.allocation * n < 1){
      stop("Need to set max.allocation higher; not enough assets to add to 1")
    }
    Amat <- cbind(Amat, -diag(n))
    bvec <- c(bvec, rep(-max.allocation, n))
  }

  # Calculate the number of loops
  loops <- risk.premium.up / risk.increment + 1
  loop <- 1

  # Initialize a matrix to contain allocation and statistics
  # This is not necessary, but speeds up processing and uses less memory
  eff <- matrix(nrow=loops, ncol=n+3)
  # Now I need to give the matrix column names
  colnames(eff) <- c(colnames(returns), "Std.Dev", "Exp.Return", "sharpe")

  # Loop through the quadratic program solver
  for (i in seq(from=0, to=risk.premium.up, by=risk.increment)){
    dvec <- colMeans(returns) * i # This moves the solution along the EF
    sol <- solve.QP(covariance, dvec=dvec, Amat=Amat, bvec=bvec, meq=meq)
    eff[loop,"Std.Dev"] <- sqrt(sum(sol$solution*colSums((covariance*sol$solution))))
    eff[loop,"Exp.Return"] <- as.numeric(sol$solution %*% colMeans(returns))
    eff[loop,"sharpe"] <- eff[loop,"Exp.Return"] / eff[loop,"Std.Dev"]
    eff[loop,1:n] <- sol$solution
    loop <- loop+1
  }

  return(as.data.frame(eff))
}

# Run the eff.frontier function based on no short and 50% alloc. restrictions
eff <- eff.frontier(returns=returns$R, short="no", max.allocation=.50,
                    risk.premium.up=1, risk.increment=.001)

##              SPY          EFA          IWM          VWO          LQD
## SPY 6.810936e-04 6.986904e-04 7.919993e-04 7.889955e-04 2.855759e-05
## EFA 6.986904e-04 8.825046e-04 8.166626e-04 9.612606e-04 4.221047e-05
## IWM 7.919993e-04 8.166626e-04 1.083438e-03 9.426306e-04 1.195884e-05
## VWO 7.889955e-04 9.612606e-04 9.426306e-04 1.339355e-03 3.826838e-05
## LQD 2.855759e-05 4.221047e-05 1.195884e-05 3.826838e-05 1.085459e-04
## HYG 3.224162e-04 3.515984e-04 3.832917e-04 3.891801e-04 8.093112e-05
##              HYG
## SPY 3.224162e-04
## EFA 3.515984e-04
## IWM 3.832917e-04
## VWO 3.891801e-04
## LQD 8.093112e-05
## HYG 2.975042e-04

# Find the optimal portfolio
eff.optimal.point <- eff[eff$sharpe==max(eff$sharpe),]

# graph efficient frontier
# Start with color scheme
ealred <- "#7D110C"
ealtan <- "#CDC4B6"
eallighttan <- "#F7F6F0"
ealdark <- "#423C30"

ggplot(eff, aes(x=Std.Dev, y=Exp.Return)) + geom_point(alpha=.1, color=ealdark) +
  geom_point(data=eff.optimal.point, aes(x=Std.Dev, y=Exp.Return, label=sharpe),
             color=ealred, size=5) +
  annotate(geom="text", x=eff.optimal.point$Std.Dev,
           y=eff.optimal.point$Exp.Return,
           label=paste("Risk: ",
                       round(eff.optimal.point$Std.Dev*100, digits=3),"%\nReturn: ",
                       round(eff.optimal.point$Exp.Return*100, digits=4),"%\nSharpe: ",
                       round(eff.optimal.point$sharpe*100, digits=2), "%", sep=""),
           hjust=0, vjust=1.2) +
  ggtitle("Efficient Frontier and Optimal Portfolio") +
  labs(x="Risk (standard deviation of portfolio)", y="Return") +
  theme(panel.background=element_rect(fill=eallighttan),
        text=element_text(color=ealdark),
        plot.title=element_text(size=24, color=ealred))

https://rpubs.com/JanpuHou/258620