library(readr)
retdata <- read_csv("m-fac9003.csv")
## Rows: 168 Columns: 14
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## dbl (14): AA, AGE, CAT, F, FDX, GM, HPQ, KMB, MEL, NYT, PG, TRB, TXN, SP500
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(retdata)
## # A tibble: 6 × 14
## AA AGE CAT F FDX GM HPQ KMB MEL NYT PG
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 -16.4 -12.2 -4.44 -0.06 -2.28 -2.12 -6.19 -11.0 -10.8 -6.3 -8.89
## 2 4.04 4.95 8.84 6.02 10.5 8.97 -4.01 -5.2 0.34 -4.62 -0.84
## 3 0.12 13.1 0.17 2.06 10.8 1.57 5.67 3.21 -0.17 -0.66 5.41
## 4 -4.28 -11.1 0.25 -5.67 -2.44 -4.19 -5.29 -0.65 -2.2 -10.6 4.26
## 5 5.81 19.7 8.52 3.89 -16.2 10.9 8.81 8.83 11.8 11.6 16.4
## 6 -4.05 -1.44 -22.1 -5.79 -2.81 -2.7 -1.47 1.55 -7.76 -0.12 4.8
## # ℹ 3 more variables: TRB <dbl>, TXN <dbl>, SP500 <dbl>
retdata <- retdata
head(retdata)
## # A tibble: 6 × 14
## AA AGE CAT F FDX GM HPQ KMB MEL NYT PG
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 -16.4 -12.2 -4.44 -0.06 -2.28 -2.12 -6.19 -11.0 -10.8 -6.3 -8.89
## 2 4.04 4.95 8.84 6.02 10.5 8.97 -4.01 -5.2 0.34 -4.62 -0.84
## 3 0.12 13.1 0.17 2.06 10.8 1.57 5.67 3.21 -0.17 -0.66 5.41
## 4 -4.28 -11.1 0.25 -5.67 -2.44 -4.19 -5.29 -0.65 -2.2 -10.6 4.26
## 5 5.81 19.7 8.52 3.89 -16.2 10.9 8.81 8.83 11.8 11.6 16.4
## 6 -4.05 -1.44 -22.1 -5.79 -2.81 -2.7 -1.47 1.55 -7.76 -0.12 4.8
## # ℹ 3 more variables: TRB <dbl>, TXN <dbl>, SP500 <dbl>
stockret <- as.matrix(retdata[, 1:13])
mktret <- as.matrix(retdata[, 14])
head(mktret,5)
## SP500
## [1,] -7.52
## [2,] 0.21
## [3,] 1.77
## [4,] -3.34
## [5,] 8.55
dim(mktret)
## [1] 168 1
TT <- dim(stockret)
TT
## [1] 168 13
# Method 2: Use formula "inv(X'X)*X'Y" ----
ones = rep(1, 168)
X = as.matrix(cbind(ones, mktret))
Y = stockret
b_hat = solve(t(X)%*%X)%*%t(X) %*% stockret
b_hat
## AA AGE CAT F FDX GM HPQ
## ones 0.549124 0.7218061 0.8393521 0.4543643 0.7995790 0.1982025 0.6835681
## SP500 1.291591 1.5141359 0.9406928 1.2192453 0.8051166 1.0457019 1.6279512
## KMB MEL NYT PG TRB TXN
## ones 0.5463020 0.8849263 0.4904120 0.8880914 0.6512465 1.438887
## SP500 0.5498052 1.1228708 0.7706495 0.4688034 0.7178808 1.796412
E_hat = Y - X %*% b_hat
head(E_hat,5)
## AA AGE CAT F FDX GM
## [1,] -7.2363588 -1.505504 1.7946575 8.6543606 2.9748981 5.5454755
## [2,] 3.2196419 3.910225 7.8031024 5.3095942 9.5013465 8.5522001
## [3,] -2.7152402 9.678173 -2.3343784 -0.5524285 8.6153645 -0.4790948
## [4,] -0.5152096 -6.724592 2.5525617 -2.0520849 -0.5504894 -0.8955583
## [5,] -5.7822280 6.032332 -0.3622754 -6.9889119 -23.8533263 1.8010466
## 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
# Excluding constant term
b_hat = as.matrix(b_hat[-1,])
b_hat
## [,1]
## AA 1.2915911
## AGE 1.5141359
## CAT 0.9406928
## F 1.2192453
## FDX 0.8051166
## GM 1.0457019
## HPQ 1.6279512
## KMB 0.5498052
## MEL 1.1228708
## NYT 0.7706495
## PG 0.4688034
## TRB 0.7178808
## TXN 1.7964117
diagD_hat = diag(t(E_hat) %*% E_hat)/(168-2) ### find it out
diagD_hat
## AA AGE CAT F FDX GM HPQ KMB
## 59.19846 60.95651 59.66741 67.91031 78.39072 66.09876 89.66712 36.84610
## MEL NYT PG TRB TXN
## 37.45483 43.43290 41.71711 52.05835 131.65239
# Covariance matrix by single factor model
cov_1f.1 = as.numeric(var(mktret))*b_hat%*%t(b_hat) + diag(diagD_hat);
var(mktret)
## SP500
## SP500 18.74401
cov_1f.1
## AA AGE CAT F FDX GM HPQ
## AA 90.46736 36.65662 22.773795 29.51744 19.491551 25.31602 39.41205
## AGE 36.65662 103.92918 26.697784 34.60338 22.850000 29.67804 46.20285
## CAT 22.77380 26.69778 76.254044 21.49817 14.196104 18.43819 28.70462
## F 29.51744 34.60338 21.498168 95.77439 18.399772 23.89800 37.20446
## FDX 19.49155 22.85000 14.196104 18.39977 90.540833 15.78081 24.56760
## GM 25.31602 29.67804 18.438188 23.89800 15.780808 86.59520 31.90890
## HPQ 39.41205 46.20285 28.704615 37.20446 24.567600 31.90890 139.34297
## KMB 13.31056 15.60401 9.694362 12.56500 8.297174 10.77654 16.77694
## MEL 27.18425 31.86817 19.798858 25.66158 16.945373 22.00899 34.26366
## NYT 18.65711 21.87179 13.588366 17.61207 11.629960 15.10523 23.51586
## PG 11.34954 13.30510 8.266108 10.71382 7.074766 9.18885 14.30522
## TRB 17.37961 20.37416 12.657930 16.40612 10.833622 14.07093 21.90566
## TXN 43.49041 50.98393 31.674972 41.05438 27.109857 35.21083 54.81631
## KMB MEL NYT PG TRB TXN
## AA 13.310564 27.184251 18.657114 11.349541 17.379606 43.49041
## AGE 15.604012 31.868174 21.871787 13.305099 20.374161 50.98393
## CAT 9.694362 19.798858 13.588366 8.266108 12.657930 31.67497
## F 12.565001 25.661582 17.612075 10.713820 16.406124 41.05438
## FDX 8.297174 16.945373 11.629960 7.074766 10.833622 27.10986
## GM 10.776539 22.008995 15.105229 9.188850 14.070929 35.21083
## HPQ 16.776942 34.263656 23.515856 14.305223 21.905656 54.81631
## KMB 42.512149 11.571807 7.941971 4.831279 7.398161 18.51302
## MEL 11.571807 61.088004 16.219938 9.866952 15.109311 37.80926
## NYT 7.941971 16.219938 54.564979 6.771894 10.369833 25.94928
## PG 4.831279 9.866952 6.771894 45.836602 6.308202 15.78553
## TRB 7.398161 15.109311 10.369833 6.308202 61.718134 24.17246
## TXN 18.513021 37.809262 25.949280 15.785529 24.172455 192.14110
#GLOBAL MINIMUM VARIANCE PORTFOLIO
library('pracma')
library('dplyr')
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
one <- rep(1,13)
one_13_1 <- matrix(one,ncol=1)
a <- inv(cov_1f.1)%*%one_13_1 # 13x13 x 13x1 = 13x1
b <- t(one_13_1)%*%inv(cov_1f.1)%*%one_13_1 # 1x13 x 13x13 x 13x1 = 1x1
mvp <- a/as.vector(b)
as_tibble(mvp)
## Warning: The `x` argument of `as_tibble.matrix()` must have unique column names if
## `.name_repair` is omitted as of tibble 2.0.0.
## ℹ Using compatibility `.name_repair`.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## # A tibble: 13 × 1
## V1
## <dbl>
## 1 0.0117
## 2 -0.0306
## 3 0.0792
## 4 0.0225
## 5 0.0802
## 6 0.0533
## 7 -0.0354
## 8 0.250
## 9 0.0703
## 10 0.154
## 11 0.243
## 12 0.140
## 13 -0.0388
colnames(mvp)="Weight"
mvp
## Weight
## AA 0.01171552
## AGE -0.03060051
## CAT 0.07924266
## F 0.02246168
## FDX 0.08020173
## GM 0.05326574
## HPQ -0.03539714
## KMB 0.25030240
## MEL 0.07031147
## NYT 0.15387826
## PG 0.24340218
## TRB 0.14003743
## TXN -0.03882144