- portfolio_noshorts.r
portfolio_noshorts.r
| # portfolio_noshorts.r # # Functions for portfolio analysis # to be used in Introduction to Computational Finance & Financial Econometrics # last updated: August 8, 2012 by Hezky Varon # November 7, 2000 by Eric Zivot # Oct 15, 2003 by Tim Hesterberg # November 18, 2003 by Eric Zivot # November 9, 2004 by Eric Zivot # November 9, 2008 by Eric Zivot # August 11, 2011 by Eric Zivot # # Functions: # 1. efficient.portfolio compute minimum variance portfolio # subject to target return # 2. globalMin.portfolio compute global minimum variance portfolio # 3. tangency.portfolio compute tangency portfolio # 4. efficient.frontier compute Markowitz bullet # 5. getPortfolio create portfolio object stopifnot("package:quadprog" %in% search() || require("quadprog",quietly = TRUE) ) getPortfolio <- function(er, cov.mat, weights) { # contruct portfolio object # # inputs: # er N x 1 vector of expected returns # cov.mat N x N covariance matrix of returns # weights N x 1 vector of portfolio weights # # output is portfolio object with the following elements # call original function call # er portfolio expected return # sd portfolio standard deviation # weights N x 1 vector of portfolio weights # call <- match.call() # # check for valid inputs # asset.names <- names(er) weights <- as.vector(weights) names(weights) = names(er) er <- as.vector(er) # assign names if none exist if(length(er) != length(weights)) stop("dimensions of er and weights do not match") cov.mat <- as.matrix(cov.mat) if(length(er) != nrow(cov.mat)) stop("dimensions of er and cov.mat do not match") if(any(diag(chol(cov.mat)) <= 0)) stop("Covariance matrix not positive definite") # # create portfolio # er.port <- crossprod(er,weights) sd.port <- sqrt(weights %*% cov.mat %*% weights) ans <- list("call" = call, "er" = as.vector(er.port), "sd" = as.vector(sd.port), "weights" = weights) class(ans) <- "portfolio" ans } efficient.portfolio <- function(er, cov.mat, target.return, shorts=TRUE) { # compute minimum variance portfolio subject to target return # # inputs: # er N x 1 vector of expected returns # cov.mat N x N covariance matrix of returns # target.return scalar, target expected return # shorts logical, allow shorts is TRUE # # output is portfolio object with the following elements # call original function call # er portfolio expected return # sd portfolio standard deviation # weights N x 1 vector of portfolio weights # call <- match.call() # # check for valid inputs # asset.names <- names(er) er <- as.vector(er) # assign names if none exist N <- length(er) cov.mat <- as.matrix(cov.mat) if(N != nrow(cov.mat)) stop("invalid inputs") if(any(diag(chol(cov.mat)) <= 0)) stop("Covariance matrix not positive definite") # remark: could use generalized inverse if cov.mat is positive semidefinite # # compute efficient portfolio # if(shorts==TRUE){ ones <- rep(1, N) top <- cbind(2*cov.mat, er, ones) bot <- cbind(rbind(er, ones), matrix(0,2,2)) A <- rbind(top, bot) b.target <- as.matrix(c(rep(0, N), target.return, 1)) x <- solve(A, b.target) w <- x[1:N] } else if(shorts==FALSE){ Dmat <- 2*cov.mat dvec <- rep.int(0, N) Amat <- cbind(rep(1,N), er, diag(1,N)) bvec <- c(1, target.return, rep(0,N)) result <- solve.QP(Dmat=Dmat,dvec=dvec,Amat=Amat,bvec=bvec,meq=2) w <- round(result$solution, 6) } else { stop("shorts needs to be logical. For no-shorts, shorts=FALSE.") } # # compute portfolio expected returns and variance # names(w) <- asset.names er.port <- crossprod(er,w) sd.port <- sqrt(w %*% cov.mat %*% w) ans <- list("call" = call, "er" = as.vector(er.port), "sd" = as.vector(sd.port), "weights" = w) class(ans) <- "portfolio" ans } globalMin.portfolio <- function(er, cov.mat, shorts=TRUE) { # Compute global minimum variance portfolio # # inputs: # er N x 1 vector of expected returns # cov.mat N x N return covariance matrix # shorts logical, allow shorts is TRUE # # output is portfolio object with the following elements # call original function call # er portfolio expected return # sd portfolio standard deviation # weights N x 1 vector of portfolio weights call <- match.call() # # check for valid inputs # asset.names <- names(er) er <- as.vector(er) # assign names if none exist cov.mat <- as.matrix(cov.mat) N <- length(er) if(N != nrow(cov.mat)) stop("invalid inputs") if(any(diag(chol(cov.mat)) <= 0)) stop("Covariance matrix not positive definite") # remark: could use generalized inverse if cov.mat is positive semi-definite # # compute global minimum portfolio # if(shorts==TRUE){ cov.mat.inv <- solve(cov.mat) one.vec <- rep(1,N) w.gmin <- rowSums(cov.mat.inv) / sum(cov.mat.inv) w.gmin <- as.vector(w.gmin) } else if(shorts==FALSE){ Dmat <- 2*cov.mat dvec <- rep.int(0, N) Amat <- cbind(rep(1,N), diag(1,N)) bvec <- c(1, rep(0,N)) result <- solve.QP(Dmat=Dmat,dvec=dvec,Amat=Amat,bvec=bvec,meq=1) w.gmin <- round(result$solution, 6) } else { stop("shorts needs to be logical. For no-shorts, shorts=FALSE.") } names(w.gmin) <- asset.names er.gmin <- crossprod(w.gmin,er) sd.gmin <- sqrt(t(w.gmin) %*% cov.mat %*% w.gmin) gmin.port <- list("call" = call, "er" = as.vector(er.gmin), "sd" = as.vector(sd.gmin), "weights" = w.gmin) class(gmin.port) <- "portfolio" gmin.port } tangency.portfolio <- function(er,cov.mat,risk.free, shorts=TRUE) { # compute tangency portfolio # # inputs: # er N x 1 vector of expected returns # cov.mat N x N return covariance matrix # risk.free scalar, risk-free rate # shorts logical, allow shorts is TRUE # # output is portfolio object with the following elements # call captures function call # er portfolio expected return # sd portfolio standard deviation # weights N x 1 vector of portfolio weights call <- match.call() # # check for valid inputs # asset.names <- names(er) if(risk.free < 0) stop("Risk-free rate must be positive") er <- as.vector(er) cov.mat <- as.matrix(cov.mat) N <- length(er) if(N != nrow(cov.mat)) stop("invalid inputs") if(any(diag(chol(cov.mat)) <= 0)) stop("Covariance matrix not positive definite") # remark: could use generalized inverse if cov.mat is positive semi-definite # # compute global minimum variance portfolio # gmin.port <- globalMin.portfolio(er, cov.mat, shorts=shorts) if(gmin.port$er < risk.free) stop("Risk-free rate greater than avg return on global minimum variance portfolio") # # compute tangency portfolio # if(shorts==TRUE){ cov.mat.inv <- solve(cov.mat) w.t <- cov.mat.inv %*% (er - risk.free) # tangency portfolio w.t <- as.vector(w.t/sum(w.t)) # normalize weights } else if(shorts==FALSE){ Dmat <- 2*cov.mat dvec <- rep.int(0, N) er.excess <- er - risk.free Amat <- cbind(er.excess, diag(1,N)) bvec <- c(1, rep(0,N)) result <- solve.QP(Dmat=Dmat,dvec=dvec,Amat=Amat,bvec=bvec,meq=1) w.t <- round(result$solution/sum(result$solution), 6) } else { stop("Shorts needs to be logical. For no-shorts, shorts=FALSE.") } names(w.t) <- asset.names er.t <- crossprod(w.t,er) sd.t <- sqrt(t(w.t) %*% cov.mat %*% w.t) tan.port <- list("call" = call, "er" = as.vector(er.t), "sd" = as.vector(sd.t), "weights" = w.t) class(tan.port) <- "portfolio" tan.port } efficient.frontier <- function(er, cov.mat, nport=20, alpha.min=-0.5, alpha.max=1.5, shorts=TRUE) { # Compute efficient frontier with no short-sales constraints # # inputs: # er N x 1 vector of expected returns # cov.mat N x N return covariance matrix # nport scalar, number of efficient portfolios to compute # shorts logical, allow shorts is TRUE # # output is a Markowitz object with the following elements # call captures function call # er nport x 1 vector of expected returns on efficient porfolios # sd nport x 1 vector of std deviations on efficient portfolios # weights nport x N matrix of weights on efficient portfolios call <- match.call() # # check for valid inputs # asset.names <- names(er) er <- as.vector(er) N <- length(er) cov.mat <- as.matrix(cov.mat) if(N != nrow(cov.mat)) stop("invalid inputs") if(any(diag(chol(cov.mat)) <= 0)) stop("Covariance matrix not positive definite") # # create portfolio names # port.names <- rep("port",nport) ns <- seq(1,nport) port.names <- paste(port.names,ns) # # compute global minimum variance portfolio # cov.mat.inv <- solve(cov.mat) one.vec <- rep(1, N) port.gmin <- globalMin.portfolio(er, cov.mat, shorts) w.gmin <- port.gmin$weights if(shorts==TRUE){ # compute efficient frontier as convex combinations of two efficient portfolios # 1st efficient port: global min var portfolio # 2nd efficient port: min var port with ER = max of ER for all assets er.max <- max(er) port.max <- efficient.portfolio(er,cov.mat,er.max) w.max <- port.max$weights a <- seq(from=alpha.min,to=alpha.max,length=nport) # convex combinations we.mat <- a %o% w.gmin + (1-a) %o% w.max # rows are efficient portfolios er.e <- we.mat %*% er # expected returns of efficient portfolios er.e <- as.vector(er.e) } else if(shorts==FALSE){ we.mat <- matrix(0, nrow=nport, ncol=N) we.mat[1,] <- w.gmin we.mat[nport, which.max(er)] <- 1 er.e <- as.vector(seq(from=port.gmin$er, to=max(er), length=nport)) for(i in 2:(nport-1)) we.mat[i,] <- efficient.portfolio(er, cov.mat, er.e[i], shorts)$weights } else { stop("shorts needs to be logical. For no-shorts, shorts=FALSE.") } names(er.e) <- port.names cov.e <- we.mat %*% cov.mat %*% t(we.mat) # cov mat of efficient portfolios sd.e <- sqrt(diag(cov.e)) # std devs of efficient portfolios sd.e <- as.vector(sd.e) names(sd.e) <- port.names dimnames(we.mat) <- list(port.names,asset.names) # # summarize results # ans <- list("call" = call, "er" = er.e, "sd" = sd.e, "weights" = we.mat) class(ans) <- "Markowitz" ans } # # print method for portfolio object print.portfolio <- function(x, ...) { cat("Call:\n") print(x$call, ...) cat("\nPortfolio expected return: ", format(x$er, ...), "\n") cat("Portfolio standard deviation: ", format(x$sd, ...), "\n") cat("Portfolio weights:\n") print(round(x$weights,4), ...) invisible(x) } # # summary method for portfolio object summary.portfolio <- function(object, risk.free=NULL, ...) # risk.free risk-free rate. If not null then # compute and print Sharpe ratio # { cat("Call:\n") print(object$call) cat("\nPortfolio expected return: ", format(object$er, ...), "\n") cat( "Portfolio standard deviation: ", format(object$sd, ...), "\n") if(!is.null(risk.free)) { SharpeRatio <- (object$er - risk.free)/object$sd cat("Portfolio Sharpe Ratio: ", format(SharpeRatio), "\n") } cat("Portfolio weights:\n") print(round(object$weights,4), ...) invisible(object) } # hard-coded 4 digits; prefer to let user control, via ... or other argument # # plot method for portfolio object plot.portfolio <- function(object, ...) { asset.names <- names(object$weights) barplot(object$weights, names=asset.names, xlab="Assets", ylab="Weight", main="Portfolio Weights", ...) invisible() } # # print method for Markowitz object print.Markowitz <- function(x, ...) { cat("Call:\n") print(x$call) xx <- rbind(x$er,x$sd) dimnames(xx)[[1]] <- c("ER","SD") cat("\nFrontier portfolios' expected returns and standard deviations\n") print(round(xx,4), ...) invisible(x) } # hard-coded 4, should let user control # # summary method for Markowitz object summary.Markowitz <- function(object, risk.free=NULL) { call <- object$call asset.names <- colnames(object$weights) port.names <- rownames(object$weights) if(!is.null(risk.free)) { # compute efficient portfolios with a risk-free asset nport <- length(object$er) sd.max <- object$sd[1] sd.e <- seq(from=0,to=sd.max,length=nport) names(sd.e) <- port.names # # get original er and cov.mat data from call er <- eval(object$call$er) cov.mat <- eval(object$call$cov.mat) # # compute tangency portfolio tan.port <- tangency.portfolio(er,cov.mat,risk.free) x.t <- sd.e/tan.port$sd # weights in tangency port rf <- 1 - x.t # weights in t-bills er.e <- risk.free + x.t*(tan.port$er - risk.free) names(er.e) <- port.names we.mat <- x.t %o% tan.port$weights # rows are efficient portfolios dimnames(we.mat) <- list(port.names, asset.names) we.mat <- cbind(rf,we.mat) } else { er.e <- object$er sd.e <- object$sd we.mat <- object$weights } ans <- list("call" = call, "er"=er.e, "sd"=sd.e, "weights"=we.mat) class(ans) <- "summary.Markowitz" ans } print.summary.Markowitz <- function(x, ...) { xx <- rbind(x$er,x$sd) port.names <- names(x$er) asset.names <- colnames(x$weights) dimnames(xx)[[1]] <- c("ER","SD") cat("Frontier portfolios' expected returns and standard deviations\n") print(round(xx,4), ...) cat("\nPortfolio weights:\n") print(round(x$weights,4), ...) invisible(x) } # hard-coded 4, should let user control # # plot efficient frontier # # things to add: plot original assets with names # tangency portfolio # global min portfolio # risk free asset and line connecting rf to tangency portfolio # plot.Markowitz <- function(object, plot.assets=FALSE, ...) # plot.assets logical. If true then plot asset sd and er { if (!plot.assets) { y.lim=c(0,max(object$er)) x.lim=c(0,max(object$sd)) plot(object$sd,object$er,type="b",xlim=x.lim, ylim=y.lim, xlab="Portfolio SD", ylab="Portfolio ER", main="Efficient Frontier", ...) } else { call = object$call mu.vals = eval(call$er) sd.vals = sqrt( diag( eval(call$cov.mat) ) ) y.lim = range(c(0,mu.vals,object$er)) x.lim = range(c(0,sd.vals,object$sd)) plot(object$sd,object$er,type="b", xlim=x.lim, ylim=y.lim, xlab="Portfolio SD", ylab="Portfolio ER", main="Efficient Frontier", ...) text(sd.vals, mu.vals, labels=names(mu.vals)) } invisible() } |
Source the r source codes above.
| source("H:/Quantitative Research/R/Optimization/portfolio_noshorts.r") |
Sample ticker code:
SPY SPDR S&P 500 ETF TRUST
EFA ISHARES MSCI EAFE ETF
LQD ISHARES IBOXX INVESTMENT GRA
HYG ISHARES IBOXX HIGH YIELD COR
Weekly basis stats (in decimal numbers):
(1) Risk-free return
| Risk-free | 0.00018631 |
(2) Average Return and Standard Deviation of Returns
| Weekly | SPY | EFA | LQD | HYG |
| Average | 0.00260290 | 0.00120427 | 0.00110844 | 0.00132166 |
| SD | 0.02029120 | 0.02491309 | 0.00788386 | 0.01151264 |
(3) Covariance matrix
| Cov | SPY | EFA | LQD | HYG |
| SPY | 0.00041173 | 0.00045214 | -0.00000957 | 0.00017600 |
| EFA | 0.00045214 | 0.00062066 | -0.00000215 | 0.00021265 |
| LQD | -0.00000957 | -0.00000215 | 0.00006216 | 0.00002826 |
| HYG | -0.00000617 | -0.00000570 | 0.00000461 | -0.00000108 |
Annual basis stats (in decimal numbers):
(1') Geometrical-Average Risk-free return
Annual
| Risk-free | 0.00973449 |
(2') Geometrical-Average Return and Standard Deviation of Returns (SD figures of (2) are multiplied by 52^0.5.)
|
| Cov | SPY | EFA | LQD | HYG |
| SPY | 0.02141010 | 0.02351144 | -0.00049749 | 0.00915216 |
| EFA | 0.02351144 | 0.03227444 | -0.00011169 | 0.01105802 |
| LQD | -0.00049749 | -0.00011169 | 0.00323208 | 0.00146939 |
| HYG | 0.00915216 | 0.01105802 | 0.00146939 | 0.00689213 |
Stats on a annual basis:
| # annual r.free <- 0.00973449 r.free [1] 0.00973449 er <- c( "SPY" = 0.14473728, "EFA" = 0.06458429, "LQD" = 0.05929884, "EEM" = 0.07109433 ) er SPY EFA LQD EEM 0.14473728 0.06458429 0.05929884 0.07109433 covmat <- matrix(c(0.02141010, 0.02351144, -0.00049749, 0.00915216, 0.02351144, 0.03227444, -0.00011169, 0.01105802, -0.00049749, -0.00011169, 0.00323208, 0.00146939, 0.00915216, 0.01105802, 0.00146939, 0.00689213), nrow = 4, ncol = 4, byrow = TRUE) rownames(covmat) <- c("SPY", "EFA", "LQD", "HYG") colnames(covmat) <- c("SPY", "EFA", "LQD", "HYG") covmat SPY EFA LQD HYG SPY 0.02141010 0.02351144 -0.00049749 0.00915216 EFA 0.02351144 0.03227444 -0.00011169 0.01105802 LQD -0.00049749 -0.00011169 0.00323208 0.00146939 HYG 0.00915216 0.01105802 0.00146939 0.00689213 |
Equal-weighted Portfolio
| # [1 1 1 1] rep(1,4) # equally-weighted ew = rep(1,4)/4 ew [1] 0.25 0.25 0.25 0.25
equalWeight.portfolio = getPortfolio(er=er,cov.mat=covmat,weights=ew)
equalWeight.portfolio Call: getPortfolio(er = er, cov.mat = covmat, weights = ew) Portfolio expected return: 0.08492869 Portfolio standard deviation: 0.09777922 Portfolio weights: SPY EFA LQD EEM 0.25 0.25 0.25 0.25 plot(equalWeight.portfolio) # If you want to change the color, then try: # plot(equalWeight.portfolio, col="blue") |
 |
Tangency Portfolio (no short, i.e., long-only)
tan.port <- tangency.portfolio(er, covmat, r.free, shorts=FALSE)
tan.port Call: tangency.portfolio(er = er, cov.mat = covmat, risk.free = r.free, shorts = FALSE) Portfolio expected return: 0.08408077 Portfolio standard deviation: 0.05679287 Portfolio weights: SPY EFA LQD EEM 0.2901 0.0000 0.7099 0.0000 print(tan.port) Call: tangency.portfolio(er = er, cov.mat = covmat, risk.free = r.free, shorts = FALSE) Portfolio expected return: 0.08408077 Portfolio standard deviation: 0.05679287 Portfolio weights: SPY EFA LQD EEM 0.2901 0.0000 0.7099 0.0000 plot(tan.port) |
 |
Minimum Volatility Portfolio (no short, i.e., long-only)
gmin.port <- globalMin.portfolio(er, covmat, shorts=FALSE)
gmin.port Call: globalMin.portfolio(er = er, cov.mat = covmat, shorts = FALSE) Portfolio expected return: 0.07045074 Portfolio standard deviation: 0.05181242 Portfolio weights: SPY EFA LQD EEM 0.1238 0.0000 0.8275 0.0487 attributes(gmin.port) $names [1] "call" "er" "sd" "weights" $class [1] "portfolio" print(gmin.port) Call: globalMin.portfolio(er = er, cov.mat = covmat, shorts = FALSE) Portfolio expected return: 0.07045074 Portfolio standard deviation: 0.05181242 Portfolio weights: SPY EFA LQD EEM 0.1238 0.0000 0.8275 0.0487 plot(gmin.port) |
 |
Target Return Portfolio (no short, i.e., long-only)
target.return = er["SPY"]
target.return SPY 0.1447373 e.port.SPY = efficient.portfolio(er, covmat, target.return, shorts=FALSE) e.port.SPY Call: efficient.portfolio(er = er, cov.mat = covmat, target.return = target.return, shorts = FALSE) Portfolio expected return: 0.1447373 Portfolio standard deviation: 0.1463219 Portfolio weights: SPY EFA LQD EEM 1 0 0 0 plot(e.port.SPY) |
 |
Efficient Frontier (no short, i.e., long-only)
ef = efficient.frontier(er, covmat, alpha.min = -2, alpha.max = 1.5, nport = 20, shorts=FALSE)
attributes(ef) $names [1] "call" "er" "sd" "weights" $class [1] "Markowitz" summary(ef) Frontier portfolios' expected returns and standard deviations port 1 port 2 port 3 port 4 port 5 port 6 port 7 port 8 port 9 port 10 ER 0.0705 0.0744 0.0783 0.0822 0.0861 0.0900 0.0939 0.0978 0.1017 0.1056 SD 0.0518 0.0521 0.0533 0.0554 0.0584 0.0621 0.0665 0.0713 0.0765 0.0820 port 11 port 12 port 13 port 14 port 15 port 16 port 17 port 18 port 19 ER 0.1095 0.1135 0.1174 0.1213 0.1252 0.1291 0.133 0.1369 0.1408 SD 0.0878 0.0938 0.1000 0.1064 0.1128 0.1194 0.126 0.1327 0.1395 port 20 ER 0.1447 SD 0.1463 Portfolio weights: SPY EFA LQD EEM port 1 0.1238 0 0.8275 0.0487 port 2 0.1763 0 0.8237 0.0000 port 3 0.2220 0 0.7780 0.0000 port 4 0.2678 0 0.7322 0.0000 port 5 0.3136 0 0.6864 0.0000 port 6 0.3593 0 0.6407 0.0000 port 7 0.4051 0 0.5949 0.0000 port 8 0.4509 0 0.5491 0.0000 port 9 0.4966 0 0.5034 0.0000 port 10 0.5424 0 0.4576 0.0000 port 11 0.5881 0 0.4119 0.0000 port 12 0.6339 0 0.3661 0.0000 port 13 0.6797 0 0.3203 0.0000 port 14 0.7254 0 0.2746 0.0000 port 15 0.7712 0 0.2288 0.0000 port 16 0.8170 0 0.1830 0.0000 port 17 0.8627 0 0.1373 0.0000 port 18 0.9085 0 0.0915 0.0000 port 19 0.9542 0 0.0458 0.0000 port 20 1.0000 0 0.0000 0.0000 plot(ef, plot.assets=T, col="blue", pch=16) points(gmin.port$sd, gmin.port$er, col="green", pch=16, cex=2) points(tan.port$sd, tan.port$er, col="red", pch=16, cex=2) text(gmin.port$sd, gmin.port$er, labels="GLOBAL MIN", pos=2) text(tan.port$sd, tan.port$er, labels="TANGENCY", pos=2) sr.tan = (tan.port$er - r.free)/tan.port$sd abline(a=r.free, b=sr.tan, col="green", lwd=2) |
 |
No comments:
Post a Comment