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
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(quantmod)
getSymbols("SSUN.F", auto.assign = TRUE)
## Warning: SSUN.F contains missing values. Some functions will not work if
## objects contain missing values in the middle of the series. Consider using
## na.omit(), na.approx(), na.fill(), etc to remove or replace them.
## [1] "SSUN.F"
tickers = c("SSUN.F", "AAPL",   "LGLG.F")
getSymbols(tickers, from = "2007-01-01", auto.assign = TRUE)
## Warning: SSUN.F contains missing values. Some functions will not work if
## objects contain missing values in the middle of the series. Consider using
## na.omit(), na.approx(), na.fill(), etc to remove or replace them.
## Warning: LGLG.F contains missing values. Some functions will not work if
## objects contain missing values in the middle of the series. Consider using
## na.omit(), na.approx(), na.fill(), etc to remove or replace them.
## [1] "SSUN.F" "AAPL"   "LGLG.F"
data = new.env()
getSymbols(tickers, from = "2007-01-01", env = data , auto.assign = TRUE)
## Warning: SSUN.F contains missing values. Some functions will not work if
## objects contain missing values in the middle of the series. Consider using
## na.omit(), na.approx(), na.fill(), etc to remove or replace them.
## Warning: LGLG.F contains missing values. Some functions will not work if
## objects contain missing values in the middle of the series. Consider using
## na.omit(), na.approx(), na.fill(), etc to remove or replace them.
## [1] "SSUN.F" "AAPL"   "LGLG.F"
ls(data)
## [1] "AAPL"   "LGLG.F" "SSUN.F"
names(data)
## [1] "AAPL"        "SSUN.F"      "LGLG.F"      ".getSymbols"
head(data$AAPL)
##            AAPL.Open AAPL.High AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
## 2007-01-03  3.081786  3.092143 2.925000   2.992857  1238319600      2.543757
## 2007-01-04  3.001786  3.069643 2.993571   3.059286   847260400      2.600218
## 2007-01-05  3.063214  3.078571 3.014286   3.037500   834741600      2.581701
## 2007-01-08  3.070000  3.090357 3.045714   3.052500   797106800      2.594451
## 2007-01-09  3.087500  3.320714 3.041071   3.306071  3349298400      2.809971
## 2007-01-10  3.383929  3.492857 3.337500   3.464286  2952880000      2.944446
for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)

names(data)
## [1] "AAPL"        "SSUN.F"      "LGLG.F"      ".getSymbols"
data.weekly <- new.env()
for(i in tickers) data.weekly[[i]] = to.weekly(data[[i]], indexAt='endof')
## Warning in to.period(x, "weeks", name = name, ...): missing values removed from
## data

## Warning in to.period(x, "weeks", name = name, ...): missing values removed from
## data
data.monthly <- new.env()
for(i in tickers) data.monthly[[i]] = to.monthly(data[[i]], indexAt='endof')
## Warning in to.period(x, "months", indexAt = indexAt, name = name, ...): missing
## values removed from data

## Warning in to.period(x, "months", indexAt = indexAt, name = name, ...): missing
## values removed from data
bt.prep(data, align='remove.na', dates='2010::2018')
bt.prep(data.monthly, align='remove.na', dates='2010::2018')
names(data)
## [1] "prices"          "AAPL"            "SSUN.F"          "dates"          
## [5] "LGLG.F"          "weight"          ".getSymbols"     "symbolnames"    
## [9] "execution.price"
#*****************************************************************
# Code Strategies
#****************************************************************** 
prices = data$prices   
n = ncol(prices)
n
## [1] 3
# find week ends
week.ends = endpoints(prices, 'weeks')
week.ends = week.ends[week.ends > 0]  
month.ends = endpoints(prices, 'months')
month.ends = month.ends[month.ends > 0]  


# Equal Weight 1/N Benchmark
data$weight[] = NA
data$weight[week.ends,] = ntop(prices[week.ends,], n)       
data$weight[month.ends,] = ntop(prices[month.ends,], n) 
#capital = 100000
#data$weight[] = (capital / prices) * data$weight
equal.weight = bt.run(data)
## Latest weights :
##             AAPL LGLG.F SSUN.F
## 2018-12-28 33.33  33.33  33.33
## 
## Performance summary :
##  CAGR    Best    Worst   
##  15.9    5.8 -6.3    
names(equal.weight)
## [1] "weight"      "type"        "ret"         "best"        "worst"      
## [6] "equity"      "cagr"        "dates.index"
head(equal.weight$equity)
##                 AAPL
## 2010-01-04 1.0000000
## 2010-01-05 1.0000000
## 2010-01-06 1.0000000
## 2010-01-07 1.0000000
## 2010-01-08 1.0000000
## 2010-01-11 0.9858916
#*****************************************************************
# 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 month.ends[month.ends >= (756 + 1)] ) {
  # one quarter is 63 days
  hist = ret[ (i- 756 +1):i, ]
  
  # create historical input assumptions
  ia = create.historical.ia(hist, 756)
  #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,10)
##              AAPL    LGLG.F    SSUN.F
## [2216,]        NA        NA        NA
## [2217,]        NA        NA        NA
## [2218,]        NA        NA        NA
## [2219,]        NA        NA        NA
## [2220,]        NA        NA        NA
## [2221,]        NA        NA        NA
## [2222,]        NA        NA        NA
## [2223,]        NA        NA        NA
## [2224,]        NA        NA        NA
## [2225,] 0.5791782 0.1608686 0.2599532
# Minimum Variance
data$weight[] = weight      
#capital = 100000
#data$weight[] = (capital / prices) * data$weight
min.var.daily = bt.run(data)
## Latest weights :
##             AAPL LGLG.F SSUN.F
## 2018-12-28 58.83  16.26  24.91
## 
## Performance summary :
##  CAGR    Best    Worst   
##  12.3    4.6 -5.3    
names(min.var.daily)
## [1] "weight"      "type"        "ret"         "best"        "worst"      
## [6] "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 month.ends[month.ends >= (756 + 1)] ) {   
  # map
  j = which(index(ret[i,]) == index(retw))
  
  # one quarter = 13 weeks
  #hist = retw[ (j- 156 +1):j, ]
  
  # create historical input assumptions
  ia = create.historical.ia(hist, n)
  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 :
##             AAPL LGLG.F SSUN.F
## 2018-12-28 56.75  16.52  26.72
## 
## Performance summary :
##  CAGR    Best    Worst   
##  11.6    4.8 -4.8    
#*****************************************************************
# Code Strategies: Monthly
#******************************************************************    
retm = data.monthly$prices / mlag(data.monthly$prices) - 1
weightm = coredata(prices)
weightm[] = NA

month.ends<-month.ends[-length(month.ends)]
  

data$weight[] = weightm     
capital = 100000
data$weight[] = (capital / prices) * data$weight
min.var.monthly = bt.run(data, type='share', capital=capital)
## Latest weights :
##            AAPL LGLG.F SSUN.F
## 2018-12-28    0      0      0
## 
## Performance summary :
##  CAGR    Best    Worst   
##  0   0   0   
#*****************************************************************
# Create Report
#****************************************************************** 
plotbt.custom.report.part3(min.var.daily, min.var.weekly, min.var.monthly, equal.weight)
#
models<-list("Min.var.daily" = min.var.daily, 
             "Min.var.weekly" = min.var.weekly,
             "Min.var.monthly" = min.var.monthly,
             "Equal.weight" = equal.weight)
#
strategy.performance.snapshoot(models, T)
## Warning in cor(y, x): the standard deviation is zero
## Warning in min(drawdown[x[1]:x[2]], na.rm = T): no non-missing arguments to
## min; returning Inf

## NULL
strategy.performance.snapshoot(models, control=list(comparison=T), 
                               sort.performance=T)
## Warning in cor(y, x): the standard deviation is zero
## Warning in min(drawdown[x[1]:x[2]], na.rm = T): no non-missing arguments to
## min; returning Inf

plotbt.strategy.sidebyside(models, return.table=T)
## Warning in cor(y, x): the standard deviation is zero
## Warning in min(drawdown[x[1]:x[2]], na.rm = T): no non-missing arguments to
## min; returning Inf

##            Min.var.daily       Min.var.weekly      Min.var.monthly    
## Period     "Jan2010 - Dec2018" "Jan2010 - Dec2018" "Jan2010 - Dec2018"
## Cagr       "12.25"             "11.58"             "0"                
## Sharpe     "0.83"              "0.79"              "NaN"              
## DVR        "0.68"              "0.66"              "NaN"              
## Volatility "15.69"             "15.55"             "0"                
## MaxDD      "-29.33"            "-28.21"            "0"                
## AvgDD      "-3.07"             "-3.29"             "NaN"              
## VaR        "-1.6"              "-1.53"             "0"                
## CVaR       "-2.33"             "-2.33"             "NaN"              
## Exposure   "65.75"             "65.75"             "0"                
##            Equal.weight       
## Period     "Jan2010 - Dec2018"
## Cagr       "15.93"            
## Sharpe     "0.79"             
## DVR        "0.69"             
## Volatility "22.04"            
## MaxDD      "-32.66"           
## AvgDD      "-3.89"            
## VaR        "-2.12"            
## CVaR       "-2.98"            
## Exposure   "99.78"
plotbt.strategy.sidebyside(min.var.daily, return.table=T)

##            min.var.daily      
## Period     "Jan2010 - Dec2018"
## Cagr       "12.25"            
## Sharpe     "0.83"             
## DVR        "0.68"             
## Volatility "15.69"            
## MaxDD      "-29.33"           
## AvgDD      "-3.07"            
## VaR        "-1.6"             
## CVaR       "-2.33"            
## Exposure   "65.75"
# plot Daily and Weekly transition maps
layout(1:3)
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')

plotbt.transition.map(min.var.monthly$weight)
legend('topright', legend = 'min.var.monthy', bty = 'n')