拓端tecdat|R語言Fama-French三因子模型實(shí)際應(yīng)用：優(yōu)化投資組合
原文鏈接：http://tecdat.cn/?p=20360?
原文出處：拓端數(shù)據(jù)部落公眾號(hào)
本文將說明金融數(shù)學(xué)中的R 語言優(yōu)化投資組合，因子模型的實(shí)現(xiàn)和使用。
具有單一市場(chǎng)因素的宏觀經(jīng)濟(jì)因素模型
我們將從一個(gè)包含單個(gè)已知因子（即市場(chǎng)指數(shù)）的簡(jiǎn)單示例開始。該模型為
其中顯式因子ft為S＆P 500指數(shù)。我們將做一個(gè)簡(jiǎn)單的最小二乘（LS）回歸來估計(jì)截距α和加載β：
大多數(shù)代碼行用于準(zhǔn)備數(shù)據(jù)，而不是執(zhí)行因子建模。讓我們開始準(zhǔn)備數(shù)據(jù)：
# 設(shè)置開始結(jié)束日期和股票名稱列表begin_date <- "2016-01-01"end_date <- "2017-12-31"# 從YahooFinance下載數(shù)據(jù)data_set <- xts()for (stock_index in 1:length(stock_namelist))data_set <- cbind(data_set, Ad(getSymbols(stock_namelist[stock_index],from = begin_date, to = end_date,head(data_set)#> ? ? ? ? ? ? ? ?AAPL ?AMD ? ? ?ADI ? ? ABBV AEZS ? ? ? ?A ? ? ? APD ? ? ? AA ? ? ? CF#> 2016-01-04 98.74225 2.77 49.99239 49.46063 4.40 39.35598 107.89010 23.00764 35.13227#> 2016-01-05 96.26781 2.75 49.62508 49.25457 4.21 39.22057 105.96097 21.96506 34.03059#> 2016-01-06 94.38389 2.51 47.51298 49.26315 3.64 39.39467 103.38042 20.40121 31.08988#> 2016-01-07 90.40047 2.28 46.30082 49.11721 3.29 37.72138 ?99.91463 19.59558 29.61520#> 2016-01-08 90.87848 2.14 45.89677 47.77789 3.29 37.32482 ?99.39687 19.12169 29.33761#> 2016-01-11 92.35001 2.34 46.98954 46.25827 3.13 36.69613 ?99.78938 18.95583 28.14919head(SP500_index)#> ? ? ? ? ? ? ?index#> 2016-01-04 2012.66#> 2016-01-05 2016.71#> 2016-01-06 1990.26#> 2016-01-07 1943.09#> 2016-01-08 1922.03#> 2016-01-11 1923.67plot(SP500_index)
# 計(jì)算股票和SP500指數(shù)的對(duì)數(shù)收益率作為顯式因子X <- diff(log(data_set), na.pad = FALSE)N <- ncol(X) ?# 股票數(shù)量T <- nrow(X) ?# 天數(shù)
現(xiàn)在我們準(zhǔn)備進(jìn)行因子模型擬合。LS擬合很容易在R中實(shí)現(xiàn)，如下所示：
beta <- cov(X,f)/as.numeric(var(f))alpha <- colMeans(X) - beta*colMeans(f)sigma2 <- rep(NA, N)print(alpha)#> ? ? ? ? ? ? ?index#> AAPL ?0.0003999086#> AMD ? 0.0013825599#> ADI ? 0.0003609968#> ABBV ?0.0006684632#> AEZS -0.0022091301#> A ? ? 0.0002810616#> APD ? 0.0001786375#> AA ? ?0.0006429140#> CF ? -0.0006029705print(beta)#> ? ? ? ? ?index#> AAPL 1.0957919#> AMD ?2.1738304#> ADI ?1.2683047#> ABBV 0.9022748#> AEZS 1.7115761#> A ? ?1.3277212#> APD ?1.0239453#> AA ? 1.8593524#> CF ? 1.5702493
或者，我們可以使用矩陣表示法進(jìn)行擬合
，我們定義
和擴(kuò)展因子
。然后最小化
t(X) %*% F_ %*% solve(t(F_) %*% F_)#> ? ? ? ? ? ? ?alpha ? ? ?beta#> AAPL ?0.0003999086 1.0957919#> AMD ? 0.0013825599 2.1738304#> ADI ? 0.0003609968 1.2683047#> ABBV ?0.0006684632 0.9022748#> AEZS -0.0022091301 1.7115761#> A ? ? 0.0002810616 1.3277212#> APD ? 0.0001786375 1.0239453#> AA ? ?0.0006429140 1.8593524#> CF ? -0.0006029705 1.5702493E <- xts(t(t(X) - Gamma %*% t(F_)), index(X)) ?# 殘差
另外，我們可以簡(jiǎn)單地使用R為我們完成工作：
cbind(alpha = factor_model$alpha, beta = factor_model$beta)#> ? ? ? ? ? ? ?alpha ? ? index#> AAPL ?0.0003999086 1.0957919#> AMD ? 0.0013825599 2.1738304#> ADI ? 0.0003609968 1.2683047#> ABBV ?0.0006684632 0.9022748#> AEZS -0.0022091301 1.7115761#> A ? ? 0.0002810616 1.3277212#> APD ? 0.0001786375 1.0239453#> AA ? ?0.0006429140 1.8593524#> CF ? -0.0006029705 1.5702493
可視化協(xié)方差矩陣
有趣的是，可視化對(duì)數(shù)收益率[算術(shù)處理誤差]
以及殘差Ψ的估計(jì)協(xié)方差矩陣。讓我們從對(duì)數(shù)收益率的協(xié)方差矩陣開始：
main = "單因子模型對(duì)數(shù)收益的協(xié)方差矩陣")
我們可以觀察到所有股票都是高度相關(guān)的，這是市場(chǎng)因素的影響。為了檢查股票相關(guān)關(guān)系，我們繪制相關(guān)圖：
plot(cov2cor(Psi),main = "殘差協(xié)方差矩陣")
cbind(stock_namelist, sector_namelist) ?# 股票的行業(yè)#> ? ? ? stock_namelist sector_namelist#> ?[1,] "AAPL" ? ? ? ? "Information Technology"#> ?[2,] "AMD" ? ? ? ? ?"Information Technology"#> ?[3,] "ADI" ? ? ? ? ?"Information Technology"#> ?[4,] "ABBV" ? ? ? ? "Health Care"#> ?[5,] "AEZS" ? ? ? ? "Health Care"#> ?[6,] "A" ? ? ? ? ? ?"Health Care"#> ?[7,] "APD" ? ? ? ? ?"Materials"#> ?[8,] "AA" ? ? ? ? ? "Materials"#> ?[9,] "CF" ? ? ? ? ? "Materials"
有趣的是，我們可以觀察到對(duì)Ψ執(zhí)行的自動(dòng)聚類可以正確識(shí)別股票的行業(yè)。
評(píng)估投資資金
在此示例中，我們將基于因子模型評(píng)估幾種投資基金的績(jī)效。我們將標(biāo)準(zhǔn)普爾500指數(shù)作為明確的市場(chǎng)因素，并假設(shè)無風(fēng)險(xiǎn)收益為零 rf = 0。特別是，我們考慮六種交易所買賣基金（ETF）：
我們首先加載數(shù)據(jù)：
# 設(shè)置開始結(jié)束日期和股票名稱列表begin_date <- "2016-10-01"end_date <- "2017-06-30"# 從YahooFinance下載數(shù)據(jù)data_set <- xts()for (stock_index in 1:length(stock_namelist))data_set <- cbind(data_set, Ad(getSymbols(stock_namelist[stock_index],head(data_set)#> ? ? ? ? ? ? ? ? SPY ? XIVH ? ? SPHB ? ? SPLV ? ? USMV ? ? ?JKD#> 2016-10-03 203.6610 29.400 31.38322 38.55683 42.88382 119.8765#> 2016-10-04 202.6228 30.160 31.29729 38.10687 42.46553 119.4081#> 2016-10-05 203.5195 30.160 31.89880 38.02249 42.37048 119.9421#> 2016-10-06 203.6610 30.160 31.83196 38.08813 42.39899 120.0826#> 2016-10-07 202.9626 30.670 31.58372 37.98500 42.35146 119.8296#> 2016-10-10 204.0197 31.394 31.87970 38.18187 42.56060 120.5978head(SP500_index)#> ? ? ? ? ? ? ?index#> 2016-10-03 2161.20#> 2016-10-04 2150.49#> 2016-10-05 2159.73#> 2016-10-06 2160.77#> 2016-10-07 2153.74#> 2016-10-10 2163.66# 計(jì)算股票和SP500指數(shù)的對(duì)數(shù)收益率作為顯式因子X <- diff(log(data_set), na.pad = FALSE)N <- ncol(X) ?# 股票數(shù)量T <- nrow(X) ?# 天數(shù)
現(xiàn)在我們可以計(jì)算所有ETF的alpha和beta：
#> ? ? ? ? ? ? ?alpha ? ? ?beta#> SPY ? 7.142225e-05 1.0071424#> XIVH ?1.810392e-03 2.4971086#> SPHB -2.422107e-04 1.5613533#> SPLV ?1.070918e-04 0.6777149#> USMV ?1.166177e-04 0.6511667#> JKD ? 2.569578e-04 0.8883843
現(xiàn)在可以進(jìn)行一些觀察：
SPY是S＆P 500的ETF，如預(yù)期的那樣，其alpha值幾乎為零，beta值幾乎為1： α= 7.142211×10-5和 β= 1.0071423。
XIVH是具有高alpha值的ETF，計(jì)算出的alpha值是ETF中最高的（高1-2個(gè)數(shù)量級(jí)）： α= 1.810392×10-3。
SPHB是一種ETF，據(jù)推測(cè)具有很高的beta，而計(jì)算出的beta卻是最高的，但不是最高的：β= 1.5613531。有趣的是，計(jì)算出的alpha為負(fù)，因此，該ETF應(yīng)謹(jǐn)慎。
SPLV是降低波動(dòng)性的ETF，實(shí)際上，計(jì)算得出的beta偏低：β= 0.6777072。
USMV還是降低波動(dòng)性的ETF，實(shí)際上，計(jì)算出的beta是最低的：β= 0.6511671。
JKD顯示出很好的折衷。
我們可以使用一些可視化：
barplot(rev(alpha), horiz = TRUE, main = "alph
我們還可以使用例如Sharpe比率，以更系統(tǒng)的比較不同的ETF。回顧一種資產(chǎn)和一個(gè)因素的因子模型
我們獲得
夏普比率如下：
假設(shè)
。因此，基于Sharpe比率對(duì)不同資產(chǎn)進(jìn)行排名的一種方法是根據(jù)α/β比率對(duì)它們進(jìn)行排名：
print(ranking)#> ? ? ? ? alpha/beta ? ? ? ? SR ? ? ? ? alpha ? ? ?beta#> XIVH ?7.249952e-04 0.13919483 ?1.810392e-03 2.4971086#> JKD ? 2.892417e-04 0.17682677 ?2.569578e-04 0.8883843#> USMV ?1.790904e-04 0.12280053 ?1.166177e-04 0.6511667#> SPLV ?1.580189e-04 0.10887903 ?1.070918e-04 0.6777149#> SPY ? 7.091574e-05 0.14170591 ?7.142225e-05 1.0071424#> SPHB -1.551287e-04 0.07401566 -2.422107e-04 1.5613533
可以看到：
就α/β而言，XIVH最佳（α最大），而SPHB最差（α負(fù)）。
就夏普比率（更確切地說，是信息比率，因?yàn)槲覀兒雎粤藷o風(fēng)險(xiǎn)利率）而言，JDK是最好的，其次是SPY。這證實(shí)了大多數(shù)投資基金的表現(xiàn)不超過市場(chǎng)的觀點(diǎn)。
顯然，無論以哪種衡量標(biāo)準(zhǔn)，SPHB都是最差的：負(fù)α，負(fù)β比率和Sharpe比率。
JDK之所以能夠取得最佳性能，是因?yàn)樗腶lpha值很好（盡管不是最好的），而同時(shí)具有0.88的中等beta值。
XIVH和SPHB有大量不同的beta，因此在市場(chǎng)上具有極端敞口。
USMV在市場(chǎng)上的曝光率最小，有可接受的alpha值，并且其Sharpe比率接近第二和第三高的位置。
Fama-French三因子模型
該示例將說明使用標(biāo)準(zhǔn)普爾500指數(shù)中的九種股票的Fama-French三因子模型。讓我們從加載數(shù)據(jù)開始：
# 設(shè)置開始結(jié)束日期和股票名稱列表begin_date <- "2013-01-01"end_date <- "2017-08-31"# 從YahooFinance下載數(shù)據(jù)data_set <- xts()for (stock_index in 1:length(stock_namelist))data_set <- cbind(data_set, Ad(getSymbols(stock_namelist[stock_index],# 下載Fama-French因子head(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.56tail(fama_lib)#> ? ? ? ? ? ?Mkt.RF ? SMB ? HML#> 2017-11-22 ?-0.05 ?0.10 -0.04#> 2017-11-24 ? 0.21 ?0.02 -0.44#> 2017-11-27 ?-0.06 -0.36 ?0.03#> 2017-11-28 ? 1.06 ?0.38 ?0.84#> 2017-11-29 ? 0.02 ?0.04 ?1.45#> 2017-11-30 ? 0.82 -0.56 -0.50# 計(jì)算股票的對(duì)數(shù)收益率和Fama-French因子X <- diff(log(data_set), na.pad = FALSE)N <- ncol(X) ?#股票數(shù)量
現(xiàn)在我們?cè)诰仃嘑中具有三個(gè)因子，并希望擬合模型
，其中現(xiàn)在的載荷是一個(gè)beta矩陣：
。我們可以做最小二乘擬合，最小化
。更方便地，我們定義
和擴(kuò)展因子?
。然后可以將LS公式寫為最小化
print(Gamma)#> ? ? ? ? ? ? ?alpha ? ? ? ?b1 ? ? ? ? ?b2 ? ? ? ? ?b3#> AAPL ?1.437845e-04 0.9657612 -0.23339130 -0.49806858#> AMD ? 6.181760e-04 1.4062105 ?0.80738336 -0.07240117#> ADI ?-2.285017e-05 1.2124008 ?0.09025928 -0.20739271#> ABBV ?1.621380e-04 1.0582340 ?0.02833584 -0.72152627#> AEZS -4.513235e-03 0.6989534 ?1.31318108 -0.25160182#> A ? ? 1.146100e-05 1.2181429 ?0.10370898 -0.20487290#> APD ? 6.281504e-05 1.0222936 -0.04394061 ?0.11060938#> AA ? -4.587722e-05 1.3391852 ?0.62590136 ?0.99858692#> CF ? -5.777426e-04 1.0387867 ?0.48430007 ?0.82014523
另外，我們可以使用R完成：
#> ? ? ? ? ? ? ?alpha ? ?Mkt.RF ? ? ? ? SMB ? ? ? ? HML#> AAPL ?1.437845e-04 0.9657612 -0.23339130 -0.49806858#> AMD ? 6.181760e-04 1.4062105 ?0.80738336 -0.07240117#> ADI ?-2.285017e-05 1.2124008 ?0.09025928 -0.20739271#> ABBV ?1.621380e-04 1.0582340 ?0.02833584 -0.72152627#> AEZS -4.513235e-03 0.6989534 ?1.31318108 -0.25160182#> A ? ? 1.146100e-05 1.2181429 ?0.10370898 -0.20487290#> APD ? 6.281504e-05 1.0222936 -0.04394061 ?0.11060938#> AA ? -4.587722e-05 1.3391852 ?0.62590136 ?0.99858692#> CF ? -5.777426e-04 1.0387867 ?0.48430007 ?0.82014523
統(tǒng)計(jì)因子模型
現(xiàn)在讓我們考慮統(tǒng)計(jì)因子模型或隱式因子模型，其中因子和載荷均不可用。調(diào)用具有 K因子的模型 XT =α1T+ BFT + ET的主成分方法：
PCA：
樣本均值：
矩陣：
樣本協(xié)方差矩陣：
特征分解：
估計(jì)：
?
更新特征分解：
重復(fù)步驟2-3，直到收斂為止。
#> ? ? ? ? ? ? ?alpha#> AAPL ?0.0007074564 0.0002732114 -0.004631647 -0.0044814226#> AMD ? 0.0013722468 0.0045782146 -0.035202146 ?0.0114549515#> ADI ? 0.0006533116 0.0004151904 -0.007379066 -0.0053058139#> ABBV ?0.0007787929 0.0017513359 -0.003967816 -0.0056000810#> AEZS -0.0041576357 0.0769496344 ?0.002935950 ?0.0006249473#> A ? ? 0.0006902482 0.0012690079 -0.005680162 -0.0061507654#> APD ? 0.0006236565 0.0005442926 -0.004229364 -0.0057976394#> AA ? ?0.0006277163 0.0027405024 -0.009796620 -0.0149177957#> CF ? -0.0000573028 0.0023108605 -0.007409061 -0.0153425661
同樣，我們可以使用R完成工作：
#> ? ? ? ? ? ? ?alpha ? ? ?factor1 ? ? ?factor2 ? ? ? factor3#> AAPL ?0.0007074564 0.0002732114 -0.004631647 -0.0044814226#> AMD ? 0.0013722468 0.0045782146 -0.035202146 ?0.0114549515#> ADI ? 0.0006533116 0.0004151904 -0.007379066 -0.0053058139#> ABBV ?0.0007787929 0.0017513359 -0.003967816 -0.0056000810#> AEZS -0.0041576357 0.0769496344 ?0.002935950 ?0.0006249473#> A ? ? 0.0006902482 0.0012690079 -0.005680162 -0.0061507654#> APD ? 0.0006236565 0.0005442926 -0.004229364 -0.0057976394#> AA ? ?0.0006277163 0.0027405024 -0.009796620 -0.0149177957#> CF ? -0.0000573028 0.0023108605 -0.007409061 -0.0153425661
通過不同因子模型進(jìn)行協(xié)方差矩陣估計(jì)的最終比較
我們最終將比較以下不同的因子模型：
樣本協(xié)方差矩陣
宏觀經(jīng)濟(jì)一因素模型
基本的三因素Fama-French模型
統(tǒng)計(jì)因素模型
我們?cè)谟?xùn)練階段估計(jì)模型，然后將估計(jì)的協(xié)方差矩陣與測(cè)試階段的樣本協(xié)方差矩陣進(jìn)行比較。估計(jì)誤差將根據(jù)PRIAL（平均損失提高百分比）進(jìn)行評(píng)估：
加載訓(xùn)練和測(cè)試集：
# 設(shè)置開始結(jié)束日期和股票名稱列表begin_date <- "2013-01-01"end_date <- "2015-12-31"# 準(zhǔn)備股票數(shù)據(jù)data_set <- xts()for (stock_index in 1:length(stock_namelist))data_set <- cbind(data_set, Ad(getSymbols(stock_namelist[stock_index],# ? Fama-French 因子mydata <- mydata[-nrow(mydata),# 準(zhǔn)備指數(shù)f_SP500 <- diff(log(SP500_index), na.pad = FALSE)# 將數(shù)據(jù)拆分為訓(xùn)練數(shù)據(jù)和測(cè)試數(shù)據(jù)T_trn <- round(0.45*T)X_trn <- X[1:T_trn, ]X_tst <- X[(T_trn+1):T, ]
?現(xiàn)在讓我們用訓(xùn)練數(shù)據(jù)估算不同的因子模型：
# 樣本協(xié)方差矩陣Sigma_SCM <- cov(X_trn)# 單因素模型Gamma <- t(solve(t(F_) %*% F_, t(F_) %*% X_trn))E <- xts(t(t(X_trn) - Gamma %*% t(F_)), index(X_trn))# Fama-French三因子模型Sigma_FamaFrench <- B %*% cov(F_FamaFrench_trn) %*% t(B) + diag(diag(Psi))# 統(tǒng)計(jì)單因子模型while (norm(Sigma - Sigma_prev, "F")/norm(Sigma, "F") > 1e-3) {B <- eigSigma$vectors[, 1:K, drop = FALSE] %*% diag(sqrt(eigSigma$values[1:K]), K, K)# 統(tǒng)計(jì)三因子模型K <- 3while (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_PCA3 <- Sigma# 統(tǒng)計(jì)五因子模型K <- 5eigSigma <- 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)))
最后，讓我們比較測(cè)試數(shù)據(jù)中的不同估計(jì)：
Sigma_true <- cov(X_tst)barplot(error, main = "協(xié)方差矩陣估計(jì)誤差",
PRIAL <- 100*(ref - error^2)/refbarplot(PRIAL, main = "協(xié)方差矩陣估計(jì)的先驗(yàn)方法",
最終可以看到使用因子模型進(jìn)行協(xié)方差矩陣估計(jì)會(huì)有所幫助。
最受歡迎的見解
1.用機(jī)器學(xué)習(xí)識(shí)別不斷變化的股市狀況—隱馬爾科夫模型(HMM)的應(yīng)用
2.R語言GARCH-DCC模型和DCC（MVT）建模估計(jì)
3.R語言實(shí)現(xiàn) Copula 算法建模依賴性案例分析報(bào)告
4.R語言COPULAS和金融時(shí)間序列數(shù)據(jù)VaR分析
5.R語言多元COPULA GARCH 模型時(shí)間序列預(yù)測(cè)
6.用R語言實(shí)現(xiàn)神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)股票實(shí)例
7.r語言預(yù)測(cè)波動(dòng)率的實(shí)現(xiàn)：ARCH模型與HAR-RV模型
8.R語言如何做馬爾科夫轉(zhuǎn)換模型markov switching model
9.matlab使用Copula仿真優(yōu)化市場(chǎng)風(fēng)險(xiǎn)
