library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.6 ✓ dplyr 1.0.8
## ✓ tidyr 1.2.0 ✓ stringr 1.4.0
## ✓ readr 2.1.2 ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(tidyr)
library(dplyr)
library(readr)
fac <- read.csv("m-fac9003(3).csv")
head(fac)
## AA AGE CAT F FDX GM HPQ KMB MEL NYT PG
## 1 -16.40 -12.17 -4.44 -0.06 -2.28 -2.12 -6.19 -11.01 -10.77 -6.30 -8.89
## 2 4.04 4.95 8.84 6.02 10.47 8.97 -4.01 -5.20 0.34 -4.62 -0.84
## 3 0.12 13.08 0.17 2.06 10.84 1.57 5.67 3.21 -0.17 -0.66 5.41
## 4 -4.28 -11.06 0.25 -5.67 -2.44 -4.19 -5.29 -0.65 -2.20 -10.60 4.26
## 5 5.81 19.70 8.52 3.89 -16.17 10.94 8.81 8.83 11.85 11.59 16.35
## 6 -4.05 -1.44 -22.10 -5.79 -2.81 -2.70 -1.47 1.55 -7.76 -0.12 4.80
## TRB TXN SP500
## 1 -13.04 -7.61 -7.52
## 2 -0.37 4.97 0.21
## 3 2.36 2.69 1.77
## 4 -7.98 -6.85 -3.34
## 5 8.82 22.88 8.55
## 6 -0.64 -5.87 -1.53
t = dim(fac)[1]
t
## [1] 168
market = fac[,14]
market
## [1] -7.52 0.21 1.77 -3.34 8.55 -1.53 -1.16 -10.05 -5.73 -1.27
## [11] 5.40 1.92 3.63 6.23 1.73 -0.44 3.40 -5.25 4.02 1.52
## [21] -2.35 0.77 -4.77 10.82 -2.31 0.64 -2.52 2.48 -0.21 -2.04
## [31] 3.67 -2.66 0.67 -0.03 2.77 0.74 0.46 0.80 1.62 -2.78
## [41] 2.03 -0.18 -0.79 3.19 -1.24 1.69 -1.55 0.75 3.00 -3.28
## [51] -4.87 0.85 0.90 -3.03 2.79 3.39 -3.08 1.67 -4.39 0.76
## [61] 1.95 3.13 2.26 2.33 3.16 1.67 2.73 -0.48 3.57 -0.94
## [71] 3.66 1.32 2.85 0.29 0.38 0.93 1.87 -0.20 -5.00 1.46
## [81] 4.99 2.20 6.92 -2.56 5.71 0.18 -4.69 5.41 5.44 3.93
## [91] 7.39 -6.17 4.90 -3.86 4.03 1.14 0.60 6.62 4.58 0.50
## [101] -2.30 3.53 -1.58 -14.99 5.86 7.70 5.55 5.27 3.74 -3.60
## [111] 3.51 3.44 -2.87 5.06 -3.58 -1.02 -3.24 5.85 1.48 5.35
## [121] -5.53 -2.47 9.20 -3.55 -2.67 1.92 -2.13 5.56 -5.85 -1.00
## [131] -8.52 -0.08 3.03 -9.64 -6.79 7.36 0.21 -2.79 -1.37 -6.69
## [141] -8.39 1.63 7.36 0.62 -1.69 -2.22 3.52 -6.29 -1.05 -7.39
## [151] -8.04 0.35 -11.14 8.51 5.60 -6.13 -2.84 -1.80 0.74 8.01
## [161] 5.00 1.06 1.55 1.71 -1.27 5.42 0.64 5.00
fac1 = fac[,c(-14)]
head(fac1)
## AA AGE CAT F FDX GM HPQ KMB MEL NYT PG
## 1 -16.40 -12.17 -4.44 -0.06 -2.28 -2.12 -6.19 -11.01 -10.77 -6.30 -8.89
## 2 4.04 4.95 8.84 6.02 10.47 8.97 -4.01 -5.20 0.34 -4.62 -0.84
## 3 0.12 13.08 0.17 2.06 10.84 1.57 5.67 3.21 -0.17 -0.66 5.41
## 4 -4.28 -11.06 0.25 -5.67 -2.44 -4.19 -5.29 -0.65 -2.20 -10.60 4.26
## 5 5.81 19.70 8.52 3.89 -16.17 10.94 8.81 8.83 11.85 11.59 16.35
## 6 -4.05 -1.44 -22.10 -5.79 -2.81 -2.70 -1.47 1.55 -7.76 -0.12 4.80
## TRB TXN
## 1 -13.04 -7.61
## 2 -0.37 4.97
## 3 2.36 2.69
## 4 -7.98 -6.85
## 5 8.82 22.88
## 6 -0.64 -5.87
fac1 = as.matrix(fac1)
head(fac1)
## AA AGE CAT F FDX GM HPQ KMB MEL NYT PG
## [1,] -16.40 -12.17 -4.44 -0.06 -2.28 -2.12 -6.19 -11.01 -10.77 -6.30 -8.89
## [2,] 4.04 4.95 8.84 6.02 10.47 8.97 -4.01 -5.20 0.34 -4.62 -0.84
## [3,] 0.12 13.08 0.17 2.06 10.84 1.57 5.67 3.21 -0.17 -0.66 5.41
## [4,] -4.28 -11.06 0.25 -5.67 -2.44 -4.19 -5.29 -0.65 -2.20 -10.60 4.26
## [5,] 5.81 19.70 8.52 3.89 -16.17 10.94 8.81 8.83 11.85 11.59 16.35
## [6,] -4.05 -1.44 -22.10 -5.79 -2.81 -2.70 -1.47 1.55 -7.76 -0.12 4.80
## TRB TXN
## [1,] -13.04 -7.61
## [2,] -0.37 4.97
## [3,] 2.36 2.69
## [4,] -7.98 -6.85
## [5,] 8.82 22.88
## [6,] -0.64 -5.87
num = rep(1,t)
num
## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [112] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [149] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
xmat = cbind(num, market)
head(xmat)
## num market
## [1,] 1 -7.52
## [2,] 1 0.21
## [3,] 1 1.77
## [4,] 1 -3.34
## [5,] 1 8.55
## [6,] 1 -1.53
betta_hat = solve(t(xmat)%*%
xmat)%*%
t(xmat)%*%
fac1
betta_hat
## AA AGE CAT F FDX GM HPQ
## num 0.549124 0.7218061 0.8393521 0.4543643 0.7995790 0.1982025 0.6835681
## market 1.291591 1.5141359 0.9406928 1.2192453 0.8051166 1.0457019 1.6279512
## KMB MEL NYT PG TRB TXN
## num 0.5463020 0.8849263 0.4904120 0.8880914 0.6512465 1.438887
## market 0.5498052 1.1228708 0.7706495 0.4688034 0.7178808 1.796412
E_hat = fac1 - xmat%*%
betta_hat
head(E_hat)
## AA AGE CAT F FDX GM
## [1,] -7.2363588 -1.5055042 1.7946575 8.6543606 2.9748981 5.5454755
## [2,] 3.2196419 3.9102254 7.8031024 5.3095942 9.5013465 8.5522001
## [3,] -2.7152402 9.6781734 -2.3343784 -0.5524285 8.6153645 -0.4790948
## [4,] -0.5152096 -6.7245922 2.5525617 -2.0520849 -0.5504894 -0.8955583
## [5,] -5.7822280 6.0323321 -0.3622754 -6.9889119 -23.8533263 1.8010466
## [6,] -2.6229896 0.1548218 -21.5000922 -4.3789190 -2.3777506 -1.2982786
## HPQ KMB MEL NYT PG TRB TXN
## [1,] 5.3686247 -7.4217666 -3.2109382 -0.9951281 -6.252690 -8.2927829 4.460129
## [2,] -5.0354379 -5.8617611 -0.7807291 -5.2722484 -1.826540 -1.1720014 3.153867
## [3,] 2.1049583 1.6905428 -3.0424075 -2.5144615 3.692127 0.4381045 -1.928535
## [4,] -0.5362112 0.6400475 0.6654621 -8.5164428 4.937712 -6.2335246 -2.288872
## [5,] -5.7925506 3.5828633 1.3645288 4.5105352 11.453640 2.0308727 6.081793
## [6,] 0.3371972 1.8449000 -6.9269340 0.5686817 4.629178 -0.1928888 -4.560377
diagD_hat = diag(t(E_hat)%*%
E_hat)/(t-2)
ret = apply(fac1, 2, var)
r2 = 1- diag(t(E_hat)%*%
E_hat)/((t-1)*ret)
rsq = sqrt(diagD_hat)
cov_factor = var(market)*t(betta_hat)%*%betta_hat + diag(diagD_hat)
one.vec = rep(1,13)
a = solve(cov_factor)%*%one.vec
b = t(one.vec)%*%a
gmv = a / as.numeric(b)
gmv
## [,1]
## AA 0.03501615
## AGE -0.03373670
## CAT 0.05314427
## F 0.05576422
## FDX 0.06451069
## GM 0.12369835
## HPQ -0.03290554
## KMB 0.28022778
## MEL 0.01901912
## NYT 0.19373120
## PG 0.19019732
## TRB 0.14314200
## TXN -0.09180887