Load some required library

Download monthly equity market capitalization of 7 countries

(see file: “Betting Against Beta Equity Factors Monthly.xlsx”):

Australia, Canada, France, Germany, Japan, United Kingdom and United States from 1989/12 to 2021/12.

Betting_Against_Beta_Equity_Factors_Monthly_1_ <- read_excel("Downloads/Betting Against Beta Equity Factors Monthly (1).xlsx", 
                                                             sheet = "ME(t-1)", range = "A19:AD1165")
data  <- Betting_Against_Beta_Equity_Factors_Monthly_1_  ########

date <- as.Date(data$DATE, "%m/%d/%Y")


Betting_Against_Beta_Equity_Factors_Monthly_1_ <- xts(coredata(Betting_Against_Beta_Equity_Factors_Monthly_1_[, -1]), order.by = date)

stockpr = subset(Betting_Against_Beta_Equity_Factors_Monthly_1_, select = c(AUS,CAN,FRA,DEU,JPN,GBR,USA) )

Question 1

stockpr1 <- with(stockpr, stockpr[(date >= "1989-11-30" & date <= "2021-12-31")])

monthly.return <- na.omit(Return.calculate(stockpr1, method = "discrete"))

head(monthly.return)
##                    AUS           CAN         FRA         DEU         JPN
## 1989-12-31 -0.02844410  0.0058178215  0.07914015  0.09570552  0.05097275
## 1990-01-31  0.03066617  0.0242062712  0.10677237  0.17590024  0.02123232
## 1990-02-28  0.01574424 -0.0631233515 -0.03985058  0.02206308 -0.04670123
## 1990-03-31 -0.06035407  0.0006305838 -0.02860599  0.00208093 -0.08398405
## 1990-04-30 -0.02950981 -0.0104092150  0.06771805  0.09567357 -0.17803619
## 1990-05-31 -0.08247848 -0.0819545270  0.04144222 -0.06148407 -0.01682893
##                    GBR         USA
## 1989-12-31  0.04729115  0.01203313
## 1990-01-31  0.09773906  0.01799629
## 1990-02-28  0.01176512 -0.07791268
## 1990-03-31 -0.03207884  0.01297222
## 1990-04-30 -0.03012199  0.02484288
## 1990-05-31 -0.06753591 -0.02987430
tail(monthly.return)
##                     AUS          CAN           FRA          DEU           JPN
## 2021-07-31 -0.007808452 -0.004923728 -0.0159977640 -0.005650794  6.743379e-05
## 2021-08-31 -0.006427044 -0.005186644  0.0148548942  0.012951049 -1.218586e-02
## 2021-09-30  0.019011675  0.001517633 -0.0004440893  0.015669145  2.783157e-02
## 2021-10-31 -0.024477429 -0.023212521 -0.0501314877 -0.052374557  1.671376e-02
## 2021-11-30  0.041072770  0.068697882  0.0495422303  0.022591338 -3.342947e-02
## 2021-12-31 -0.058206115 -0.046620684 -0.0319030400 -0.047759221 -3.347900e-02
##                     GBR          USA
## 2021-07-31 -0.033224039  0.026827519
## 2021-08-31  0.025475826  0.009958025
## 2021-09-30  0.006104884  0.026998156
## 2021-10-31 -0.036415768 -0.040433011
## 2021-11-30  0.023093370  0.066567980
## 2021-12-31 -0.054600574 -0.016133549

Question 2

stockpr3 <- with(stockpr, stockpr[(date >= "1993-01-31" & date <= "2021-12-31")])
data <- new.env()
#create 4 required input elements in data
data$prices <- stockpr3
data$weight <- stockpr3
data$execution.price <- stockpr3
data$execution.price[] <- NA
data$symbolnames <- colnames(data$prices)
prices <- data$prices
n = ncol(prices)

data$weight = ntop(prices, n)

model <-list()
model$equal.weight <- bt.run(data, trade.summary=T)
## Latest weights :
##              AUS   CAN   FRA   DEU   JPN   GBR   USA
## 2021-12-31 14.29 14.29 14.29 14.29 14.29 14.29 14.29
## 
## Performance summary :
##  CAGR    Best    Worst   
##  8.2 14.1    -21.7   
capital = 100000
data$weight[] = (capital / prices) * data$weight

####2.1 Compute the equal-weighted portfolio returns EACH month starting from 1993/01 to 2021/12. ####Denote this strategy as the Benchmark portfolio and create its backtesting report using SIT package.

equal.weight = bt.run(data, type='share')
## Latest weights :
##              AUS   CAN   FRA   DEU   JPN   GBR   USA
## 2021-12-31 14.29 14.29 14.29 14.29 14.29 14.29 14.29
## 
## Performance summary :
##  CAGR    Best    Worst   
##  8.2 14.1    -21.7   
head(equal.weight$ret)
##                    AUS
## 1993-01-31 0.000000000
## 1993-02-28 0.003976969
## 1993-03-31 0.044861517
## 1993-04-30 0.065207441
## 1993-05-31 0.035969320
## 1993-06-30 0.024251503
bt.detail.summary(model$equal.weight)
## $System
## $System$Period
## [1] "Jan1993 - Dec2021"
## 
## $System$Cagr
## [1] 8.2
## 
## $System$Sharpe
## [1] 0.57
## 
## $System$DVR
##     [,1]
## AUS 0.51
## 
## $System$Volatility
## [1] 16.02
## 
## $System$MaxDD
## [1] -57.92
## 
## $System$AvgDD
## [1] -7.69
## 
## $System$VaR
##    5% 
## -7.26 
## 
## $System$CVaR
## [1] -11.07
## 
## $System$Exposure
## [1] 99.71
## 
## 
## $Trade
## $Trade$Win.Percent
## [1] 100
## 
## $Trade$Avg.Trade
## [1] 130.2
## 
## $Trade$Avg.Win
## [1] 130.2
## 
## $Trade$Avg.Loss
## [1] NaN
## 
## $Trade$Best.Trade
## [1] 213.68
## 
## $Trade$Worst.Trade
## [1] 25.81
## 
## $Trade$WinLoss.Ratio
## [1] NaN
## 
## $Trade$Avg.Len
## [1] 347
## 
## $Trade$Num.Trades
## [1] 7
## 
## 
## $Period
## $Period$Win.Percent.Day
## [1] 60.3
## 
## $Period$Best.Day
## [1] 14.1
## 
## $Period$Worst.Day
## [1] -21.7
## 
## $Period$Win.Percent.Month
## [1] 60.3
## 
## $Period$Best.Month
## [1] 14.1
## 
## $Period$Worst.Month
## [1] -21.7
## 
## $Period$Win.Percent.Year
## [1] 69
## 
## $Period$Best.Year
## [1] 43
## 
## $Period$Worst.Year
## [1] -48.3
strategy.performance.snapshoot(model, T)

## NULL
data$weight = ntop(prices, n)
head(data$weight)
##                  AUS       CAN       FRA       DEU       JPN       GBR
## 1993-01-31 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
## 1993-02-28 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
## 1993-03-31 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
## 1993-04-30 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
## 1993-05-31 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
## 1993-06-30 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
##                  USA
## 1993-01-31 0.1428571
## 1993-02-28 0.1428571
## 1993-03-31 0.1428571
## 1993-04-30 0.1428571
## 1993-05-31 0.1428571
## 1993-06-30 0.1428571

Question 2.2: Compute MVP portfolio returns by rebalancing EACH month starting from 1993/01 to 2021/12.

Use in-sample data range of 36 months to compute covariance matrix.

Denote this strategy as the MVP portfolio and create its backtesting report using SIT.

data$prices <- stockpr3
data$weight <- stockpr3
data$execution.price <- stockpr3
data$execution.price[] <- NA
prices <- data$prices

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
i = 36
for (i in 36:dim(weight)[1]) {
  hist = ret[ (i- 36 +1):i, ]
  hist = na.omit(hist)
  ia = create.historical.ia(hist, 12)
  ia$cov = cov(coredata(hist))
  weight[i,] = min.risk.portfolio(ia, constraints)
}
## Loading required package: kernlab
## 
## Attaching package: 'kernlab'
## The following object is masked from 'package:SIT':
## 
##     cross
## The following object is masked from 'package:purrr':
## 
##     cross
## The following object is masked from 'package:ggplot2':
## 
##     alpha
data$weight[] = weight   
capital = 100000
data$weight[] = (capital / prices) * data$weight
model$min.var.monthly <- bt.run(data, type = 'share', capital = capital)
## Latest weights :
##               AUS   CAN   FRA   DEU    JPN   GBR  USA
## 2021-12-31 -33.98 21.66 -8.95 15.53 107.47 -4.46 2.72
## 
## Performance summary :
##  CAGR    Best    Worst   
##  4.5 11.4    -15.9   
sum(as.numeric(weight[36,])*as.numeric(ret[37,]))
## [1] 0.0318602
model$min.var.monthly$ret[37, ]
##                  AUS
## 1996-01-31 0.0318602

Question 2.3: Plot both strategies side by side and compare their performance and comment.

plotbt.custom.report.part1(model$min.var.monthly, model$equal.weight)

layout(1:2)
plotbt.transition.map(model$min.var.monthly$weight)
legend('topright', legend = 'min.var.monthly', bty = 'n')
plotbt.transition.map(model$equal.weight$weight)
legend('topright', legend = 'equal weight', bty = 'n')

strategy.performance.snapshoot(model, T)

## NULL
model <- rev(model)
plotbt.custom.report(model)

Question 3

stockpr4 <- with(stockpr, stockpr[(date >= "1988-01-31" & date <= "2021-12-31")])

hist.caps = stockpr4
hist.caps.weight = hist.caps/rowSums(hist.caps)
#Load up aa.test.create.ia.country()
aa.test.create.ia.country <- function(dates = '1990::2021')
{
  # load.packages('quantmod,quadprog')
  symbols = spl('EWA,EWC,EWQ,EWG,EWJ,EWU,SPY')
  symbol.names = spl('Australia, Canada, France, Germany, Japan, UK, USA')
  getSymbols(symbols, from = '1980-01-01', auto.assign = TRUE)
  hist.prices = merge(EWA,EWC,EWQ,EWG,EWJ,EWU,SPY)
  period.ends = endpoints(hist.prices, 'months')
  hist.prices = Ad(hist.prices)[period.ends, ]
  colnames(hist.prices) = symbol.names
  annual.factor = 12
  hist.prices = na.omit(hist.prices[dates])
  hist.returns = na.omit( ROC(hist.prices, type = 'discrete') )
  ia = create.historical.ia(hist.returns, annual.factor)
  return(ia)
}

# 3. Load up efficient frontier plotting function:
plot.ef <- function(
  ia,
  efs,
  portfolio.risk.fn = portfolio.risk,
  transition.map = TRUE,
  layout = NULL
)
{
  risk.label = as.character(substitute(portfolio.risk.fn))
  n = ia$n
  x = match.fun(portfolio.risk.fn)(diag(n), ia)
  y = ia$expected.return
  xlim = range(c(0, x,
                 max( sapply(efs, function(x) max(match.fun(portfolio.risk.fn)(x$weight,ia))) )
  ), na.rm = T)
  ylim = range(c(0, y,
                 min( sapply(efs, function(x) min(portfolio.return(x$weight,ia))) ),
                 max( sapply(efs, function(x) max(portfolio.return(x$weight,ia))) )
  ), na.rm = T)
  x = 100 * x
  y = 100 * y
  xlim = 100 * xlim
  ylim = 100 * ylim
  if( !transition.map ) layout = T
  if( is.null(layout) ) layout(1:2)
  par(mar = c(4,3,2,1), cex = 0.8)
  plot(x, y, xlim = xlim, ylim = ylim,
       xlab='', ylab='', main=paste(risk.label, 'vs Return'), col='black')
  mtext('Return', side = 2,line = 2, cex = par('cex'))
  mtext(risk.label, side = 1,line = 2, cex = par('cex'))
  grid();
  text(x, y, ia$symbols,    col = 'blue', adj = c(1,1), cex = 0.8)
  for(i in len(efs):1) {
    ef = efs[[ i ]]
    x = 100 * match.fun(portfolio.risk.fn)(ef$weight, ia)
    y = 100 * ef$return
    lines(x, y, col=i)
  }
  plota.legend(sapply(efs, function(x) x$name), 1:len(efs))
  if(transition.map) {
    plot.transition.map(efs[[i]]$weight, x, risk.label, efs[[i]]$name)
  }
}
# Use reverse optimization to compute the vector of equilibrium returns
bl.compute.eqret <- function(
  risk.aversion,  # Risk Aversion
  cov,        # Covariance matrix
  cap.weight,     # Market Capitalization Weights
  risk.free = 0   # Rsik Free Interest Rate
)
{
  return( risk.aversion * cov %*% cap.weight +  risk.free)    
}
#  3.1 Visualize Market Capitalization History: 
#  a. Plot Transition of Market Cap Weights in time 
#b. Plot History for each Country's Market Cap
stockpr4 <- with(stockpr, stockpr[(date >= "1988-01-31" & date <= "2021-12-31")])

hist.caps = stockpr4
hist.caps.weight = hist.caps/rowSums(hist.caps)

# Plot Transition of Market Cap Weights in time
plot.transition.map(hist.caps.weight, index(hist.caps.weight), xlab='', name='Market Capitalization Weight History')

# Plot History for each Country's Market Cap
layout( matrix(1:9, nrow = 3, byrow=T) )
col = plota.colors(ncol(hist.caps))
for(i in 1:ncol(hist.caps)) {
  plota(hist.caps[,i], type='l', lwd=5, col=col[i], main=colnames(hist.caps)[i])
}
#3.2 Compute Risk Aversion, prepare Black-Litterman input assumptions
#--------------------------------------------------------------------------
# Compute Risk Aversion, prepare Black-Litterman input assumptions
#--------------------------------------------------------------------------
ia = aa.test.create.ia.country()

# compute Risk Aversion
risk.aversion = bl.compute.risk.aversion( ia$hist.returns$` USA` )
risk.aversion
## [1] 4.663067
# the latest market capitalization weights
cap.weight = last(hist.caps.weight) 

# create Black-Litterman input assumptions  
ia.bl = ia
ia.bl$expected.return = bl.compute.eqret( risk.aversion, ia$cov, as.vector(cap.weight) )

#3.2 Plot market capitalization weights and implied equilibrium returns
layout( matrix(c(1,1,2,3), nrow=2, byrow=T) )

pie(coredata(cap.weight), paste(colnames(cap.weight), round(100*cap.weight), '%'), 
    main = paste('Country Market Capitalization Weights for', format(index(cap.weight),'%b %Y'))
    , col=plota.colors(ia$n))

plot.ia(ia.bl, T)

#--------------------------------------------------------------------------
# Create Efficient Frontier(s)
#--------------------------------------------------------------------------
n = ia$n

# -1 <= x.i <= 1
constraints = new.constraints(n, lb = 0, ub = 1)

# SUM x.i = 1
constraints = add.constraints(rep(1, n), 1, type = '=', constraints)        

# create efficient frontier(s)
ef.risk = portopt(ia, constraints, 50, 'Historical', equally.spaced.risk = T)       
## 
## Attaching package: 'corpcor'
## The following object is masked from 'package:SIT':
## 
##     cov.shrink
ef.risk.bl = portopt(ia.bl, constraints, 50, 'Black-Litterman', equally.spaced.risk = T)    

# Plot multiple Efficient Frontiers and Transition Maps
layout( matrix(1:4, nrow = 2) )
plot.ef(ia, list(ef.risk), portfolio.risk, T, T)            
plot.ef(ia.bl, list(ef.risk.bl), portfolio.risk, T, T) 

##

bl.compute.posterior <- function(
  mu,         # Equilibrium returns
  cov,        # Covariance matrix
  pmat=NULL,  # Views pick matrix
  qmat=NULL,  # Views mean vector
  tau=0.025   # Measure of uncertainty of the prior estimate of the mean returns
)
{
  out = list()    
  omega = diag(c(1,diag(tau * pmat %*% cov %*% t(pmat))))[-1,-1]
  
  temp = solve(solve(tau * cov) + t(pmat) %*% solve(omega) %*% pmat)  
  out$cov = cov + temp
  
  out$expected.return = temp %*% (solve(tau * cov) %*% mu + t(pmat) %*% solve(omega) %*% qmat)
  return(out)
}
#--------------------------------------------------------------------------
# Create Views
#--------------------------------------------------------------------------
temp = matrix(rep(0, n), nrow = 1)
colnames(temp) = ia$symbols

# Relative View
# Japan will outperform UK by 2%
temp[,' Japan'] = 1
temp[,' UK'] = -1


pmat = temp
qmat = c(0.02)

# Absolute View
# Australia's expected return is 12%
temp[] = 0
# temp[,'Australia'] = 1
temp[,1] = 1
pmat = rbind(pmat, temp)    
qmat = c(qmat, 0.12)

# compute posterior distribution parameters
post = bl.compute.posterior(ia.bl$expected.return, ia$cov, pmat, qmat, tau = 0.025 )

# create Black-Litterman input assumptions with Views   
ia.bl.view = ia.bl
ia.bl.view$expected.return = post$expected.return
ia.bl.view$cov = post$cov
ia.bl.view$risk = sqrt(diag(ia.bl.view$cov))

# create efficient frontier(s)
ef.risk.bl.view = portopt(ia.bl.view, constraints, 50, 'Black-Litterman + View(s)', equally.spaced.risk = T)    

# Plot multiple Efficient Frontiers and Transition Maps
layout( matrix(1:4, nrow = 2) )
plot.ef(ia.bl, list(ef.risk.bl), portfolio.risk, T, T)          
plot.ef(ia.bl.view, list(ef.risk.bl.view), portfolio.risk, T, T)