How to get the weights on the tangency portfolio using R package PortfolioAnalytics

Question

I am trying to use the R PortfolioAnalytics package to compute the weights of the tangency portfolio for the efficient frontier when there is access to a risk free asset.

An interesting and replicable experiment is shown below using the edhec dataset:

library(ROI)
library(PortfolioAnalytics)

data(edhec)
R <- edhec[, 1:5]
# change the column names for better legends in plotting
colnames(R) <- c("CA", "CTAG", "DS", "EM", "EQM")
funds <- colnames(R)

# Mean returns for each asset
(assetMeanReturns <- colMeans(R))

(rf <- min(assetMeanReturns) / 2.0)

# Initialise the portfolio specification using names of assets
portfolio <- portfolio.spec(assets=funds)
portfolio <- add.constraint(portfolio, type="full_investment")
portfolio <- add.constraint(portfolio=portfolio, type="box", min=-1, max=1)
portfolio <- add.objective(portfolio=portfolio, type="risk", name="StdDev")

iterations <-20
results <- matrix(nrow = iterations, ncol = 8)

for (i in 1:iterations) {
  efficientFrontier <- create.EfficientFrontier(R=R, portfolio=portfolio, type="mean-sd", n.portfolios=4+i)
  (means <- efficientFrontier$frontier[,c("mean")] * 100)
  (stddevs <- efficientFrontier$frontier[,c("StdDev")] * 100)
  (sharpeRatios <- (means - rf*100) / stddevs)
  (tangencyPortfolio <- efficientFrontier$frontier[sharpeRatios == max(sharpeRatios),])
  tangencyPortfolioDetails <- efficientFrontier$frontier[sharpeRatios == max(sharpeRatios),1:ncol(efficientFrontier$frontier)]
  results[i,] <- efficientFrontier$frontier[sharpeRatios == max(sharpeRatios),1:ncol(efficientFrontier$frontier)]  
}

results

This code computes portfolios along the efficient frontier with equal distance between steps along the vertical expected return axis, and then choose the one that maximises the sharpe ratio based on the mean and the standard deviation of each efficient frontier portfolio.

It puzzles me that there is not a more direct way to get this portfolio and

This approach does not find the actual tangency portfolio. Instead, it flags the portfolio along the efficient frontier that is closest to being the tangent portfolio as the actual tangent portfolio, given the predetermined E(return) step sizes chosen when computing the efficient frontier.

The variation in the tangent portfolio weights, as you increase the E(Return) step size along the efficient portfolio shows how imprecise this approach is.

Right now I am just going to optimise the portfolio given an expected return and maximising the sharpe ratio directly given the risk free rate of return as per the following code snippet:

# Define the sharpe ratio calculation function - that nests a portfolio optimisation problem.
frontierPortfolioSharpeRatio <- function(expectedReturn) {
  testPortfolio <- add.constraint(portfolio=portfolio, type="return", return_target=expectedReturn)
  optimisedTestPortfolio <- optimize.portfolio(R=R, portfolio=testPortfolio, optimize_method="ROI", trace=TRUE)
  return <- t(optimisedTestPortfolio$weights) %*% assetMeanReturns
  sharpeRatio <- (return - rf) /optimisedTestPortfolio$objective_measures$StdDev
  -sharpeRatio
}

# Get the sharpe ratio of portfolio that minimises variance given a required rate of return
# Initial guess of tangency portfolio required return is minimum variance portfolio expected return. Try others.
results <- optim(mean(assetMeanReturns), frontierPortfolioSharpeRatio, method="BFGS")

tangencyPortfolioRequiredReturn <- results$par
tangencyPortfolioDefinition <- add.constraint(portfolio=portfolio, type="return", return_target=tangencyPortfolioRequiredReturn)
tangencyPortfolio <- optimize.portfolio(R=R, portfolio=tangencyPortfolioDefinition, optimize_method="ROI", trace=TRUE)
(tangencyPortfolioReturn <- t(tangencyPortfolio$weights) %*% assetMeanReturns)
tangencyPortfolio$objective_measures
(sharpeRatio <- (tangencyPortfolioReturn - rf) /tangencyPortfolio$objective_measures$StdDev)
tangencyPortfolio$weights
plot(tangencyPortfolio, risk.col="StdDev", return.col="mean", main="Tangency Portfolio", chart.assets=TRUE)

Any insight on how to be doing this more directly with PortfolioAnalytics would be appreciated.

Geoff S