globMin<-function(cov,e=NULL){
call<-match.call()
if(is.null(e)){e<-rep(0.05,dim(cov)[2])}
#label the assets
if(is.null(names(e))){
asset.names <- rep("asset",length(e))
pnumb <- c(1:length(e))
asset.names <- paste(asset.names,pnumb)
}else{asset.names<-names(e)}
#inverse matrix
cov.inv<-solve(cov)
#solve the global min variance portfolio
w<-as.vector(rowSums(cov.inv)/sum(cov.inv))
names(w)<-asset.names
eR<-crossprod(w,e)
sd<-sqrt(t(w)%*%cov%*%w)
globMin<-list("call"=call,"eR"=eR,"sd"=sd,"weights"=w)
class(globMin)<-"portfolio"
return(globMin)
}
getPortf<-function(cov,e,w){
if(!matrixcalc::is.positive.semi.definite(cov)){stop("Covariance matrix has to be PSD")}
n<-dim(cov)[1]
#label the assets
if(is.null(names(e))){
asset.names <- rep("asset",length(e))
pnumb <- c(1:length(e))
asset.names <- paste(asset.names,pnumb)
}else{asset.names<-names(e)}
e<-matrix(e,nrow=n)
w<-matrix(w,nrow=n)
names(w)<-asset.names
er<-crossprod(e,w)
sd<-sqrt(t(w)%*%cov%*%w)
pfolio<-list("eR"=er,"sd"=sd,"weights"=w)
class(pfolio)<-"portfolio"
return(pfolio)
}
effPortf<-function(cov=NULL,e=NULL,target=NULL){
call<-match.call()
if(is.null(cov)){cov<-simCov(20)$cov}
if(is.null(e)){e<-(runif(nrow(cov))*0.1)}
if(is.null(target)){target<-max(e)}
n<-dim(cov)[1]
e<-matrix(e,ncol=1)
#label the assets
if(is.null(names(e))){
asset.names <- rep("asset",length(e))
pnumb <- c(1:length(e))
asset.names <- paste(asset.names,pnumb)
}else{asset.names<-names(e)}
unit<-rep(1,n)
#compose matrix A then solve system
top<-cbind(2*cov,e,unit)
bottom<-cbind(rbind(t(e),t(unit)),matrix(0,ncol=2,nrow=2))
A<-rbind(top,bottom)
#prepare right handside
b0<-matrix(c(rep(0,n),target,1),ncol=1)
w<-solve(A,b0)[1:n]
#compute characteristics
er<-crossprod(e,w)
sd<-sqrt(t(w)%*%cov%*%w)
pfolio<-list("call"=call,"eR"=er,"sd"=sd,"weights"=w,"cov"=cov)
class(pfolio)<-"portfolio"
return(pfolio)
}
tanPortf <-
function(cov,e,risk.free)
{
call <- match.call()
#
#label the assets
if(is.null(names(e))){
asset.names <- rep("asset",length(e))
pnumb <- c(1:length(e))
asset.names <- paste(asset.names,pnumb)
}else{asset.names<-names(e)}
# remark: could use generalized inverse if cov.mat is positive semi-definite
#
# compute global minimum variance portfolio
#
gmin.port <- globMin(cov,e)
if(gmin.port$eR < risk.free)
stop("Risk-free rate cannot be greater than expected return on global minimum variance portfolio")
#
# compute tangency portfolio
#
cov.inv <- solve(cov)
w.t <- cov.inv %*% (e - risk.free) # tangency portfolio
w.t <- as.vector(w.t/sum(w.t)) # normalize weights
names(w.t) <- asset.names
er.t <- crossprod(w.t,e)
sd.t <- sqrt(t(w.t) %*% cov %*% 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
}
effFront <-
function(cov=NULL,e=NULL, nport=20, alpha.min=-0.1, alpha.max=1.1)
{
call<-match.call()
#
#introduce for "Ser class"
if(is.null(cov)){cov<-simCov(20)$cov}
if(is.null(e)){e<-(runif(nrow(cov))*0.1)}
#
n<-dim(cov)[1]
e<-matrix(as.vector(e),nrow=n)
#
# enumarate portfolio names
#
port.names <- rep("port",nport)
pnumb <- c(1:nport)
port.names <- paste(port.names,pnumb)
if(is.null(names(e))){
asset.names <- rep("asset",length(e))
pnumb <- c(1:length(e))
asset.names <- paste(asset.names,pnumb)
}else{asset.names<-names(e)}
#
# compute global minimum varian
port.gmin <- globMin(cov,e)
w.gmin <- port.gmin$weights
#
# 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
#
e.max <- max(e)
port.max <- effPortf(cov,e,e.max)
w.max <- port.max$weights
#create convex combination sequence
a <- seq(from=alpha.min,to=alpha.max,length=nport)
#can spare transpose by using outter product
# all(a%*%t(w.gmin)==a%o%w.gmin) is TRUE
we.mat <- a %o% w.gmin + (1-a) %o% w.max # rows are efficient portfolios
er.efp <- c(we.mat %*% e) # expected returns of efficient portfolios
#we need var-covar matrix of efficient portofilos!
#the diagonal are variances of portfolios.
cov.efp <- we.mat %*% cov %*% t(we.mat)
sd.efp <- c(sqrt(diag(cov.efp))) # std devs of efficient frontier portfolios
dimnames(we.mat) <- list(port.names,asset.names)
#
# summarize results
#
ans <- list("call" = call,
"eR" = er.efp,
"sd" = sd.efp,
"cov"=cov,
"weights" = we.mat)
class(ans) <- "EF"
ans
}
#
# print method for portfolio object
print.portfolio <- function(x, ...)
{
cat("Calucalated:\n")
print(x$call, ...)
cat("\nExpected return: ", format(x$eR, ...), "\n")
cat("Standard deviation: ", format(x$sd, ...), "\n")
cat("Weights:\n")
print(round(x$weights,4), ...)
invisible(x)
}
#
# summary method
summary.portfolio <- function(object, risk.free=NULL, ...)
# risk-free rate. If not null then
# compute and print Sharpe ratio
#
{
cat("Calculated:\n")
print(object$call)
cat("\nExpected return: ", format(object$eR, ...), "\n")
cat( "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("Weights:\n")
print(round(object$weights,4), ...)
invisible(object)
}
# hard-coded 4 digits; prefer to let user control, via ... or other argument
#
# plot weights for portfolio object
plot.portfolio<-function(object, ...){
a<-data.frame(weight=object$weights,x=names(object$weights))
ggplot(a,aes(x,weight))+geom_bar(stat="identity")
}
#
# print method for Markowitz object
print.EF <- function(x, ...)
{
cat("Calculated:\n")
print(x$call)
xx <- rbind(x$eR,x$sd)
row.names(xx) <- c("ER","SD")
colnames(xx)<-dimnames(x$weights)[[1]]
cat("\nFrontier expected returns and standard deviations\n")
print(round(xx,4), ...)
invisible(x)
}
#
# summary method for Efficient frontier
summary.EF <- 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 <- eval(object$call$cov)
#
# compute tangency portfolio
tan.port <- tanPortf(cov,e,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.EF"
ans
}
print.summary.EF <- 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 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.EF<-function(object,details=FALSE){
cols<-c("#00C8E0", "#EE0049")
a<-data.frame(sd=object$sd,eR=object$eR)
if(details){
min.p<-as.data.frame(globMin(eval(object$cov),eval(object$call$e)))
max.p<-as.data.frame(effPortf(eval(object$call$cov),eval(object$call$e),max(eval(object$call$e))))
p<-ggplot(a,aes(sd,eR))+geom_point()+geom_point(data=min.p,aes(sd,eR,color="Min Variance Portfolio"),color=cols[2],show.legend= TRUE,size=2.5)+geom_point(data=max.p,aes(sd,eR,color="Max eR"),color=cols[1],show.legend = T,size=2.5)
p<-p+scale_colour_manual(name="legend", values=c("blue","red"))+ggtitle("Portfolio frontier")
return(p)
}
ggplot(a,aes(sd,eR))+geom_point()
}
as.data.frame.portfolio<-function(object){
data.frame(sd=object$sd,eR=object$eR)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.