Walk Forward Optimization

#*****************************************************************
# Load historical data
#*****************************************************************
library(SIT)
load.packages('quantmod')

# load saved Proxies Raw Data, data.proxy.raw
# please see http://systematicinvestor.github.io/Data-Proxy/ for more details
load('data/data.proxy.raw.Rdata')

tickers = '
EQ = SPY + VFINX # S&P 500
FI = TLT + VUSTX # 20 Year Treasury
'

# uncomment if you want to use same data as in the source
#tickers = 'EQ=SPY, FI=TLT'

data <- new.env()
getSymbols.extra(tickers, src = 'yahoo', from = '1970-01-01', env = data, raw.data = data.proxy.raw, set.symbolnames = T, auto.assign = T)
  for(i in data$symbolnames) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
bt.prep(data, align='remove.na')

# Check data
plota.matplot(scale.one(data$prices),main='Asset Perfromance')

#*****************************************************************
# Setup
#*****************************************************************
prices = data$prices

frequency = 'months'
period.ends = endpoints(prices, frequency)
  period.ends = period.ends[period.ends > 0]

# all possible combinations
choices = expand.grid(
  EQ = seq(0, 100, by = 5),
  FI = seq(0, 100, by = 5),
KEEP.OUT.ATTRS=F)

# only select ones that sum up to 1
choices = choices[choices$EQ + choices$FI == 100,]
  choices = choices[sort.list(choices$EQ),]
  index = 1:nrow(choices)

# run back test over all combinations
result = rep.col(prices[,1], nrow(choices))
  colnames(result) = index

for(i in 1:nrow(choices)) {
  data$weight[] = NA
    data$weight$EQ[period.ends] = choices$EQ[i]/100
    data$weight$FI[period.ends] = choices$FI[i]/100
  model = bt.run.weight.fast(data)
  # uncomment if you want to get same results as in the source
  #model = bt.run.share(data, clean.signal=F, silent=T)
  result[,i] = model$equity
}

plota.matplot(result,main='Strategy Perfromance')

#*****************************************************************
# Pick Strategy based on modified Sharpe over 72 days
#*****************************************************************
sd.factor = 2.5
lookback.len = 72

lookback.return = (result / mlag(result,lookback.len))^(252/lookback.len) - 1

lookback.sd = bt.apply.matrix(result / mlag(result)-1, runSD, lookback.len)*sqrt(252)

mod.sharpe = lookback.return / lookback.sd ^ sd.factor
  mod.sharpe = mod.sharpe[period.ends,]

# pick best one
best.sharpe = ntop(mod.sharpe, 1)

# map back to original weights
weight = t(apply(best.sharpe, 1, function(x) 
  colMeans(choices[index[x!=0],,drop=F])
)) / 100

weight = make.xts(weight, data$dates[period.ends])
#*****************************************************************
# Test strategy
#*****************************************************************
commission = list(cps = 0.01, fixed = 10.0, percentage = 0.0)

models = list()

data$weight[] = NA
	data$weight$EQ = 1
models$SP500 = bt.run.share(data, clean.signal=T, commission = commission, trade.summary=T, silent=T)


data$weight[] = NA
  data$weight[period.ends,] = as.matrix(weight)
models$UIS = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)


#*****************************************************************
# Create Report
#*****************************************************************
plotbt(models, plotX = T, log = 'y', LeftMargin = 3, main = NULL)
	mtext('Cumulative Performance', side = 2, line = 1)

print(plotbt.strategy.sidebyside(models, make.plot=F, return.table=T, perfromance.fn=engineering.returns.kpi))

print(last.trades(models$UIS, make.plot=F, return.table=T))

print(plotbt.monthly.table(models$UIS$equity, make.plot = F))

plota(weight$EQ, type='s', main='SP500 Allocation in UIS Model')

#*****************************************************************
# Modify strategy
#*****************************************************************
# 1. pick top 5
best.sharpe = ntop(mod.sharpe, 5)

  # 2. only consider portfolios with sharpe > 0
  best.sharpe = iif(mod.sharpe > 1, best.sharpe, 0)

# map back to original weights
weight = t(apply(best.sharpe, 1, function(x) 
  colMeans(choices[index[x!=0],,drop=F])
)) / 100

weight = make.xts(weight, data$dates[period.ends])
  # 3. if no portfolio is selected, invest into 100% TLT
  weight = ifna(weight,0)
  weight$FI = weight$FI + 1 - rowSums(weight)


data$weight[] = NA
  data$weight[period.ends,] = as.matrix(weight)
models$UIS5 = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)


#*****************************************************************
# Create Report
#*****************************************************************
plotbt(models, plotX = T, log = 'y', LeftMargin = 3, main = NULL)
	mtext('Cumulative Performance', side = 2, line = 1)

print(plotbt.strategy.sidebyside(models, make.plot=F, return.table=T, perfromance.fn=engineering.returns.kpi))

#*****************************************************************
# make lookback period and the F-Factor adaptive
#*****************************************************************
lookbacks = round(21 * seq(2, 6, by = 0.5))  # 2-6 months range
ffactors  = seq(1.5, 3.5, by = 0.5) # 1.5-3.5 range

# best sharpe across all parameters
best.weight = weight
  best.weight[] = 0
# average of best sharpes for each parameter 
avg.weight = weight
  avg.weight[] = 0
best.mod.sharpe = rep(-10e10, nrow(weight))

for(lookback.len in lookbacks) {
  lookback.return = (result / mlag(result,lookback.len))^(252/lookback.len) - 1
  lookback.sd = bt.apply.matrix(result / mlag(result)-1, runSD, lookback.len)*sqrt(252)
  
  for(sd.factor in ffactors) {
    mod.sharpe = lookback.return / lookback.sd ^ sd.factor
      mod.sharpe = mod.sharpe[period.ends,]

    best.sharpe = ntop(mod.sharpe, 1)

    # map back to original weights
    iweight = t(apply(best.sharpe, 1, function(x) 
      colMeans(choices[index[x!=0],,drop=F])
    )) / 100

    ibest.mod.sharpe = rowSums(best.sharpe * mod.sharpe)

    avg.weight = avg.weight +  iweight
    
    select.index =  ifna(ibest.mod.sharpe > best.mod.sharpe, F)
    best.mod.sharpe[select.index] = ibest.mod.sharpe[select.index]
    best.weight[select.index,] = iweight[select.index,]
  }
}

avg.weight = avg.weight / rowSums(avg.weight)
best.weight = best.weight / rowSums(best.weight)

data$weight[] = NA
  data$weight[period.ends,] = as.matrix(avg.weight)
models$UIS.A = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)


data$weight[] = NA
  data$weight[period.ends,] = as.matrix(best.weight)
models$UIS.B = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)



#*****************************************************************
# Create Report
#*****************************************************************
plotbt(models, plotX = T, log = 'y', LeftMargin = 3, main = NULL)
	mtext('Cumulative Performance', side = 2, line = 1)

print(plotbt.strategy.sidebyside(models, make.plot=F, return.table=T, perfromance.fn=engineering.returns.kpi))

# simple monthly switching, based on prior three month returns, works well
position.score = prices / mlag(prices, 3 * 21)	

data$weight[] = NA
	data$weight[period.ends,] = ntop(position.score[period.ends,], 1)
models$TOP.3M = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)

plotbt(models, plotX = T, log = 'y', LeftMargin = 3, main = NULL)
	mtext('Cumulative Performance', side = 2, line = 1)

print(plotbt.strategy.sidebyside(models, make.plot=F, return.table=T, perfromance.fn=engineering.returns.kpi))

test.models = list()

for(i in 1:12) {
	position.score = prices / mlag(prices, i * 21)	
		
	data$weight[] = NA
		data$weight[period.ends,] = ntop(position.score[period.ends,], 1)
	test.models[[paste0('TOP.', i, 'M')]] = bt.run.share(data, clean.signal=F, commission = commission, trade.summary=T, silent=T)
}

stats = plotbt.strategy.sidebyside(test.models, make.plot=F, return.table=T, perfromance.fn=engineering.returns.kpi)
performance.barchart.helper(stats, 'Sharpe,Cagr,Win.Percent,MaxDD,Volatility,Profit.Factor', c(T,T,T,T,F,T), sort.performance = F)