#======================================================
# Use systematic investors toolbox (SIT)
# Load Systematic Investor Toolbox (SIT)
# http://www.r-bloggers.com/backtesting-minimum-variance-portfolios/
# https://systematicinvestor.wordpress.com/2011/12/13/backtesting-minimum-variance-portfolios/
# https://systematicinvestor.wordpress.com/2013/03/22/maximum-sharpe-portfolio/
#=========================================================
#setInternet2(TRUE)
# Load Systematic Investor Toolbox (SIT)
rm(list=ls())
#setInternet2(TRUE)
con = gzcon(url('https://github.com/systematicinvestor/SIT/raw/master/sit.gz', 'rb'))
source(con)
close(con)
#*****************************************************************
# Load historical data
#******************************************************************
load.packages('quantmod,quadprog,lpSolve')
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Attaching package: 'TTR'
## The following object is masked _by_ '.GlobalEnv':
##
## DVI
## Version 0.4-0 included new data defaults. See ?getSymbols.
tickers = spl('SPY,QQQ,EEM,IWM,EFA,TLT,IYR,GLD')
data <- new.env()
getSymbols(tickers, src = 'yahoo', from = '1980-01-01', env = data, auto.assign = T)
## 'getSymbols' currently uses auto.assign=TRUE by default, but will
## use auto.assign=FALSE in 0.5-0. You will still be able to use
## 'loadSymbols' to automatically load data. getOption("getSymbols.env")
## and getOption("getSymbols.auto.assign") will still be checked for
## alternate defaults.
##
## This message is shown once per session and may be disabled by setting
## options("getSymbols.warning4.0"=FALSE). See ?getSymbols for details.
## pausing 1 second between requests for more than 5 symbols
## pausing 1 second between requests for more than 5 symbols
## pausing 1 second between requests for more than 5 symbols
## pausing 1 second between requests for more than 5 symbols
## [1] "SPY" "QQQ" "EEM" "IWM" "EFA" "TLT" "IYR" "GLD"
for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
names(data)
## [1] "GLD" "IWM" "IYR" "TLT" "QQQ"
## [6] ".getSymbols" "EEM" "SPY" "EFA"
data.weekly <- new.env()
for(i in tickers) data.weekly[[i]] = to.weekly(data[[i]], indexAt='endof')
bt.prep(data, align='remove.na', dates='1990::2018')
bt.prep(data.weekly, align='remove.na', dates='1990::2018')
names(data)
## [1] "GLD" "prices" "IWM"
## [4] "dates" "IYR" "TLT"
## [7] "weight" "QQQ" ".getSymbols"
## [10] "symbolnames" "execution.price" "EEM"
## [13] "SPY" "EFA"
#*****************************************************************
# Code Strategies
#******************************************************************
prices = data$prices
n = ncol(prices)
n
## [1] 8
# find week ends
week.ends = endpoints(prices, 'weeks')
week.ends = week.ends[week.ends > 0]
# Equal Weight 1/N Benchmark
data$weight[] = NA
data$weight[week.ends,] = ntop(prices[week.ends,], n)
#capital = 100000
#data$weight[] = (capital / prices) * data$weight
equal.weight = bt.run(data)
## Latest weights :
## EEM EFA GLD IWM IYR QQQ SPY TLT
## 2018-12-31 08:00:00 12.5 12.5 12.5 12.5 12.5 12.5 12.5 12.5
##
## Performance summary :
## CAGR Best Worst
## 8.3 10.3 -8.5
names(equal.weight)
## [1] "weight" "type" "ret" "best" "worst"
## [6] "equity" "cagr" "dates.index"
head(equal.weight$equity)
## EEM
## 2004-11-18 08:00:00 1.000000
## 2004-11-19 08:00:00 1.000000
## 2004-11-22 08:00:00 1.006032
## 2004-11-23 08:00:00 1.008216
## 2004-11-24 08:00:00 1.014851
## 2004-11-26 08:00:00 1.016392
#*****************************************************************
# Create Constraints
#*****************************************************************
constraints = new.constraints(n, lb = -Inf, ub = +Inf)
# SUM x.i = 1
constraints = add.constraints(rep(1, n), 1, type = '=', constraints)
ret = prices / mlag(prices) - 1
weight = coredata(prices)
weight[] = NA
for( i in week.ends[week.ends >= (63 + 1)] ) {
# one quarter is 63 days
hist = ret[ (i- 63 +1):i, ]
# create historical input assumptions
ia = create.historical.ia(hist, 252)
#s0 = apply(coredata(hist),2,sd)
#ia$cov = cor(coredata(hist), use='complete.obs',method='pearson') * (s0 %*% t(s0))
weight[i,] = min.risk.portfolio(ia, constraints)
}
tail(weight,1)
## EEM EFA GLD IWM IYR QQQ
## [3553,] -0.02227524 0.1746793 0.3294391 -0.07523821 0.03104374 -0.244077
## SPY TLT
## [3553,] 0.4325373 0.373891
# Minimum Variance
data$weight[] = weight
capital = 100000
data$weight[] = (capital / prices) * data$weight
min.var.daily = bt.run(data, type='share', capital=capital)
## Latest weights :
## EEM EFA GLD IWM IYR QQQ SPY TLT
## 2018-12-31 08:00:00 -2.79 16.47 32.64 -6.75 2.35 -24.88 44.89 38.07
##
## Performance summary :
## CAGR Best Worst
## 8.8 2.8 -2.6
names(min.var.daily)
## [1] "weight" "type" "share" "capital" "ret"
## [6] "best" "worst" "equity" "cagr" "dates.index"
#
#*****************************************************************
# Code Strategies: Weekly
#******************************************************************
retw = data.weekly$prices / mlag(data.weekly$prices) - 1
weightw = coredata(prices)
weightw[] = NA
week.ends<-week.ends[ -length(week.ends) ]
for( i in week.ends[week.ends >= (63 + 1)] ) {
# map
j = which(index(ret[i,]) == index(retw))
# one quarter = 13 weeks
hist = retw[ (j- 13 +1):j, ]
# create historical input assumptions
ia = create.historical.ia(hist, 52)
#s0 = apply(coredata(hist),2,sd)
#ia$cov = cor(coredata(hist), use='complete.obs',method='pearson') * (s0 %*% t(s0))
weightw[i,] = min.risk.portfolio(ia, constraints)
}
data$weight[] = weightw
capital = 100000
data$weight[] = (capital / prices) * data$weight
min.var.weekly = bt.run(data, type='share', capital=capital)
## Latest weights :
## EEM EFA GLD IWM IYR QQQ SPY TLT
## 2018-12-31 08:00:00 -6.45 50.97 44.24 -26.24 5.46 -10.94 22.15 20.82
##
## Performance summary :
## CAGR Best Worst
## 10.8 4.3 -4.8
#*****************************************************************
# Create Report
#******************************************************************
plotbt.custom.report.part1(min.var.daily, equal.weight)
#
models<-list(min.var.daily, equal.weight)
models<-list("Min.var.daily" = min.var.daily,
"Min.var.weekly" = min.var.weekly,
"Equal.weight" = equal.weight)
#
strategy.performance.snapshoot(models,T)
## NULL
strategy.performance.snapshoot(models, control=list(comparison=T),
sort.performance = T )
plotbt.strategy.sidebyside(models,return.table = T )
## Min.var.daily Min.var.weekly Equal.weight
## Period "Nov2004 - Dec2018" "Nov2004 - Dec2018" "Nov2004 - Dec2018"
## Cagr "8.84" "10.81" "8.33"
## Sharpe "1.21" "1" "0.57"
## DVR "1.17" "0.95" "0.54"
## Volatility "7.22" "10.94" "16.35"
## MaxDD "-13.99" "-16.12" "-44.48"
## AvgDD "-1.1" "-1.74" "-1.67"
## VaR "-0.68" "-1.06" "-1.48"
## CVaR "-1.03" "-1.54" "-2.44"
## Exposure "98.2" "98.2" "99.94"
plotbt.strategy.sidebyside(min.var.daily,return.table = T)
## min.var.daily
## Period "Nov2004 - Dec2018"
## Cagr "8.84"
## Sharpe "1.21"
## DVR "1.17"
## Volatility "7.22"
## MaxDD "-13.99"
## AvgDD "-1.1"
## VaR "-0.68"
## CVaR "-1.03"
## Exposure "98.2"
# plot Daily and Weekly transition maps
layout(1:2)
plotbt.transition.map(min.var.daily$weight)
legend('topright', legend = 'min.var.daily', bty = 'n')
plotbt.transition.map(min.var.weekly$weight)
legend('topright', legend = 'min.var.weekly', bty = 'n')
