R Markdown
#homwework08
#-----------------------------
library(ellipse)
## Warning: package 'ellipse' was built under R version 3.5.3
##
## Attaching package: 'ellipse'
## The following object is masked from 'package:graphics':
##
## pairs
#library(fEcofin)
#library(PerformanceAnalytics)
library(zoo)
## Warning: package 'zoo' was built under R version 3.5.3
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(readxl)
## Warning: package 'readxl' was built under R version 3.5.3
library(writexl)
## Warning: package 'writexl' was built under R version 3.5.3
?read_excel
## starting httpd help server ...
## done
fdata = read_excel("fama_french.xls")
fdata
## # A tibble: 360 x 9
## date `Mkt-RF` SMB HML RF ge ibm mobil
## <dttm> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1969-01-01 00:00:00 -1.2 -0.8 1.57 0.53 -1.20 -5.95 -1.40
## 2 1969-02-01 00:00:00 -5.82 -3.9 0.93 0.46 -6.04 -0.700 -7.84
## 3 1969-03-01 00:00:00 2.59 -0.28 -0.45 0.46 6.65 7.03 21.5
## 4 1969-04-01 00:00:00 1.52 -0.85 0.06 0.53 5.96 4.46 3.00
## 5 1969-05-01 00:00:00 0.02 -0.27 0.74 0.48 -3.58 -2.5 2.67
## 6 1969-06-01 00:00:00 -7.25 -5.31 -1.15 0.51 -3.82 5.88 -13.0
## 7 1969-07-01 00:00:00 -7.05 -3.27 1.36 0.53 -4.31 -3.92 -6.10
## 8 1969-08-01 00:00:00 4.65 0.89 -3.83 0.5 -2.76 6.63 10.8
## 9 1969-09-01 00:00:00 -2.88 1.2 -3.24 0.62 2.27 0.0725 -10.4
## 10 1969-10-01 00:00:00 4.96 3.78 -3.19 0.6 -1.03 4.42 -4.43
## # ... with 350 more rows, and 1 more variable: CRSP <dbl>
str(fdata)
## Classes 'tbl_df', 'tbl' and 'data.frame': 360 obs. of 9 variables:
## $ date : POSIXct, format: "1969-01-01" "1969-02-01" ...
## $ Mkt-RF: num -1.2 -5.82 2.59 1.52 0.02 -7.25 -7.05 4.65 -2.88 4.96 ...
## $ SMB : num -0.8 -3.9 -0.28 -0.85 -0.27 -5.31 -3.27 0.89 1.2 3.78 ...
## $ HML : num 1.57 0.93 -0.45 0.06 0.74 -1.15 1.36 -3.83 -3.24 -3.19 ...
## $ RF : num 0.53 0.46 0.46 0.53 0.48 0.51 0.53 0.5 0.62 0.6 ...
## $ ge : num -1.2 -6.04 6.65 5.96 -3.58 ...
## $ ibm : num -5.95 -0.7 7.03 4.46 -2.5 ...
## $ mobil : num -1.4 -7.84 21.51 3 2.67 ...
## $ CRSP : num -0.671 -5.364 3.05 2.053 0.504 ...
returns.mat = as.matrix(fdata[, c(-1,-2, -3, -4, -5, -9)])
market.mat = as.matrix(fdata[,2, drop=F])
n.obs = nrow(returns.mat)
X.mat = cbind(rep(1,n.obs),market.mat)
colnames(X.mat)[1] = "intercept"
XX.mat = crossprod(X.mat)
G.hat = solve(XX.mat)%*%crossprod(X.mat,returns.mat)
beta.hat = G.hat[2,]
E.hat = returns.mat - X.mat%*%G.hat
diagD.hat = diag(crossprod(E.hat)/(n.obs-2))
sumSquares = apply(returns.mat, 2, function(x) {sum( (x - mean(x))^2 )})
R.square = 1 - (n.obs-2)*diagD.hat/sumSquares
R.square
## ge ibm mobil
## 0.5824167 0.3039684 0.3263774
cbind(beta.hat, diagD.hat, R.square)
## beta.hat diagD.hat R.square
## ge 1.0580825 17.02494 0.5824167
## ibm 0.8149949 32.25874 0.3039684
## mobil 0.8158072 29.13458 0.3263774
par(mfrow=c(1,2))
barplot(beta.hat, horiz=T, main="Beta values", col="blue", cex.names = 0.75, las=1)
barplot(R.square, horiz=T, main="R-square values", col="blue", cex.names = 0.75, las=1)
par(mfrow=c(1,1))
cov.si = as.numeric(var(market.mat))*beta.hat%*%t(beta.hat) + diag(diagD.hat)
cor.si = cov2cor(cov.si)
rownames(cor.si) = colnames(cor.si)
ord <- order(cor.si[1,])
ordered.cor.si <- cor.si[ord, ord]
plotcorr(ordered.cor.si, col=cm.colors(11)[5*ordered.cor.si + 6])
w.gmin.si = solve(cov.si)%*%rep(1,nrow(cov.si))
w.gmin.si = w.gmin.si/sum(w.gmin.si)
colnames(w.gmin.si) = "single.index"
w.gmin.sample = solve(var(returns.mat))%*%rep(1,nrow(cov.si))
w.gmin.sample = w.gmin.sample/sum(w.gmin.sample)
colnames(w.gmin.sample) = "sample"
cbind(w.gmin.si, sample = w.gmin.sample)
## single.index sample
## ge 0.3312630 0.3204386
## ibm 0.3176039 0.2924230
## mobil 0.3511331 0.3871383
par(mfrow=c(2,1))
barplot(t(w.gmin.si), horiz=F, main="Single Index Weights", col="blue", cex.names = 0.75, las=2)
barplot(t(w.gmin.sample), horiz=F, main="Sample Weights", col="blue", cex.names = 0.75, las=2)
par(mfrow=c(1,1))
# compare means and sd values on global min variance portfolios
mu.vals = colMeans(returns.mat)
mu.gmin.si = as.numeric(crossprod(w.gmin.si, mu.vals))
sd.gmin.si = as.numeric(sqrt(t(w.gmin.si)%*%cov.si%*%w.gmin.si))
mu.gmin.sample = as.numeric(crossprod(w.gmin.sample, mu.vals))
sd.gmin.sample = as.numeric(sqrt(t(w.gmin.sample)%*%var(returns.mat)%*%w.gmin.sample))
cbind(mu.gmin.si,mu.gmin.sample, sd.gmin.si, sd.gmin.sample)
## mu.gmin.si mu.gmin.sample sd.gmin.si sd.gmin.sample
## [1,] 1.274778 1.283081 5.068253 4.948165
asset.names = colnames(returns.mat)
asset.names
## [1] "ge" "ibm" "mobil"