FINAL PROJECT
Nicole Jancelle Echague
6/18/2021
# https://palomar.home.ece.ust.hk/MAFS6010R_lectures/Rsession_factor_models.html#final_comparison_of_covariance_matrix_estimations_via_different_factor_models
# https://systematicinvestor.wordpress.com/2012/03/19/backtesting-asset-allocation-portfolios/rm(list = ls())library(lubridate)##
## Attaching package: 'lubridate'## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, unionlibrary(timetk)
library(purrr)
library(quantmod)## Loading required package: xts## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric## Loading required package: TTR## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(lpSolve)#devtools::install_github(‘joshuaulrich/xts’, force = T) #devtools::install_github(‘joshuaulrich/quantmod’, force = T) # #devtools::install_github(‘systematicinvestor/SIT.date’, force = T) #curl::curl_download(‘https://github.com/systematicinvestor/SIT/raw/master/SIT.tar.gz’, #‘SIT.tar.gz’,mode = ‘wb’,quiet=T) #install.packages(‘SIT.tar.gz’, repos = NULL, type=‘source’)library(SIT)## Loading required package: SIT.date##
## Attaching package: 'SIT'## The following object is masked from 'package:TTR':
##
## DVI## The following object is masked from 'package:purrr':
##
## cross## The following object is masked from 'package:base':
##
## closelibrary(PerformanceAnalytics)##
## Attaching package: 'PerformanceAnalytics'## The following object is masked from 'package:graphics':
##
## legendlibrary(xts)
library(tidyr)#Import Datatickers = spl('MA,HAS,JNJ,VRTX,AMZN,NFLX,MSFT,MNST')
datas <- new.env()
getSymbols(tickers, src = 'yahoo', from = '2005-01-01', env = datas, 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] "MA" "HAS" "JNJ" "VRTX" "AMZN" "NFLX" "MSFT" "MNST"for(i in ls(datas)) datas[[i]] = adjustOHLC(datas[[i]], use.Adjusted=T)
names(datas)## [1] "VRTX" "HAS" ".getSymbols" "NFLX" "MA"
## [6] "AMZN" "MSFT" "MNST" "JNJ"datas$EEM## NULLhead(datas$EEM)## NULLbt.prep(datas, align='remove.na')
names(datas)## [1] "VRTX" "prices" "dates" "HAS"
## [5] "weight" ".getSymbols" "symbolnames" "execution.price"
## [9] "NFLX" "MA" "AMZN" "MSFT"
## [13] "MNST" "JNJ"#
head(datas$SPY)## NULLhead(datas$prices)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2006-05-25 35.63 12.31116 38.98571 4.289672 7.580000 17.29757 4.107143 33.23
## 2006-05-26 36.07 12.39033 39.14048 4.189889 7.632500 17.28300 4.118571 33.75
## 2006-05-30 34.64 12.20560 38.59875 4.103165 7.471250 16.86767 3.928571 32.62
## 2006-05-31 34.61 12.23199 38.83738 4.190823 7.702917 16.50336 3.955714 34.50
## 2006-06-01 35.07 12.33095 39.11469 4.430486 7.656667 16.62723 4.072857 34.35
## 2006-06-02 34.76 12.42991 39.17275 4.381993 7.678333 16.58352 3.992857 34.67prices_monthly <- to.monthly(datas$prices, indexAt = "last", OHLC = FALSE) # indexAt = 'endof'
head(prices_monthly)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2006-05-31 34.61 12.23199 38.83738 4.190823 7.702917 16.50336 3.955714 34.50
## 2006-06-30 38.68 11.94828 38.64388 4.476179 7.932083 16.97697 3.887143 36.71
## 2006-07-31 26.89 12.41766 40.34005 4.277548 7.665000 17.53073 2.955714 33.52
## 2006-08-31 30.83 13.48013 41.94282 5.212882 4.590000 18.79462 2.860000 34.45
## 2006-09-29 32.12 15.10705 42.12446 6.560400 5.413333 20.00129 3.254286 33.65
## 2006-10-31 38.09 17.29476 43.72017 6.919045 5.291667 20.99586 3.951429 40.60monthly.ret <- na.omit(Return.calculate(prices_monthly, method = "discrete"))
head(monthly.ret)## AMZN HAS JNJ MA MNST
## 2006-06-30 0.11759604 -0.02319345 -0.004982313 0.06809068 0.02975055
## 2006-07-31 -0.30480871 0.03928422 0.043892354 -0.04437512 -0.03367123
## 2006-08-31 0.14652291 0.08556095 0.039731504 0.21866125 -0.40117417
## 2006-09-29 0.04184233 0.12069001 0.004330586 0.25849770 0.17937538
## 2006-10-31 0.18586554 0.14481394 0.037880769 0.05466816 -0.02247525
## 2006-11-30 0.05907062 0.03202161 -0.016583650 0.37314442 -0.11401587
## MSFT NFLX VRTX
## 2006-06-30 0.02869766 -0.01733467 0.06405794
## 2006-07-31 0.03261806 -0.23961789 -0.08689728
## 2006-08-31 0.07209598 -0.03238270 0.02774466
## 2006-09-29 0.06420299 0.13786224 -0.02322203
## 2006-10-31 0.04972508 0.21422303 0.20653776
## 2006-11-30 0.02613658 0.05856818 0.09113303#download Fama-French factors from website
url <- "http://mba.tuck.dartmouth.edu/pages/faculty/ken.french/ftp/F-F_Research_Data_Factors_daily_CSV.zip"
temp <- tempfile()
download.file(url, temp, method = "libcurl", mode = "wb")
unzip(temp, "F-F_Research_Data_Factors_daily.CSV")
unlink(temp)mydata <- read.csv("F-F_Research_Data_Factors_daily.CSV", skip = 4)
mydata <- mydata[-nrow(mydata), ] # remove last row
fama_lib <- xts(x = mydata[, c(2,3,4)], order.by = as.Date(paste(mydata[, 1]), "%Y%m%d"))
str(fama_lib)## An 'xts' object on 1926-07-01/2021-04-30 containing:
## Data: num [1:24978, 1:3] 0.1 0.45 0.17 0.09 0.21 -0.71 0.62 0.04 0.48 0.04 ...
## - attr(*, "dimnames")=List of 2
## ..$ : NULL
## ..$ : chr [1:3] "Mkt.RF" "SMB" "HML"
## Indexed by objects of class: [Date] TZ: UTC
## xts Attributes:
## NULLhead(fama_lib)## Mkt.RF SMB HML
## 1926-07-01 0.10 -0.24 -0.28
## 1926-07-02 0.45 -0.32 -0.08
## 1926-07-06 0.17 0.27 -0.35
## 1926-07-07 0.09 -0.59 0.03
## 1926-07-08 0.21 -0.36 0.15
## 1926-07-09 -0.71 0.44 0.56dates <- '2005-01::2021-03'
X <- monthly.ret[dates]
head(X)## AMZN HAS JNJ MA MNST
## 2006-06-30 0.11759604 -0.02319345 -0.004982313 0.06809068 0.02975055
## 2006-07-31 -0.30480871 0.03928422 0.043892354 -0.04437512 -0.03367123
## 2006-08-31 0.14652291 0.08556095 0.039731504 0.21866125 -0.40117417
## 2006-09-29 0.04184233 0.12069001 0.004330586 0.25849770 0.17937538
## 2006-10-31 0.18586554 0.14481394 0.037880769 0.05466816 -0.02247525
## 2006-11-30 0.05907062 0.03202161 -0.016583650 0.37314442 -0.11401587
## MSFT NFLX VRTX
## 2006-06-30 0.02869766 -0.01733467 0.06405794
## 2006-07-31 0.03261806 -0.23961789 -0.08689728
## 2006-08-31 0.07209598 -0.03238270 0.02774466
## 2006-09-29 0.06420299 0.13786224 -0.02322203
## 2006-10-31 0.04972508 0.21422303 0.20653776
## 2006-11-30 0.02613658 0.05856818 0.09113303dim(X)## [1] 178 8#
fama_lib_month <- to.monthly(fama_lib, indexAt = "last", OHLC = FALSE) # indexAt = 'endof'
head(fama_lib_month)## Mkt.RF SMB HML
## 1926-07-31 0.46 -0.12 -0.17
## 1926-08-31 0.58 -0.52 -0.04
## 1926-09-30 0.30 -0.22 0.06
## 1926-10-30 -0.25 0.20 0.25
## 1926-11-30 0.26 0.03 -0.50
## 1926-12-31 0.33 0.30 0.67f <- fama_lib_month[dates]/100
head(f)## Mkt.RF SMB HML
## 2005-01-31 0.0097 0.0076 0.0014
## 2005-02-28 -0.0061 0.0027 0.0007
## 2005-03-31 0.0000 -0.0006 0.0065
## 2005-04-29 0.0104 -0.0040 -0.0003
## 2005-05-31 -0.0049 0.0046 0.0033
## 2005-06-30 -0.0063 0.0021 0.0025dim(f)## [1] 195 3#======================================================================
# Use Professor Palomar's package to compute estimated coefficients
# Based on CAPM model, compute covariance matrix for the 8-asset portfolio by using past 60-month returns from 2010/01 - 2014/12.
#======================================================================
devtools::install_github("dppalomar/covFactorModel")## Error in get(genname, envir = envir) : object 'testthat_print' not found## Skipping install of 'covFactorModel' from a github remote, the SHA1 (1e7cc3df) has not changed since last install.
## Use `force = TRUE` to force installationlibrary(covFactorModel)
#
insample_range <- '2010-01::2014-12'
one_factor_model <- factorModel(X[insample_range], type = "M", econ_fact = f$Mkt.RF[insample_range])
names(one_factor_model)## [1] "alpha" "beta" "factors" "residual"cbind(alpha = one_factor_model$alpha, beta = one_factor_model$beta)## alpha Mkt.RF
## AMZN 0.01697563 -0.2085693
## HAS 0.01356192 -0.2763406
## JNJ 0.01239375 1.2967969
## MA 0.02358281 1.2943668
## MNST 0.03401129 1.6321305
## MSFT 0.01124666 0.3790205
## NFLX 0.05205511 1.2534143
## VRTX 0.02456127 -0.5415507#======================================================================
# we can do the fitting using a more compact matrix notation
F_ <- cbind(ones = 1, f$Mkt.RF[insample_range])
Gamma <- t(X[insample_range]) %*% F_ %*% solve(t(F_) %*% F_) # better: Gamma <- t(solve(t(F_) %*% F_, t(F_) %*% X))
colnames(Gamma) <- c("alpha", "beta")
alpha <- Gamma[, 1] # or alpha <- Gamma[, "alpha"]
beta <- Gamma[, 2] # or beta <- Gamma[, "beta"]
print(Gamma)## alpha beta
## AMZN 0.01697563 -0.2085693
## HAS 0.01356192 -0.2763406
## JNJ 0.01239375 1.2967969
## MA 0.02358281 1.2943668
## MNST 0.03401129 1.6321305
## MSFT 0.01124666 0.3790205
## NFLX 0.05205511 1.2534143
## VRTX 0.02456127 -0.5415507# compute Sigma
X <- X[insample_range]
f <- f$Mkt.RF[insample_range]
T <- dim(X)[1]
E <- xts(t(t(X) - Gamma %*% t(F_)), index(X)) # residuals
Psi <- (1/(T-2)) * t(E) %*% E
Sigma <- as.numeric(var(f)) * beta %o% beta + diag(diag(Psi))
Sigma # This is covariance matrix computed from single factor model## AMZN HAS JNJ MA MNST
## AMZN 6.539575e-03 6.295385e-06 -2.954265e-05 -2.948729e-05 -3.718197e-05
## HAS 6.295385e-06 4.333560e-03 -3.914207e-05 -3.906872e-05 -4.926367e-05
## JNJ -2.954265e-05 -3.914207e-05 1.513699e-03 1.833397e-04 2.311820e-04
## MA -2.948729e-05 -3.906872e-05 1.833397e-04 4.267586e-03 2.307487e-04
## MNST -3.718197e-05 -4.926367e-05 2.311820e-04 2.307487e-04 8.168392e-03
## MSFT -8.634561e-06 -1.144022e-05 5.368609e-05 5.358548e-05 6.756856e-05
## NFLX -2.855434e-05 -3.783263e-05 1.775390e-04 1.772063e-04 2.234481e-04
## VRTX 1.233720e-05 1.634598e-05 -7.670757e-05 -7.656382e-05 -9.654308e-05
## MSFT NFLX VRTX
## AMZN -8.634561e-06 -2.855434e-05 1.233720e-05
## HAS -1.144022e-05 -3.783263e-05 1.634598e-05
## JNJ 5.368609e-05 1.775390e-04 -7.670757e-05
## MA 5.358548e-05 1.772063e-04 -7.656382e-05
## MNST 6.756856e-05 2.234481e-04 -9.654308e-05
## MSFT 3.640926e-03 5.189009e-05 -2.241966e-05
## NFLX 5.189009e-05 4.530811e-02 -7.414142e-05
## VRTX -2.241966e-05 -7.414142e-05 1.777782e-02#=======================================================================
# We can also use lm() to compute estimated coefficients
fit = lm(formula = X~f)
sigF = as.numeric(var(f))
beta_ = as.matrix(fit$coefficients)
beta_ = as.matrix(beta_[-1,])
beta_## [,1]
## AMZN -0.2085693
## HAS -0.2763406
## JNJ 1.2967969
## MA 1.2943668
## MNST 1.6321305
## MSFT 0.3790205
## NFLX 1.2534143
## VRTX -0.5415507sigeps = crossprod(fit$residuals)/(T-2)
# sigeps = as.matrix(var(fit$residuals)) # you can use this way too
sigeps = diag(diag(sigeps))
sigeps## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.006534823 0.000000000 0.000000000 0.00000000 0.000000000 0.000000000
## [2,] 0.000000000 0.004325219 0.000000000 0.00000000 0.000000000 0.000000000
## [3,] 0.000000000 0.000000000 0.001330015 0.00000000 0.000000000 0.000000000
## [4,] 0.000000000 0.000000000 0.000000000 0.00408459 0.000000000 0.000000000
## [5,] 0.000000000 0.000000000 0.000000000 0.00000000 0.007877429 0.000000000
## [6,] 0.000000000 0.000000000 0.000000000 0.00000000 0.000000000 0.003625235
## [7,] 0.000000000 0.000000000 0.000000000 0.00000000 0.000000000 0.000000000
## [8,] 0.000000000 0.000000000 0.000000000 0.00000000 0.000000000 0.000000000
## [,7] [,8]
## [1,] 0.00000000 0.00000000
## [2,] 0.00000000 0.00000000
## [3,] 0.00000000 0.00000000
## [4,] 0.00000000 0.00000000
## [5,] 0.00000000 0.00000000
## [6,] 0.00000000 0.00000000
## [7,] 0.04513651 0.00000000
## [8,] 0.00000000 0.01774579cov_1f = sigF*beta_%*%t(beta_)+sigeps
cov_1f # This is covariance matrix computed from single factor model (from lm())## AMZN HAS JNJ MA MNST
## AMZN 6.539575e-03 6.295385e-06 -2.954265e-05 -2.948729e-05 -3.718197e-05
## HAS 6.295385e-06 4.333560e-03 -3.914207e-05 -3.906872e-05 -4.926367e-05
## JNJ -2.954265e-05 -3.914207e-05 1.513699e-03 1.833397e-04 2.311820e-04
## MA -2.948729e-05 -3.906872e-05 1.833397e-04 4.267586e-03 2.307487e-04
## MNST -3.718197e-05 -4.926367e-05 2.311820e-04 2.307487e-04 8.168392e-03
## MSFT -8.634561e-06 -1.144022e-05 5.368609e-05 5.358548e-05 6.756856e-05
## NFLX -2.855434e-05 -3.783263e-05 1.775390e-04 1.772063e-04 2.234481e-04
## VRTX 1.233720e-05 1.634598e-05 -7.670757e-05 -7.656382e-05 -9.654308e-05
## MSFT NFLX VRTX
## AMZN -8.634561e-06 -2.855434e-05 1.233720e-05
## HAS -1.144022e-05 -3.783263e-05 1.634598e-05
## JNJ 5.368609e-05 1.775390e-04 -7.670757e-05
## MA 5.358548e-05 1.772063e-04 -7.656382e-05
## MNST 6.756856e-05 2.234481e-04 -9.654308e-05
## MSFT 3.640926e-03 5.189009e-05 -2.241966e-05
## NFLX 5.189009e-05 4.530811e-02 -7.414142e-05
## VRTX -2.241966e-05 -7.414142e-05 1.777782e-02#=======================================================================
# Backtesting portfolio using SIT package
#=======================================================================
dates <- '2005-01::2021-03'
X <- monthly.ret[dates]
head(X)## AMZN HAS JNJ MA MNST
## 2006-06-30 0.11759604 -0.02319345 -0.004982313 0.06809068 0.02975055
## 2006-07-31 -0.30480871 0.03928422 0.043892354 -0.04437512 -0.03367123
## 2006-08-31 0.14652291 0.08556095 0.039731504 0.21866125 -0.40117417
## 2006-09-29 0.04184233 0.12069001 0.004330586 0.25849770 0.17937538
## 2006-10-31 0.18586554 0.14481394 0.037880769 0.05466816 -0.02247525
## 2006-11-30 0.05907062 0.03202161 -0.016583650 0.37314442 -0.11401587
## MSFT NFLX VRTX
## 2006-06-30 0.02869766 -0.01733467 0.06405794
## 2006-07-31 0.03261806 -0.23961789 -0.08689728
## 2006-08-31 0.07209598 -0.03238270 0.02774466
## 2006-09-29 0.06420299 0.13786224 -0.02322203
## 2006-10-31 0.04972508 0.21422303 0.20653776
## 2006-11-30 0.02613658 0.05856818 0.09113303dim(X)## [1] 178 8tail(X)## AMZN HAS JNJ MA MNST
## 2020-10-30 -0.0357541227 0.008173986 -0.07905702 -0.145470847 -0.04526181
## 2020-11-30 0.0434399293 0.124637388 0.06254355 0.165846635 0.10722214
## 2020-12-31 0.0280583237 0.005482124 0.08777995 0.060711424 0.09082336
## 2021-01-29 -0.0155760124 0.002993361 0.03653580 -0.112761444 -0.06109430
## 2021-02-26 -0.0353284326 0.006119783 -0.02257026 0.118751799 0.01048020
## 2021-03-31 0.0003718629 0.025717667 0.03717029 0.006217229 0.03818097
## MSFT NFLX VRTX
## 2020-10-30 -0.037369745 -0.04857710 -0.23430838
## 2020-11-30 0.060060517 0.03144579 0.09306008
## 2020-12-31 0.039005834 0.10195632 0.03771678
## 2021-01-29 0.042891914 -0.01542353 -0.03071843
## 2021-02-26 0.004117931 0.01213389 -0.07215819
## 2021-03-31 0.014588161 -0.03190128 0.01100916#
prices_monthly <- to.monthly(datas$prices, indexAt = "last", OHLC = FALSE) # indexAt = 'endof'
head(prices_monthly)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2006-05-31 34.61 12.23199 38.83738 4.190823 7.702917 16.50336 3.955714 34.50
## 2006-06-30 38.68 11.94828 38.64388 4.476179 7.932083 16.97697 3.887143 36.71
## 2006-07-31 26.89 12.41766 40.34005 4.277548 7.665000 17.53073 2.955714 33.52
## 2006-08-31 30.83 13.48013 41.94282 5.212882 4.590000 18.79462 2.860000 34.45
## 2006-09-29 32.12 15.10705 42.12446 6.560400 5.413333 20.00129 3.254286 33.65
## 2006-10-31 38.09 17.29476 43.72017 6.919045 5.291667 20.99586 3.951429 40.60prices_monthly <- prices_monthly[dates]
head(prices_monthly)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2006-05-31 34.61 12.23199 38.83738 4.190823 7.702917 16.50336 3.955714 34.50
## 2006-06-30 38.68 11.94828 38.64388 4.476179 7.932083 16.97697 3.887143 36.71
## 2006-07-31 26.89 12.41766 40.34005 4.277548 7.665000 17.53073 2.955714 33.52
## 2006-08-31 30.83 13.48013 41.94282 5.212882 4.590000 18.79462 2.860000 34.45
## 2006-09-29 32.12 15.10705 42.12446 6.560400 5.413333 20.00129 3.254286 33.65
## 2006-10-31 38.09 17.29476 43.72017 6.919045 5.291667 20.99586 3.951429 40.60tail(prices_monthly)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2020-10-30 3036.15 81.55690 134.4810 287.9328 76.57 201.0131 475.74 208.36
## 2020-11-30 3168.04 91.72194 142.8919 335.6855 84.78 213.0861 490.70 227.75
## 2020-12-31 3256.93 92.22477 155.4350 356.0654 92.48 221.3977 540.73 236.34
## 2021-01-29 3206.20 92.50083 161.1139 315.9150 86.83 230.8938 532.39 229.08
## 2021-02-26 3092.93 93.06692 157.4775 353.4305 87.74 231.8447 538.85 212.55
## 2021-03-31 3094.08 95.46038 163.3310 355.6278 91.09 235.2268 521.66 214.89#
#*****************************************************************
# Code Strategies
#******************************************************************
# We need to create environment variable required by SIT
# data_m <- new.env()
#
# data_m$EEM <- prices_monthly$EEM
# colnames(data_m$EEM) <- 'Close'
# head(data_m$EEM)
# Change column names to Close which will be consistent with the requirement of SIT package
# data_m$SPY <- prices_monthly$SPY %>% `colnames<-` (c("Close"))
# data_m$QQQ <- prices_monthly$QQQ %>% `colnames<-` (c("Close"))
# data_m$EEM <- prices_monthly$EEM %>% `colnames<-` (c("Close"))
# data_m$IWM <- prices_monthly$IWM %>% `colnames<-` (c("Close"))
# data_m$EFA <- prices_monthly$EFA %>% `colnames<-` (c("Close"))
# data_m$TLT <- prices_monthly$TLT %>% `colnames<-` (c("Close"))
# data_m$IYR <- prices_monthly$IYR %>% `colnames<-` (c("Close"))
# data_m$GLD <- prices_monthly$GLD %>% `colnames<-` (c("Close"))
# names(data_m)
# head(data_m$EEM)
# Using loop to save time in naming
data_m <- list()
#i = 1
for ( i in 1:length(tickers)){
data_m[[tickers[i]]] <- prices_monthly[,i] %>% `colnames<-` (c("Close"))
}
# data_env <- new.env()
# list2env(data_m, envir = data_env)
# convert from list to environment variable because SIT package requires
# the input data to be environmental variable which can be processed using bt.prep()
data_m <- list2env(data_m)
names(data_m)## [1] "MNST" "MSFT" "NFLX" "AMZN" "VRTX" "JNJ" "HAS" "MA"# Check if the column name is 'Close'
head(data_m$EEM)## NULL#
bt.prep(data_m, align='remove.na', dates = dates)
names(data_m)## [1] "prices" "execution.price" "weight" "dates"
## [5] "symbolnames" "MNST" "MSFT" "NFLX"
## [9] "AMZN" "VRTX" "JNJ" "HAS"
## [13] "MA"#===============================================================
# Equal weighting portfolio
# Equal Weight 1/N Benchmark
#===============================================================
models <- list()
prices <- data_m$prices
# data_m$weight[] = NA
N <- length(tickers)
data_m$weight = ntop(prices, N)
head(data_m$weight)## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2006-05-31 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## 2006-06-30 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## 2006-07-31 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## 2006-08-31 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## 2006-09-29 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## 2006-10-31 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125data_m$weight[1:59, ] <- NA
models$equal.weight = bt.run.share(data_m)## Latest weights :
## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2021-03-31 10.69 3.58 4.33 21.21 5.18 21.77 15.09 18.16
##
## Performance summary :
## CAGR Best Worst
## 16.2 16.5 -13.7 # head(models$equal.weight$ret, 62)
# Slightly difference between bt.run.share() and bt.run()
capital = 100000
data_m$weight[] = (capital / prices) * data_m$weight
models$equal.weight.share = bt.run(data_m, type='share')## Latest weights :
## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2021-03-31 12.5 12.5 12.5 12.5 12.5 12.5 12.5 12.5
##
## Performance summary :
## CAGR Best Worst
## 17.2 16.3 -11.2 # head(models$equal.weight.share$ret, 62)
# head(equal.weight$ret)
#================================================================
# MVP portfolio
#*****************************************************************
# 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(ret)
weight[] = NA
#
# head(data_m$weight)
# compute covariance matrix based on historical 60 months returns
# i = 60
# To make covariance matrix estimate more stable, use the Ledoit-Wolf covariance shrinkage estimator from tawny package
# 1. ia$cov = tawny::cov.shrink(hist) or
# 2. ia$cov = cor(coredata(hist), use='complete.obs', method='spearman') * (s0 %*% t(s0)) or
# 3. ia$cov = cor(coredata(hist), use='complete.obs', method='kendall') * (s0 %*% t(s0)
# i = 60
for( i in 60:dim(weight)[1]) {
hist = ret[ (i- 60 + 1):i, ]
# create historical input assumptions
ia = create.historical.ia(hist, 12) # 12 is annulized factor for monthly returns
s0 = apply(na.omit(coredata(hist)), 2, sd)
ia$cov = cor(coredata(hist), use='complete.obs',method='pearson') * (s0 %*% t(s0))
weight[i,] = min.risk.portfolio(ia, constraints) # use SIT's function min.risk.portfolio()
}## Loading required package: kernlab##
## Attaching package: 'kernlab'## The following object is masked from 'package:SIT':
##
## cross## The following object is masked from 'package:purrr':
##
## crossdim(weight)## [1] 179 8# 195 8
head(weight, 70)## AMZN HAS JNJ MA MNST MSFT
## [1,] NA NA NA NA NA NA
## [2,] NA NA NA NA NA NA
## [3,] NA NA NA NA NA NA
## [4,] NA NA NA NA NA NA
## [5,] NA NA NA NA NA NA
## [6,] NA NA NA NA NA NA
## [7,] NA NA NA NA NA NA
## [8,] NA NA NA NA NA NA
## [9,] NA NA NA NA NA NA
## [10,] NA NA NA NA NA NA
## [11,] NA NA NA NA NA NA
## [12,] NA NA NA NA NA NA
## [13,] NA NA NA NA NA NA
## [14,] NA NA NA NA NA NA
## [15,] NA NA NA NA NA NA
## [16,] NA NA NA NA NA NA
## [17,] NA NA NA NA NA NA
## [18,] NA NA NA NA NA NA
## [19,] NA NA NA NA NA NA
## [20,] NA NA NA NA NA NA
## [21,] NA NA NA NA NA NA
## [22,] NA NA NA NA NA NA
## [23,] NA NA NA NA NA NA
## [24,] NA NA NA NA NA NA
## [25,] NA NA NA NA NA NA
## [26,] NA NA NA NA NA NA
## [27,] NA NA NA NA NA NA
## [28,] NA NA NA NA NA NA
## [29,] NA NA NA NA NA NA
## [30,] NA NA NA NA NA NA
## [31,] NA NA NA NA NA NA
## [32,] NA NA NA NA NA NA
## [33,] NA NA NA NA NA NA
## [34,] NA NA NA NA NA NA
## [35,] NA NA NA NA NA NA
## [36,] NA NA NA NA NA NA
## [37,] NA NA NA NA NA NA
## [38,] NA NA NA NA NA NA
## [39,] NA NA NA NA NA NA
## [40,] NA NA NA NA NA NA
## [41,] NA NA NA NA NA NA
## [42,] NA NA NA NA NA NA
## [43,] NA NA NA NA NA NA
## [44,] NA NA NA NA NA NA
## [45,] NA NA NA NA NA NA
## [46,] NA NA NA NA NA NA
## [47,] NA NA NA NA NA NA
## [48,] NA NA NA NA NA NA
## [49,] NA NA NA NA NA NA
## [50,] NA NA NA NA NA NA
## [51,] NA NA NA NA NA NA
## [52,] NA NA NA NA NA NA
## [53,] NA NA NA NA NA NA
## [54,] NA NA NA NA NA NA
## [55,] NA NA NA NA NA NA
## [56,] NA NA NA NA NA NA
## [57,] NA NA NA NA NA NA
## [58,] NA NA NA NA NA NA
## [59,] NA NA NA NA NA NA
## [60,] 0.07782913 0.03137722 0.7903689 -0.05434008 0.009173707 0.06541928
## [61,] 0.07639972 0.04071648 0.7785839 -0.05206156 0.010367126 0.06178763
## [62,] 0.07712825 0.04118901 0.7793875 -0.05308010 0.009052984 0.06201955
## [63,] 0.07155715 0.03160720 0.7952658 -0.06284437 0.005566458 0.05733215
## [64,] 0.06945856 0.03392369 0.7881626 -0.06246487 0.007740626 0.06024423
## [65,] 0.08870090 0.02079686 0.8024643 -0.05434820 0.011109044 0.03091883
## [66,] 0.08527225 0.02606438 0.7903214 -0.05087071 0.016958742 0.03537449
## [67,] 0.08809357 0.02847967 0.7938247 -0.04857773 0.014312228 0.03470435
## [68,] 0.09177453 0.03527649 0.7917262 -0.04481854 0.009960164 0.03226565
## [69,] 0.09746619 0.03301324 0.8090799 -0.04601894 0.009831585 0.01639175
## [70,] 0.09970935 0.02194186 0.8092013 -0.05133660 0.016537753 0.01773513
## NFLX VRTX
## [1,] NA NA
## [2,] NA NA
## [3,] NA NA
## [4,] NA NA
## [5,] NA NA
## [6,] NA NA
## [7,] NA NA
## [8,] NA NA
## [9,] NA NA
## [10,] NA NA
## [11,] NA NA
## [12,] NA NA
## [13,] NA NA
## [14,] NA NA
## [15,] NA NA
## [16,] NA NA
## [17,] NA NA
## [18,] NA NA
## [19,] NA NA
## [20,] NA NA
## [21,] NA NA
## [22,] NA NA
## [23,] NA NA
## [24,] NA NA
## [25,] NA NA
## [26,] NA NA
## [27,] NA NA
## [28,] NA NA
## [29,] NA NA
## [30,] NA NA
## [31,] NA NA
## [32,] NA NA
## [33,] NA NA
## [34,] NA NA
## [35,] NA NA
## [36,] NA NA
## [37,] NA NA
## [38,] NA NA
## [39,] NA NA
## [40,] NA NA
## [41,] NA NA
## [42,] NA NA
## [43,] NA NA
## [44,] NA NA
## [45,] NA NA
## [46,] NA NA
## [47,] NA NA
## [48,] NA NA
## [49,] NA NA
## [50,] NA NA
## [51,] NA NA
## [52,] NA NA
## [53,] NA NA
## [54,] NA NA
## [55,] NA NA
## [56,] NA NA
## [57,] NA NA
## [58,] NA NA
## [59,] NA NA
## [60,] 0.04248265 0.03768918
## [61,] 0.04918436 0.03502229
## [62,] 0.04971010 0.03459274
## [63,] 0.07105093 0.03046468
## [64,] 0.07436125 0.02857395
## [65,] 0.06079157 0.03956668
## [66,] 0.05952938 0.03735004
## [67,] 0.06092575 0.02823740
## [68,] 0.05948049 0.02433502
## [69,] 0.04406403 0.03617223
## [70,] 0.04879693 0.03741424tail(weight)## AMZN HAS JNJ MA MNST MSFT
## [174,] -0.04746791 0.1341246 0.5971092 -0.1829207 0.07514950 0.09825748
## [175,] -0.05060838 0.1292616 0.6116818 -0.1712008 0.07331142 0.10873543
## [176,] -0.05536037 0.1436823 0.5768308 -0.1610075 0.08083710 0.09914714
## [177,] -0.05061918 0.1381929 0.5615687 -0.1553686 0.08266347 0.09549654
## [178,] -0.04331664 0.1321348 0.5386311 -0.1591294 0.07824848 0.08660776
## [179,] -0.04016588 0.1311905 0.5316657 -0.1588674 0.06826865 0.08785905
## NFLX VRTX
## [174,] 0.3383856 -0.012637731
## [175,] 0.3612724 -0.062453482
## [176,] 0.3807184 -0.064847947
## [177,] 0.3594550 -0.031388810
## [178,] 0.3737855 -0.006961676
## [179,] 0.3822735 -0.002224169# Assign minimum variance weights to data_m$weight
capital = 100000
data_m$weight <- data_m$prices
data_m$weight[] <- NA
data_m$weight[] <- weight
data_m$weight[] = (capital / prices) * data_m$weight
models$mvp.hist.cov = bt.run(data_m, type='share')## Latest weights :
## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2021-03-31 -4.33 13.21 53.86 -15.91 7.82 8.66 37.38 -0.7
##
## Performance summary :
## CAGR Best Worst
## 11 11.5 -10.2 #=================================================================
# CAPM model (single factor model)
#*****************************************************************
# 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_1 = coredata(ret)
weight_1[] = NA
#
# head(data_m$weight)
# compute covariance matrix based on historical 60 months returns
# i = 60
f <- fama_lib_month[dates]/100
for( i in 60:dim(weight_1)[1]) {
hist = ret[ (i- 60 + 1):i, ]
Xi <- hist
fi <- f$Mkt.RF[(i - 60 + 1):i, ]
fiti = lm(formula = Xi ~ fi)
sigF = as.numeric(var(fi))
beta_ = as.matrix(fiti$coefficients)
beta_ = as.matrix(beta_[-1,])
sigeps = crossprod(fiti$residuals)/(T-2)
# sigeps = as.matrix(var(fit$residuals)) # you can use this way too
sigeps = diag(diag(sigeps))
cov_1f = sigF*beta_%*%t(beta_)+sigeps
# cov_1f
ia$cov = cov_1f
weight_1[i, ] = min.risk.portfolio(ia, constraints)
}
dim(weight_1)## [1] 179 8# 195 8
head(weight_1, 70)## AMZN HAS JNJ MA MNST MSFT NFLX
## [1,] NA NA NA NA NA NA NA
## [2,] NA NA NA NA NA NA NA
## [3,] NA NA NA NA NA NA NA
## [4,] NA NA NA NA NA NA NA
## [5,] NA NA NA NA NA NA NA
## [6,] NA NA NA NA NA NA NA
## [7,] NA NA NA NA NA NA NA
## [8,] NA NA NA NA NA NA NA
## [9,] NA NA NA NA NA NA NA
## [10,] NA NA NA NA NA NA NA
## [11,] NA NA NA NA NA NA NA
## [12,] NA NA NA NA NA NA NA
## [13,] NA NA NA NA NA NA NA
## [14,] NA NA NA NA NA NA NA
## [15,] NA NA NA NA NA NA NA
## [16,] NA NA NA NA NA NA NA
## [17,] NA NA NA NA NA NA NA
## [18,] NA NA NA NA NA NA NA
## [19,] NA NA NA NA NA NA NA
## [20,] NA NA NA NA NA NA NA
## [21,] NA NA NA NA NA NA NA
## [22,] NA NA NA NA NA NA NA
## [23,] NA NA NA NA NA NA NA
## [24,] NA NA NA NA NA NA NA
## [25,] NA NA NA NA NA NA NA
## [26,] NA NA NA NA NA NA NA
## [27,] NA NA NA NA NA NA NA
## [28,] NA NA NA NA NA NA NA
## [29,] NA NA NA NA NA NA NA
## [30,] NA NA NA NA NA NA NA
## [31,] NA NA NA NA NA NA NA
## [32,] NA NA NA NA NA NA NA
## [33,] NA NA NA NA NA NA NA
## [34,] NA NA NA NA NA NA NA
## [35,] NA NA NA NA NA NA NA
## [36,] NA NA NA NA NA NA NA
## [37,] NA NA NA NA NA NA NA
## [38,] NA NA NA NA NA NA NA
## [39,] NA NA NA NA NA NA NA
## [40,] NA NA NA NA NA NA NA
## [41,] NA NA NA NA NA NA NA
## [42,] NA NA NA NA NA NA NA
## [43,] NA NA NA NA NA NA NA
## [44,] NA NA NA NA NA NA NA
## [45,] NA NA NA NA NA NA NA
## [46,] NA NA NA NA NA NA NA
## [47,] NA NA NA NA NA NA NA
## [48,] NA NA NA NA NA NA NA
## [49,] NA NA NA NA NA NA NA
## [50,] NA NA NA NA NA NA NA
## [51,] NA NA NA NA NA NA NA
## [52,] NA NA NA NA NA NA NA
## [53,] NA NA NA NA NA NA NA
## [54,] NA NA NA NA NA NA NA
## [55,] NA NA NA NA NA NA NA
## [56,] NA NA NA NA NA NA NA
## [57,] NA NA NA NA NA NA NA
## [58,] NA NA NA NA NA NA NA
## [59,] NA NA NA NA NA NA NA
## [60,] 0.05686594 0.1510929 0.4305811 0.05015941 0.05779332 0.05014609 0.1448397
## [61,] 0.05687191 0.1511282 0.4302659 0.05041776 0.05774930 0.04963676 0.1451084
## [62,] 0.05622834 0.1506700 0.4306342 0.05089247 0.05776042 0.04961835 0.1451982
## [63,] 0.05561291 0.1444594 0.4294557 0.05680025 0.05829546 0.05308205 0.1433751
## [64,] 0.06680349 0.1429472 0.4252998 0.05601470 0.05542776 0.05160215 0.1419638
## [65,] 0.06976960 0.1382716 0.4298696 0.05647285 0.05789554 0.03914050 0.1440393
## [66,] 0.06989978 0.1347698 0.4339531 0.05757943 0.05895068 0.03731231 0.1431713
## [67,] 0.07146909 0.1325493 0.4331714 0.05616649 0.05471233 0.03541899 0.1421959
## [68,] 0.07486064 0.1281545 0.4338631 0.05529816 0.05467167 0.03598601 0.1425410
## [69,] 0.07626443 0.1282828 0.4409171 0.05656702 0.05469353 0.02792388 0.1386400
## [70,] 0.07599612 0.1268518 0.4487767 0.05482404 0.05560280 0.02631219 0.1379022
## VRTX
## [1,] NA
## [2,] NA
## [3,] NA
## [4,] NA
## [5,] NA
## [6,] NA
## [7,] NA
## [8,] NA
## [9,] NA
## [10,] NA
## [11,] NA
## [12,] NA
## [13,] NA
## [14,] NA
## [15,] NA
## [16,] NA
## [17,] NA
## [18,] NA
## [19,] NA
## [20,] NA
## [21,] NA
## [22,] NA
## [23,] NA
## [24,] NA
## [25,] NA
## [26,] NA
## [27,] NA
## [28,] NA
## [29,] NA
## [30,] NA
## [31,] NA
## [32,] NA
## [33,] NA
## [34,] NA
## [35,] NA
## [36,] NA
## [37,] NA
## [38,] NA
## [39,] NA
## [40,] NA
## [41,] NA
## [42,] NA
## [43,] NA
## [44,] NA
## [45,] NA
## [46,] NA
## [47,] NA
## [48,] NA
## [49,] NA
## [50,] NA
## [51,] NA
## [52,] NA
## [53,] NA
## [54,] NA
## [55,] NA
## [56,] NA
## [57,] NA
## [58,] NA
## [59,] NA
## [60,] 0.05852149
## [61,] 0.05882179
## [62,] 0.05899806
## [63,] 0.05891918
## [64,] 0.05994116
## [65,] 0.06454099
## [66,] 0.06436360
## [67,] 0.07431648
## [68,] 0.07462496
## [69,] 0.07671122
## [70,] 0.07373423tail(weight_1)## AMZN HAS JNJ MA MNST MSFT
## [174,] 0.09246307 0.09035330 0.2684925 0.08387234 0.06446031 0.04891578
## [175,] 0.09538089 0.08969731 0.2665227 0.08676005 0.06446391 0.05262093
## [176,] 0.09536256 0.09276303 0.2578677 0.08738041 0.06541591 0.05332276
## [177,] 0.09497650 0.09200393 0.2497972 0.09164210 0.07597399 0.05751137
## [178,] 0.09554648 0.09310171 0.2472471 0.09027778 0.07543603 0.05445638
## [179,] 0.09512519 0.09314546 0.2451374 0.08941992 0.07648533 0.05360068
## NFLX VRTX
## [174,] 0.2076641 0.1437786
## [175,] 0.2096846 0.1348697
## [176,] 0.2117629 0.1361247
## [177,] 0.2085731 0.1295218
## [178,] 0.2181206 0.1258139
## [179,] 0.2207411 0.1263448# Assign minimum variance weights to data_m$weight
capital = 100000
data_m$weight <- data_m$prices
data_m$weight[] <- NA
data_m$weight[] <- weight_1
data_m$weight[] = (capital / prices) * data_m$weight
models$mvp.capm.cov = bt.run(data_m, type='share')## Latest weights :
## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2021-03-31 9.55 9.31 24.72 9.03 7.54 5.45 21.81 12.58
##
## Performance summary :
## CAGR Best Worst
## 13.9 13.2 -10.4 models$mvp.capm.cov$cagr## [1] 0.139108models$mvp.hist.cov$cagr## [1] 0.1096295#=================================================================
# FF 3 factor model
#*****************************************************************
# 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_3 = coredata(ret)
weight_3[] = NA
#
# head(data_m$weight)
# compute covariance matrix based on historical 60 months returns
# i = 60
f <- fama_lib_month[dates]/100
for( i in 60:dim(weight_1)[1]) {
hist = ret[ (i- 60 + 1):i, ]
Xi <- hist
fi <- f[(i - 60 + 1):i, ]
fiti = lm(formula = Xi ~ fi)
sigF = as.matrix(var(fi))
beta_ = as.matrix(fiti$coefficients)
beta_ = as.matrix(beta_[-1,])
sigeps = crossprod(fiti$residuals)/(T-4) # note (T - 4)
# sigeps = as.matrix(var(fit$residuals)) # you can use this way too
sigeps = diag(diag(sigeps))
cov_3f = t(beta_)%*% sigF %*% beta_ + sigeps
# cov_1f
ia$cov = cov_3f
weight_3[i, ] = min.risk.portfolio(ia, constraints)
}
dim(weight_3)## [1] 179 8# 195 8
head(weight_3, 70)## AMZN HAS JNJ MA MNST MSFT NFLX
## [1,] NA NA NA NA NA NA NA
## [2,] NA NA NA NA NA NA NA
## [3,] NA NA NA NA NA NA NA
## [4,] NA NA NA NA NA NA NA
## [5,] NA NA NA NA NA NA NA
## [6,] NA NA NA NA NA NA NA
## [7,] NA NA NA NA NA NA NA
## [8,] NA NA NA NA NA NA NA
## [9,] NA NA NA NA NA NA NA
## [10,] NA NA NA NA NA NA NA
## [11,] NA NA NA NA NA NA NA
## [12,] NA NA NA NA NA NA NA
## [13,] NA NA NA NA NA NA NA
## [14,] NA NA NA NA NA NA NA
## [15,] NA NA NA NA NA NA NA
## [16,] NA NA NA NA NA NA NA
## [17,] NA NA NA NA NA NA NA
## [18,] NA NA NA NA NA NA NA
## [19,] NA NA NA NA NA NA NA
## [20,] NA NA NA NA NA NA NA
## [21,] NA NA NA NA NA NA NA
## [22,] NA NA NA NA NA NA NA
## [23,] NA NA NA NA NA NA NA
## [24,] NA NA NA NA NA NA NA
## [25,] NA NA NA NA NA NA NA
## [26,] NA NA NA NA NA NA NA
## [27,] NA NA NA NA NA NA NA
## [28,] NA NA NA NA NA NA NA
## [29,] NA NA NA NA NA NA NA
## [30,] NA NA NA NA NA NA NA
## [31,] NA NA NA NA NA NA NA
## [32,] NA NA NA NA NA NA NA
## [33,] NA NA NA NA NA NA NA
## [34,] NA NA NA NA NA NA NA
## [35,] NA NA NA NA NA NA NA
## [36,] NA NA NA NA NA NA NA
## [37,] NA NA NA NA NA NA NA
## [38,] NA NA NA NA NA NA NA
## [39,] NA NA NA NA NA NA NA
## [40,] NA NA NA NA NA NA NA
## [41,] NA NA NA NA NA NA NA
## [42,] NA NA NA NA NA NA NA
## [43,] NA NA NA NA NA NA NA
## [44,] NA NA NA NA NA NA NA
## [45,] NA NA NA NA NA NA NA
## [46,] NA NA NA NA NA NA NA
## [47,] NA NA NA NA NA NA NA
## [48,] NA NA NA NA NA NA NA
## [49,] NA NA NA NA NA NA NA
## [50,] NA NA NA NA NA NA NA
## [51,] NA NA NA NA NA NA NA
## [52,] NA NA NA NA NA NA NA
## [53,] NA NA NA NA NA NA NA
## [54,] NA NA NA NA NA NA NA
## [55,] NA NA NA NA NA NA NA
## [56,] NA NA NA NA NA NA NA
## [57,] NA NA NA NA NA NA NA
## [58,] NA NA NA NA NA NA NA
## [59,] NA NA NA NA NA NA NA
## [60,] 0.04842715 0.1357863 0.4592521 0.05613880 0.05975124 0.05517867 0.1268433
## [61,] 0.04895882 0.1356565 0.4587749 0.05610711 0.05925775 0.05534075 0.1270160
## [62,] 0.04871341 0.1345362 0.4587999 0.05709223 0.05929200 0.05485767 0.1274504
## [63,] 0.04755795 0.1284438 0.4619051 0.06154855 0.05913337 0.05658082 0.1261895
## [64,] 0.06834662 0.1237107 0.4621842 0.05822498 0.05291668 0.05525798 0.1214689
## [65,] 0.07395036 0.1197986 0.4646616 0.05825459 0.05424676 0.04042605 0.1238027
## [66,] 0.07252312 0.1186421 0.4644371 0.06020387 0.05651412 0.03914572 0.1242230
## [67,] 0.07304868 0.1160514 0.4643233 0.05791766 0.05106528 0.03626438 0.1228793
## [68,] 0.07534593 0.1138328 0.4620345 0.05847413 0.05028882 0.03699772 0.1237610
## [69,] 0.07552435 0.1135430 0.4730067 0.06002339 0.05039039 0.02764751 0.1182869
## [70,] 0.07577531 0.1104433 0.4802681 0.05871723 0.05124038 0.02530697 0.1191462
## VRTX
## [1,] NA
## [2,] NA
## [3,] NA
## [4,] NA
## [5,] NA
## [6,] NA
## [7,] NA
## [8,] NA
## [9,] NA
## [10,] NA
## [11,] NA
## [12,] NA
## [13,] NA
## [14,] NA
## [15,] NA
## [16,] NA
## [17,] NA
## [18,] NA
## [19,] NA
## [20,] NA
## [21,] NA
## [22,] NA
## [23,] NA
## [24,] NA
## [25,] NA
## [26,] NA
## [27,] NA
## [28,] NA
## [29,] NA
## [30,] NA
## [31,] NA
## [32,] NA
## [33,] NA
## [34,] NA
## [35,] NA
## [36,] NA
## [37,] NA
## [38,] NA
## [39,] NA
## [40,] NA
## [41,] NA
## [42,] NA
## [43,] NA
## [44,] NA
## [45,] NA
## [46,] NA
## [47,] NA
## [48,] NA
## [49,] NA
## [50,] NA
## [51,] NA
## [52,] NA
## [53,] NA
## [54,] NA
## [55,] NA
## [56,] NA
## [57,] NA
## [58,] NA
## [59,] NA
## [60,] 0.05862239
## [61,] 0.05888821
## [62,] 0.05925820
## [63,] 0.05864089
## [64,] 0.05788988
## [65,] 0.06485939
## [66,] 0.06431092
## [67,] 0.07844995
## [68,] 0.07926503
## [69,] 0.08157785
## [70,] 0.07910255tail(weight_3)## AMZN HAS JNJ MA MNST MSFT
## [174,] 0.09297260 0.08611206 0.2738790 0.08635669 0.06048155 0.04969762
## [175,] 0.09562278 0.08720020 0.2703185 0.08898391 0.06123324 0.05296350
## [176,] 0.09584257 0.09056380 0.2611364 0.08899016 0.06280090 0.05381083
## [177,] 0.09513075 0.08991215 0.2540937 0.09206149 0.07353568 0.05740522
## [178,] 0.09576737 0.09184274 0.2494079 0.09133644 0.07389897 0.05456011
## [179,] 0.09497598 0.09300893 0.2446081 0.09083962 0.07683494 0.05410113
## NFLX VRTX
## [174,] 0.2086505 0.1418501
## [175,] 0.2097610 0.1339168
## [176,] 0.2115070 0.1353484
## [177,] 0.2083110 0.1295499
## [178,] 0.2180522 0.1251343
## [179,] 0.2201700 0.1254614# Assign minimum variance weights to data_m$weight
capital = 100000
data_m$weight <- data_m$prices
data_m$weight[] <- NA
data_m$weight[] <- weight_3
data_m$weight[] = (capital / prices) * data_m$weight
models$mvp.ff3.cov = bt.run(data_m, type='share')## Latest weights :
## AMZN HAS JNJ MA MNST MSFT NFLX VRTX
## 2021-03-31 9.58 9.18 24.94 9.13 7.39 5.46 21.81 12.51
##
## Performance summary :
## CAGR Best Worst
## 13.8 13.2 -10.4 #============================================================
# Principal Component Analysis (PCA)
#=============================================================
library(covFactorModel)
X_i <-X[1:60,]
factor_pca <- factorModel(X_i, type = "S", K = 3, max_iter = 10)
cbind(alpha = factor_pca$alpha, beta = factor_pca$beta)## alpha factor1 factor2 factor3
## AMZN 0.037739835 -0.08360209 0.009611094 -0.009244520
## HAS 0.020044733 -0.03915069 0.013163221 0.005104457
## JNJ 0.005388372 -0.02306895 0.012034167 -0.008177582
## MA 0.038714604 -0.07436711 0.037833409 0.062646553
## MNST 0.015800446 -0.01908231 -0.030757167 0.009161600
## MSFT 0.006234520 -0.04432063 0.019191767 0.010594500
## NFLX 0.047374756 -0.07180474 -0.089119226 -0.004363538
## VRTX 0.015123667 -0.05337934 0.032324795 -0.079211047#
# Statistical 3-factor model
K <- 3
X_trn <- X_i
T_trn <- dim(X_i)[1]
alpha <- colMeans(X_trn)
X_trn_ <- X_trn - matrix(alpha, T_trn, N, byrow = TRUE)
Sigma_prev <- matrix(0, N, N)
Sigma <- (1/(T_trn-1)) * t(X_trn_) %*% X_trn_
eigSigma <- eigen(Sigma)
while (norm(Sigma - Sigma_prev, "F")/norm(Sigma, "F") > 1e-3) {
B <- eigSigma$vectors[, 1:K] %*% diag(sqrt(eigSigma$values[1:K]), K, K)
Psi <- diag(diag(Sigma - B %*% t(B)))
Sigma_prev <- Sigma
Sigma <- B %*% t(B) + Psi
eigSigma <- eigen(Sigma - Psi)
}
# B: factor loadings
B## [,1] [,2] [,3]
## [1,] 0.10929581 -0.030483746 0.003133449
## [2,] 0.03879183 -0.009906642 -0.012893169
## [3,] 0.02119561 -0.010404282 0.002155805
## [4,] 0.07268224 -0.023307620 -0.084581828
## [5,] 0.02962411 0.090106764 -0.006842262
## [6,] 0.04416130 -0.008559109 -0.021442576
## [7,] 0.07897783 0.081583122 0.033852929
## [8,] 0.05839933 -0.053151139 0.081089571# fiti = lm(formula = X_i ~ t(B))
Sigma_PCA3 <- Sigma
Sigma_PCA3## [,1] [,2] [,3] [,4] [,5]
## [1,] 0.017603623 0.0045013757 0.0026405076 0.0083893351 0.0004695590
## [2,] 0.004501376 0.0058612250 0.0008974926 0.0041409048 0.0003447362
## [3,] 0.002640508 0.0008974926 0.0020597991 0.0016007012 -0.0003243459
## [4,] 0.008389335 0.0041409048 0.0016007012 0.0150401900 0.0006317032
## [5,] 0.000469559 0.0003447362 -0.0003243459 0.0006317032 0.0155581163
## [6,] 0.005020369 0.0020743522 0.0009788509 0.0052228867 0.0006837210
## [7,] 0.006251063 0.0018190079 0.0008981495 0.0009754344 0.0094592080
## [8,] 0.008257138 0.0017464646 0.0019656221 -0.0013752833 -0.0036140853
## [,6] [,7] [,8]
## [1,] 0.0050203694 0.0062510631 0.008257138
## [2,] 0.0020743522 0.0018190079 0.001746465
## [3,] 0.0009788509 0.0008981495 0.001965622
## [4,] 0.0052228867 0.0009754344 -0.001375283
## [5,] 0.0006837210 0.0094592080 -0.003614085
## [6,] 0.0060980190 0.0020635905 0.001295147
## [7,] 0.0020635905 0.0183744713 0.003021136
## [8,] 0.0012951474 0.0030211361 0.015900071diag(Sigma_PCA3)## [1] 0.017603623 0.005861225 0.002059799 0.015040190 0.015558116 0.006098019
## [7] 0.018374471 0.015900071# error
norm(Sigma_PCA3 - cov(X_i), "F")## [1] 0.01159327#---------------------------------------------------------------
# By Eric Zivot
# use R princomp() function for principal component analysis
#---------------------------------------------------------------
pc.fit = princomp(X_i)
class(pc.fit)## [1] "princomp"names(pc.fit)## [1] "sdev" "loadings" "center" "scale" "n.obs" "scores" "call"pc.fit## Call:
## princomp(x = X_i)
##
## Standard deviations:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 0.17643836 0.13794816 0.12370749 0.11255611 0.08953133 0.06534236 0.05699564
## Comp.8
## 0.03493793
##
## 8 variables and 60 observations.summary(pc.fit)## Importance of components:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5
## Standard deviation 0.1764384 0.1379482 0.1237075 0.1125561 0.08953133
## Proportion of Variance 0.3280788 0.2005506 0.1612813 0.1335151 0.08447772
## Cumulative Proportion 0.3280788 0.5286293 0.6899107 0.8234257 0.90790347
## Comp.6 Comp.7 Comp.8
## Standard deviation 0.06534236 0.05699564 0.03493793
## Proportion of Variance 0.04499682 0.03423540 0.01286431
## Cumulative Proportion 0.95290029 0.98713569 1.00000000#eigenvalues & eigenvectors
eig.value <- eigen(cov(X_i))$values
pc.fit## Call:
## princomp(x = X_i)
##
## Standard deviations:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 0.17643836 0.13794816 0.12370749 0.11255611 0.08953133 0.06534236 0.05699564
## Comp.8
## 0.03493793
##
## 8 variables and 60 observations.eig.vec <- eigen(cov(X_i))$vectors
#
B_hat <- eig.vec[, 1:3]%*%diag(sqrt(eig.value[1:3]), 3, 3)
B_hat## [,1] [,2] [,3]
## [1,] -0.10929581 -0.030483746 0.003133449
## [2,] -0.03879183 -0.009906642 -0.012893169
## [3,] -0.02119561 -0.010404282 0.002155805
## [4,] -0.07268224 -0.023307620 -0.084581828
## [5,] -0.02962411 0.090106764 -0.006842262
## [6,] -0.04416130 -0.008559109 -0.021442576
## [7,] -0.07897783 0.081583122 0.033852929
## [8,] -0.05839933 -0.053151139 0.081089571#
#
plot(pc.fit)
loadings(pc.fit)##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 Comp.8
## AMZN 0.614 0.219 0.186 0.718 0.147
## HAS 0.218 -0.103 -0.101 -0.288 0.764 0.454 0.231
## JNJ 0.119 0.256 -0.952
## MA 0.408 0.168 -0.678 -0.389 -0.423 0.109
## MNST 0.166 -0.648 -0.685 0.201 -0.105 0.171
## MSFT 0.248 -0.172 -0.248 0.314 -0.844 0.179
## NFLX 0.444 -0.586 0.271 0.518 -0.318 -0.125
## VRTX 0.328 0.382 0.650 -0.387 -0.329 -0.249
##
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 Comp.8
## SS loadings 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000
## Proportion Var 0.125 0.125 0.125 0.125 0.125 0.125 0.125 0.125
## Cumulative Var 0.125 0.250 0.375 0.500 0.625 0.750 0.875 1.000pc.fit$loadings[, 1:3]## Comp.1 Comp.2 Comp.3
## AMZN 0.6142723 0.21913049 0.02511753
## HAS 0.2180206 0.07121327 -0.10335085
## JNJ 0.1191251 0.07479053 0.01728080
## MA 0.4084940 0.16754536 -0.67800274
## MNST 0.1664955 -0.64772681 -0.05484715
## MSFT 0.2481985 0.06152662 -0.17188237
## NFLX 0.4438769 -0.58645514 0.27136300
## VRTX 0.3282202 0.38207362 0.65000903# pc factors are in the scores component. Note these scores are based on
# centered data
head(pc.fit$scores[, 1:3])## Comp.1 Comp.2 Comp.3
## 2006-06-30 0.04562732 0.06755789 -0.003999946
## 2006-07-31 -0.39814153 0.07826193 -0.099607154
## 2006-08-31 0.07436684 0.38695527 -0.127339653
## 2006-09-29 0.18331852 -0.12529221 -0.178636022
## 2006-10-31 0.26988699 0.04920042 0.144891876
## 2006-11-30 0.16295140 0.16770275 -0.171682568pc.fit$scores[, 1:3]## Comp.1 Comp.2 Comp.3
## 2006-06-30 0.045627319 0.067557891 -0.003999946
## 2006-07-31 -0.398141527 0.078261931 -0.099607154
## 2006-08-31 0.074366836 0.386955268 -0.127339653
## 2006-09-29 0.183318520 -0.125292209 -0.178636022
## 2006-10-31 0.269886987 0.049200420 0.144891876
## 2006-11-30 0.162951400 0.167702749 -0.171682568
## 2006-12-29 -0.162138197 -0.111009282 -0.120590843
## 2007-01-31 -0.076315192 -0.002810830 -0.170678204
## 2007-02-28 -0.156530036 0.017217457 -0.036491478
## 2007-03-30 -0.074197071 -0.088977909 -0.039977776
## 2007-04-30 0.343624393 0.217157013 0.011602793
## 2007-05-31 0.144932574 0.071897250 -0.250327170
## 2007-06-29 -0.105946135 0.023711402 -0.124261516
## 2007-07-31 -0.040207805 0.184896970 0.102077697
## 2007-08-31 -0.027429284 -0.004199007 0.240790002
## 2007-09-28 0.192467578 -0.212801545 -0.019237009
## 2007-10-31 0.197120482 -0.277654356 -0.277360570
## 2007-11-30 -0.259237493 0.248284814 -0.163353995
## 2007-12-31 0.008392794 -0.105650031 -0.059462145
## 2008-01-31 -0.297615527 0.035029859 -0.049075185
## 2008-02-29 -0.175684605 -0.300246815 0.060324018
## 2008-03-31 0.234330487 0.258911642 0.148202495
## 2008-04-30 0.145732269 0.171926944 -0.165274974
## 2008-05-30 -0.004506148 0.195480772 0.003698828
## 2008-06-30 -0.228295855 0.161764204 0.180290271
## 2008-07-31 -0.012882649 0.055318047 0.146336388
## 2008-08-29 -0.067668224 -0.183335297 -0.150108908
## 2008-09-30 -0.172917390 -0.034587810 0.346460649
## 2008-10-31 -0.545658790 0.052010033 -0.026594878
## 2008-11-28 -0.303754659 -0.151846113 -0.034888591
## 2008-12-31 0.286744780 -0.106593362 0.248177861
## 2009-01-30 0.047196534 -0.071133737 0.194715057
## 2009-02-27 -0.012067246 0.017843675 -0.146576046
## 2009-03-31 0.173405677 -0.105240214 -0.048102436
## 2009-04-30 0.135510170 -0.026913534 -0.026343266
## 2009-05-29 -0.188282604 0.133609172 -0.020413133
## 2009-06-30 0.040861816 0.183838909 0.181908369
## 2009-07-31 0.064131809 0.020075204 -0.091778011
## 2009-08-31 -0.041305553 0.003991303 -0.021050290
## 2009-09-30 0.074614372 -0.060712606 0.024359619
## 2009-10-30 0.171040347 -0.034863371 -0.087438896
## 2009-11-30 0.187867393 0.103367572 0.053521044
## 2009-12-31 -0.001488515 0.044982248 -0.004794693
## 2010-01-29 -0.131965427 -0.129561398 0.005032901
## 2010-02-26 -0.047307341 -0.067037389 0.101388762
## 2010-03-31 0.156201898 -0.014748954 -0.055144022
## 2010-04-30 0.070836693 -0.215445677 0.072060797
## 2010-05-28 -0.239705767 -0.085952146 0.120471123
## 2010-06-30 -0.203925251 -0.025351981 -0.010763861
## 2010-07-30 0.025845258 0.045765327 -0.056335333
## 2010-08-31 0.016147461 -0.174176544 0.114591645
## 2010-09-30 0.328572932 -0.069480103 0.010653846
## 2010-10-29 0.105672903 -0.014197225 0.023135044
## 2010-11-30 -0.004752261 -0.162264530 -0.018229748
## 2010-12-31 -0.105135444 0.137354616 0.025684694
## 2011-01-31 0.041614166 -0.135184918 0.100189466
## 2011-02-28 -0.001441982 0.111069651 0.120122107
## 2011-03-31 0.046510814 -0.078779681 0.034525275
## 2011-04-29 0.093966962 0.066147915 0.026152109
## 2011-05-31 0.017010353 -0.135281680 0.014553584# time series plot of principal component factors
# library(PerformanceAnalytics)
# chart.TimeSeries(pc.fit$scores[, 1, drop=FALSE], colorset="purple")
# compare with direct eigen-value analysis
# notice the sign change on the first set of loadings
eigen.fit = eigen(var(X_i))
names(eigen.fit)## [1] "values" "vectors"names(eigen.fit$values) = rownames(eigen.fit$vectors) = colnames(X_i)
cbind(pc.fit$loadings[,1:3], eigen.fit$vectors[, 1:3])## Comp.1 Comp.2 Comp.3
## AMZN 0.6142723 0.21913049 0.02511753 -0.6142723 -0.21913049 0.02511753
## HAS 0.2180206 0.07121327 -0.10335085 -0.2180206 -0.07121327 -0.10335085
## JNJ 0.1191251 0.07479053 0.01728080 -0.1191251 -0.07479053 0.01728080
## MA 0.4084940 0.16754536 -0.67800274 -0.4084940 -0.16754536 -0.67800274
## MNST 0.1664955 -0.64772681 -0.05484715 -0.1664955 0.64772681 -0.05484715
## MSFT 0.2481985 0.06152662 -0.17188237 -0.2481985 -0.06152662 -0.17188237
## NFLX 0.4438769 -0.58645514 0.27136300 -0.4438769 0.58645514 0.27136300
## VRTX 0.3282202 0.38207362 0.65000903 -0.3282202 -0.38207362 0.65000903# compute uncentered pc factors from eigenvectors and return data
pc.factors.uc = X_i %*% eigen.fit$vectors
colnames(pc.factors.uc) = paste(colnames(pc.fit$scores),".uc",sep="")
# compare centered and uncentered scores. Note sign change on first factor
# We can treat centered scores as unobservable factor values (F)
cbind(pc.fit$scores[,1,drop=F], -pc.factors.uc[, 1, drop=F])## Comp.1 Comp.1.uc
## 2006-06-30 0.045627319 1.198071e-01
## 2006-07-31 -0.398141527 -3.239617e-01
## 2006-08-31 0.074366836 1.485467e-01
## 2006-09-29 0.183318520 2.574983e-01
## 2006-10-31 0.269886987 3.440668e-01
## 2006-11-30 0.162951400 2.371312e-01
## 2006-12-29 -0.162138197 -8.795837e-02
## 2007-01-31 -0.076315192 -2.135363e-03
## 2007-02-28 -0.156530036 -8.235021e-02
## 2007-03-30 -0.074197071 -1.724226e-05
## 2007-04-30 0.343624393 4.178042e-01
## 2007-05-31 0.144932574 2.191124e-01
## 2007-06-29 -0.105946135 -3.176631e-02
## 2007-07-31 -0.040207805 3.397202e-02
## 2007-08-31 -0.027429284 4.675054e-02
## 2007-09-28 0.192467578 2.666474e-01
## 2007-10-31 0.197120482 2.713003e-01
## 2007-11-30 -0.259237493 -1.850577e-01
## 2007-12-31 0.008392794 8.257262e-02
## 2008-01-31 -0.297615527 -2.234357e-01
## 2008-02-29 -0.175684605 -1.015048e-01
## 2008-03-31 0.234330487 3.085103e-01
## 2008-04-30 0.145732269 2.199121e-01
## 2008-05-30 -0.004506148 6.967368e-02
## 2008-06-30 -0.228295855 -1.541160e-01
## 2008-07-31 -0.012882649 6.129718e-02
## 2008-08-29 -0.067668224 6.511604e-03
## 2008-09-30 -0.172917390 -9.873756e-02
## 2008-10-31 -0.545658790 -4.714790e-01
## 2008-11-28 -0.303754659 -2.295748e-01
## 2008-12-31 0.286744780 3.609246e-01
## 2009-01-30 0.047196534 1.213764e-01
## 2009-02-27 -0.012067246 6.211258e-02
## 2009-03-31 0.173405677 2.475855e-01
## 2009-04-30 0.135510170 2.096900e-01
## 2009-05-29 -0.188282604 -1.141028e-01
## 2009-06-30 0.040861816 1.150416e-01
## 2009-07-31 0.064131809 1.383116e-01
## 2009-08-31 -0.041305553 3.287428e-02
## 2009-09-30 0.074614372 1.487942e-01
## 2009-10-30 0.171040347 2.452202e-01
## 2009-11-30 0.187867393 2.620472e-01
## 2009-12-31 -0.001488515 7.269131e-02
## 2010-01-29 -0.131965427 -5.778560e-02
## 2010-02-26 -0.047307341 2.687249e-02
## 2010-03-31 0.156201898 2.303817e-01
## 2010-04-30 0.070836693 1.450165e-01
## 2010-05-28 -0.239705767 -1.655259e-01
## 2010-06-30 -0.203925251 -1.297454e-01
## 2010-07-30 0.025845258 1.000251e-01
## 2010-08-31 0.016147461 9.032729e-02
## 2010-09-30 0.328572932 4.027528e-01
## 2010-10-29 0.105672903 1.798527e-01
## 2010-11-30 -0.004752261 6.942757e-02
## 2010-12-31 -0.105135444 -3.095562e-02
## 2011-01-31 0.041614166 1.157940e-01
## 2011-02-28 -0.001441982 7.273785e-02
## 2011-03-31 0.046510814 1.206906e-01
## 2011-04-29 0.093966962 1.681468e-01
## 2011-05-31 0.017010353 9.119018e-02#
# use first 3 eigen-vectors to compue three factor (with normalization to have pos correlation with market)
# note: cannot treat pc as a portfolio b/c weights do not sum to unity
p3 = pc.fit$loadings[, 1:3]
p3 # refers to B ## Comp.1 Comp.2 Comp.3
## AMZN 0.6142723 0.21913049 0.02511753
## HAS 0.2180206 0.07121327 -0.10335085
## JNJ 0.1191251 0.07479053 0.01728080
## MA 0.4084940 0.16754536 -0.67800274
## MNST 0.1664955 -0.64772681 -0.05484715
## MSFT 0.2481985 0.06152662 -0.17188237
## NFLX 0.4438769 -0.58645514 0.27136300
## VRTX 0.3282202 0.38207362 0.65000903colSums(p3)## Comp.1 Comp.2 Comp.3
## 2.54670313 -0.25790205 -0.04431275# create factor mimicking portfolio by normalizing weights to unity
p3 = p3/colSums(p3)
p3## Comp.1 Comp.2 Comp.3
## AMZN 0.2412029 -4.94508901 -0.097391758
## HAS -0.8453622 0.02796293 2.332304996
## JNJ -2.6882801 -0.28999587 0.006785559
## MA 0.1604011 -3.78097406 2.628915647
## MNST -0.6455767 -0.25433935 1.237728455
## MSFT -5.6010635 -0.23856583 -0.067492112
## NFLX 0.1742947 13.23445593 -1.052194032
## VRTX -1.2726543 0.15002676 -14.668668093barplot(p3[,1], horiz=F, main="Factor mimicking weights", col="pink", cex.names = 0.75, las=2)
# create first 3 factors
f3 = X_i %*% p3
head(f3)## Comp.1 Comp.2 Comp.3
## [1,] -0.19220160 -1.0723932 -0.7730953
## [2,] -0.32397533 -1.5200096 1.4878626
## [3,] -0.29450308 -1.9000210 -0.1139192
## [4,] -0.48393608 0.5779408 1.3702721
## [5,] -0.66017127 1.7272087 -2.8225888
## [6,] -0.08694696 -0.8857004 -0.4915262head(pc.fit$scores[, 1:3])## Comp.1 Comp.2 Comp.3
## 2006-06-30 0.04562732 0.06755789 -0.003999946
## 2006-07-31 -0.39814153 0.07826193 -0.099607154
## 2006-08-31 0.07436684 0.38695527 -0.127339653
## 2006-09-29 0.18331852 -0.12529221 -0.178636022
## 2006-10-31 0.26988699 0.04920042 0.144891876
## 2006-11-30 0.16295140 0.16770275 -0.171682568# chart.TimeSeries(f1, main="First principal component factor", colorset="green")
# estimate factor betas by multivariate regression
n.obs <- dim(X_i)[1]
X.mat = cbind(rep(1, n.obs), f3)
colnames(X.mat) = c("intercept", "Factor 1", "Factor 2", "Factor 3")
XX.mat = crossprod(X.mat)
# multivariate least squares
G.hat = solve(XX.mat)%*%crossprod(X.mat, X_i)
t(G.hat)## intercept Factor 1 Factor 2 Factor 3
## AMZN 0.033633161 -0.07351432 -0.0079492779 -0.010217222
## HAS 0.014518425 -0.07452909 0.0026532762 0.005809860
## JNJ 0.001978040 -0.04685116 -0.0004155822 -0.001363209
## MA 0.037835912 -0.09858280 -0.0130707619 0.022799102
## MNST 0.004204255 -0.09857632 0.0229773464 0.020257954
## MSFT -0.002083790 -0.13394453 -0.0002948429 0.011803247
## NFLX 0.023346763 -0.06044369 0.0693329756 0.003955469
## VRTX 0.007556629 -0.03241744 -0.0049012245 -0.061724186# can also use solve(qr(X.mat), returns.mat)
beta.hat = G.hat[2:4,]
beta.hat## AMZN HAS JNJ MA MNST
## Factor 1 -0.073514322 -0.074529094 -0.0468511563 -0.09858280 -0.09857632
## Factor 2 -0.007949278 0.002653276 -0.0004155822 -0.01307076 0.02297735
## Factor 3 -0.010217222 0.005809860 -0.0013632093 0.02279910 0.02025795
## MSFT NFLX VRTX
## Factor 1 -0.1339445266 -0.060443693 -0.032417440
## Factor 2 -0.0002948429 0.069332976 -0.004901224
## Factor 3 0.0118032473 0.003955469 -0.061724186B## [,1] [,2] [,3]
## [1,] 0.10929581 -0.030483746 0.003133449
## [2,] 0.03879183 -0.009906642 -0.012893169
## [3,] 0.02119561 -0.010404282 0.002155805
## [4,] 0.07268224 -0.023307620 -0.084581828
## [5,] 0.02962411 0.090106764 -0.006842262
## [6,] 0.04416130 -0.008559109 -0.021442576
## [7,] 0.07897783 0.081583122 0.033852929
## [8,] 0.05839933 -0.053151139 0.081089571E.hat = X_i - X.mat%*%G.hat
diagD.hat = diag(crossprod(E.hat)/(n.obs-4))
# compute covariance/correlation matrices with three pc factor
cov.pc3 = t(beta.hat)%*%var(f3)%*%(beta.hat) + diag(diagD.hat)
cov.pc3## AMZN HAS JNJ MA MNST
## AMZN 0.0183709596 0.0018822098 0.0015848218 0.002128499 0.0009535913
## HAS 0.0018822098 0.0060862053 0.0011124362 0.002151545 0.0019287963
## JNJ 0.0015848218 0.0011124362 0.0021251297 0.001274140 0.0009761364
## MA 0.0021284992 0.0021515450 0.0012741404 0.015601418 0.0021188736
## MNST 0.0009535913 0.0019287963 0.0009761364 0.002118874 0.0161849615
## MSFT 0.0034827756 0.0030067848 0.0020170356 0.004252950 0.0032487444
## NFLX -0.0009199986 0.0008902887 0.0002516180 -0.002642138 0.0049085665
## VRTX 0.0051970799 0.0014497783 0.0020845785 -0.001405968 -0.0011556113
## MSFT NFLX VRTX
## AMZN 0.0034827756 -0.0009199986 0.0051970799
## HAS 0.0030067848 0.0008902887 0.0014497783
## JNJ 0.0020170356 0.0002516180 0.0020845785
## MA 0.0042529504 -0.0026421382 -0.0014059682
## MNST 0.0032487444 0.0049085665 -0.0011556113
## MSFT 0.0061283715 0.0005144678 0.0023391923
## NFLX 0.0005144678 0.0185715836 0.0007272709
## VRTX 0.0023391923 0.0007272709 0.0159231690Sigma_PCA3## [,1] [,2] [,3] [,4] [,5]
## [1,] 0.017603623 0.0045013757 0.0026405076 0.0083893351 0.0004695590
## [2,] 0.004501376 0.0058612250 0.0008974926 0.0041409048 0.0003447362
## [3,] 0.002640508 0.0008974926 0.0020597991 0.0016007012 -0.0003243459
## [4,] 0.008389335 0.0041409048 0.0016007012 0.0150401900 0.0006317032
## [5,] 0.000469559 0.0003447362 -0.0003243459 0.0006317032 0.0155581163
## [6,] 0.005020369 0.0020743522 0.0009788509 0.0052228867 0.0006837210
## [7,] 0.006251063 0.0018190079 0.0008981495 0.0009754344 0.0094592080
## [8,] 0.008257138 0.0017464646 0.0019656221 -0.0013752833 -0.0036140853
## [,6] [,7] [,8]
## [1,] 0.0050203694 0.0062510631 0.008257138
## [2,] 0.0020743522 0.0018190079 0.001746465
## [3,] 0.0009788509 0.0008981495 0.001965622
## [4,] 0.0052228867 0.0009754344 -0.001375283
## [5,] 0.0006837210 0.0094592080 -0.003614085
## [6,] 0.0060980190 0.0020635905 0.001295147
## [7,] 0.0020635905 0.0183744713 0.003021136
## [8,] 0.0012951474 0.0030211361 0.015900071diag(cov.pc3)## AMZN HAS JNJ MA MNST MSFT
## 0.018370960 0.006086205 0.002125130 0.015601418 0.016184961 0.006128372
## NFLX VRTX
## 0.018571584 0.015923169# error difference between pca and empirical covariance matrix
norm(cov.pc3 - cov(X_i), "F")## [1] 0.01382605norm(Sigma_PCA3 - cov(X_i), "F")## [1] 0.01159327#==========================================================
# Plot perfromance
plotbt(models, plotX = T, log = 'y', LeftMargin = 3)
mtext('Cumulative Performance', side = 2, line = 1)
# Plot Strategy Statistics Side by Side
layout(1:1)
plotbt.strategy.sidebyside(models)
# Plot transition maps
layout(1:len(models))
for(m in names(models)) {
plotbt.transition.map(models[[m]]$weight, name=m)
legend('topright', legend = m, bty = 'n')
}
#```