rm(list=ls())
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
library(tidyr)
library(tidyquant)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## ── Attaching core tidyquant packages ─────────────────────── tidyquant 1.0.11 ──
## ✔ PerformanceAnalytics 2.0.8 ✔ TTR 0.24.4
## ✔ quantmod 0.4.26 ✔ xts 0.14.1
## ── Conflicts ────────────────────────────────────────── tidyquant_conflicts() ──
## ✖ zoo::as.Date() masks base::as.Date()
## ✖ zoo::as.Date.numeric() masks base::as.Date.numeric()
## ✖ dplyr::filter() masks stats::filter()
## ✖ xts::first() masks dplyr::first()
## ✖ dplyr::lag() masks stats::lag()
## ✖ xts::last() masks dplyr::last()
## ✖ PerformanceAnalytics::legend() masks graphics::legend()
## ✖ quantmod::summary() masks base::summary()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(readr)
data <- read_csv("capm.csv")
## Rows: 2363 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): Date
## dbl (5): Close-tbill, Close-sp500, Close-msft, Close-ge, Close-ford
##
## ℹ 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(data)
## # A tibble: 6 × 6
## Date `Close-tbill` `Close-sp500` `Close-msft` `Close-ge` `Close-ford`
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1993/11/1 3.06 469. 2.52 6.72 8.16
## 2 1993/11/2 3.12 468. 2.5 6.69 8.31
## 3 1993/11/3 3.08 463. 2.45 6.64 8.24
## 4 1993/11/4 3.07 457. 2.38 6.46 7.98
## 5 1993/11/5 3.07 460. 2.45 6.49 7.94
## 6 1993/11/8 3.06 460. 2.45 6.46 8.02
dat <- read_csv("capm.csv") %>%
mutate(Date = as.character(Date) %>% as.Date(., "%Y/%m/%d")) %>%
filter(Date>=as.Date("1993-11-01") & Date<= as.Date("1998-11-30")) %>%
rename(rf = `Close-tbill`,
sp500 = `Close-sp500`,
msft = `Close-msft`,
ge = `Close-ge`,
ford = `Close-ford`) %>%
mutate(rf = rf/(100*360))
## Rows: 2363 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): Date
## dbl (5): Close-tbill, Close-sp500, Close-msft, Close-ge, Close-ford
##
## ℹ 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(dat)
## # A tibble: 6 × 6
## Date rf sp500 msft ge ford
## <date> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1993-11-01 0.000085 469. 2.52 6.72 8.16
## 2 1993-11-02 0.0000867 468. 2.5 6.69 8.31
## 3 1993-11-03 0.0000856 463. 2.45 6.64 8.24
## 4 1993-11-04 0.0000853 457. 2.38 6.46 7.98
## 5 1993-11-05 0.0000853 460. 2.45 6.49 7.94
## 6 1993-11-08 0.000085 460. 2.45 6.46 8.02
glimpse(dat)
## Rows: 1,278
## Columns: 6
## $ Date <date> 1993-11-01, 1993-11-02, 1993-11-03, 1993-11-04, 1993-11-05, 199…
## $ rf <dbl> 8.500000e-05, 8.666667e-05, 8.555556e-05, 8.527778e-05, 8.527778…
## $ sp500 <dbl> 469.10, 468.44, 463.02, 457.49, 459.57, 460.21, 460.33, 463.72, …
## $ msft <dbl> 2.52, 2.50, 2.45, 2.38, 2.45, 2.45, 2.44, 2.53, 2.55, 2.51, 2.56…
## $ ge <dbl> 6.72, 6.69, 6.64, 6.46, 6.49, 6.46, 6.51, 6.52, 6.47, 6.50, 6.65…
## $ ford <dbl> 8.16, 8.31, 8.24, 7.98, 7.94, 8.02, 7.94, 8.05, 8.07, 8.10, 8.10…
tail(dat)
## # A tibble: 6 × 6
## Date rf sp500 msft ge ford
## <date> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1998-11-20 0.000121 1164. 28.3 28.6 25.9
## 2 1998-11-23 0.000124 1188. 29.7 29.2 26.1
## 3 1998-11-24 0.000125 1183. 30.3 29.1 25.9
## 4 1998-11-25 0.000124 1187. 31.0 28.8 25.7
## 5 1998-11-27 0.000123 1192. 31.9 28.5 25.8
## 6 1998-11-30 0.000123 1164. 30.4 27.9 25.8
ret4 <- dat %>% select(-rf) %>%
gather(key = stock, value = price, -Date) %>%
group_by(stock) %>%
tq_transmute( mutate_fun = periodReturn,
period = "daily",
type = "arithmetic",
col_rename = "daily.returns") %>%
ungroup() %>%
spread(stock, daily.returns) %>%
bind_cols(., rf = dat$rf) %>%
# subract each returns by risk-free rate
mutate(ford_rf = ford - rf,
ge_rf = ge - rf,
msft_rf = msft - rf,
sp500_rf = sp500 - rf) %>%
# Delete the first row with 0 data
slice(-1) %>%
select(Date, ends_with("_rf"))
head(ret4)
## # A tibble: 6 × 5
## Date ford_rf ge_rf msft_rf sp500_rf
## <date> <dbl> <dbl> <dbl> <dbl>
## 1 1993-11-02 0.0183 -0.00455 -0.00802 -0.00149
## 2 1993-11-03 -0.00851 -0.00756 -0.0201 -0.0117
## 3 1993-11-04 -0.0316 -0.0272 -0.0287 -0.0120
## 4 1993-11-05 -0.00510 0.00456 0.0293 0.00446
## 5 1993-11-08 0.00999 -0.00471 -0.000085 0.00131
## 6 1993-11-09 -0.0101 0.00765 -0.00417 0.000174
ret4.reg <- ret4 %>% lm(formula = cbind(msft_rf, ge_rf, ford_rf) ~ sp500_rf, data = .)
head(cbind(ret4$msft_rf, ret4$ge_rf, ret4$ford_rf))
## [,1] [,2] [,3]
## [1,] -0.008023175 -0.004550952 0.018295686
## [2,] -0.020085556 -0.007559397 -0.008509142
## [3,] -0.028656706 -0.027193712 -0.031638676
## [4,] 0.029326487 0.004558685 -0.005097809
## [5,] -0.000085000 -0.004707496 0.009990567
## [6,] -0.004168022 0.007653549 -0.010061451
b_hat <- ret4.reg$coefficients
b_hat
## msft_rf ge_rf ford_rf
## (Intercept) 0.001233949 0.000344806 0.0002995505
## sp500_rf 1.282382154 1.198199332 1.0250994574
diagD_hat <- ret4.reg$residuals %>% cov() %>% diag() %>% diag(nrow = length(.))
diagD_hat
## [,1] [,2] [,3]
## [1,] 0.0002861862 0.000000e+00 0.0000000000
## [2,] 0.0000000000 9.148714e-05 0.0000000000
## [3,] 0.0000000000 0.000000e+00 0.0002484053
sigm2 = var(ret4$sp500_rf)
omega = sigm2*t(b_hat)%*%b_hat + diagD_hat
omega
## msft_rf ge_rf ford_rf
## msft_rf 0.0004170745 1.222959e-04 1.046283e-04
## ge_rf 0.0001222959 2.057549e-04 9.775984e-05
## ford_rf 0.0001046283 9.775984e-05 3.320421e-04