hw8

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

# Load Systematic Investor Toolbox (SIT)
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

## 
## ################################### WARNING ###################################
## # We noticed you have dplyr installed. The dplyr lag() function breaks how    #
## # base R's lag() function is supposed to work, which breaks lag(my_xts).      #
## #                                                                             #
## # If you call library(dplyr) later in this session, then calls to lag(my_xts) #
## # that you enter or source() into this session won't work correctly.          #
## #                                                                             #
## # All package code is unaffected because it is protected by the R namespace   #
## # mechanism.                                                                  #
## #                                                                             #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning.  #
## #                                                                             #
## # You can use stats::lag() to make sure you're not using dplyr::lag(), or you #
## # can add conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop   #
## # dplyr from breaking base R's lag() function.                                #
## ################################### WARNING ###################################

## 
## Attaching package: 'TTR'

## The following object is masked _by_ '.GlobalEnv':
## 
##     DVI

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

tickers = spl('SPY,QQQ,EEM,IWM,EFA,TLT,IYR,GLD')


data <- new.env()
getSymbols(tickers, src = 'yahoo', from = '2000-01-01', env = data, auto.assign = T)

## [1] "SPY" "QQQ" "EEM" "IWM" "EFA" "TLT" "IYR" "GLD"

for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)

data.weekly <- new.env()
for(i in tickers) data.weekly[[i]] = to.weekly(data[[i]], indexAt='endof')

bt.prep(data, align='remove.na', dates='2000::2022')
bt.prep(data.weekly, align='remove.na', dates='2000::2022')

#*****************************************************************
# Code Strategies
#****************************************************************** 
prices = data$prices   
n = ncol(prices)

# 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, type='share')

## Latest weights :
##                       EEM   EFA   GLD  IWM   IYR   QQQ   SPY   TLT
## 2022-12-30 08:00:00 12.67 12.55 12.58 12.5 12.52 12.42 12.48 12.28
## 
## Performance summary :
##  CAGR    Best    Worst   
##  8   10.2    -8.8

#*****************************************************************
# 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)
}

# 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
## 2022-12-30 08:00:00 34.25 -51.19 48.19 -82.82 7.65 -93.89 217.12 20.69
## 
## Performance summary :
##  CAGR    Best    Worst   
##  7.9 4.1 -5.2

#*****************************************************************
# Code Strategies: Weekly
#******************************************************************     
retw = data.weekly$prices / mlag(data.weekly$prices) - 1
weightw = coredata(prices)
weightw[] = NA

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
## 2022-12-30 08:00:00 -25.91 -93.66 138.66 -70.47 17.57 -108.08 238.22 3.65
## 
## Performance summary :
##  CAGR    Best    Worst   
##  9.8 7.1 -13.9

#*****************************************************************
# Create Report
#****************************************************************** 
plotbt.custom.report.part1(min.var.weekly, min.var.daily, equal.weight)

# 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')

hw8

2023-05-20

R Markdown

Including Plots