m.barra.9003 <- read.csv("~/Homework 4/m-barra-9003.txt", sep="")
da <- m.barra.9003
colnames(da)=c("AGE","C","MWD","MER","DELL","HPQ","IBM","AA","CAT","PG")
round(apply(da,2,mean),2)
## AGE C MWD MER DELL HPQ IBM AA CAT PG
## 1.36 2.08 1.87 2.08 4.82 1.37 1.06 1.09 1.23 1.08
round(apply(da,2,sd),2)
## AGE C MWD MER DELL HPQ IBM AA CAT PG
## 10.18 9.60 11.17 10.39 16.44 11.78 9.47 9.49 8.71 6.75
rmean <- apply(da,2,mean)
R.rm <- da-rmean
fin <- c(rep(1,4),rep(0,6))
tech <- c(rep(0,4),rep(1,3),rep(0,3))
oth <- c(rep(0,7),rep(1,3))
ind.dum <- cbind(fin,tech,oth)
ind.dum
## fin tech oth
## [1,] 1 0 0
## [2,] 1 0 0
## [3,] 1 0 0
## [4,] 1 0 0
## [5,] 0 1 0
## [6,] 0 1 0
## [7,] 0 1 0
## [8,] 0 0 1
## [9,] 0 0 1
## [10,] 0 0 1
cov.R <- var(R.rm)
sd.R <- sqrt(diag(cov.R))
corr.R <- cov.R/outer(sd.R,sd.R)
print(corr.R)
## AGE C MWD MER DELL HPQ IBM
## AGE 1.0000000 0.6182390 0.6208611 0.6451494 0.2723693 0.2991778 0.2695563
## C 0.6182390 1.0000000 0.7103460 0.6606422 0.2419522 0.3653988 0.4029239
## MWD 0.6208611 0.7103460 1.0000000 0.7945580 0.2552549 0.4640199 0.3846602
## MER 0.6451494 0.6606422 0.7945580 1.0000000 0.2293966 0.4827828 0.3043901
## DELL 0.2723693 0.2419522 0.2552549 0.2293966 1.0000000 0.4380318 0.3428415
## HPQ 0.2991778 0.3653988 0.4640199 0.4827828 0.4380318 1.0000000 0.4563985
## IBM 0.2695563 0.4029239 0.3846602 0.3043901 0.3428415 0.4563985 1.0000000
## AA 0.2937591 0.3728564 0.4039369 0.3806761 0.3214074 0.5013769 0.4113859
## CAT 0.2852101 0.3667355 0.2807931 0.2968603 0.1264958 0.2376912 0.3263783
## PG 0.1466319 0.2721835 0.2358181 0.2257318 0.1024763 0.0632573 -0.0339072
## AA CAT PG
## AGE 0.2937591 0.2852101 0.1466319
## C 0.3728564 0.3667355 0.2721835
## MWD 0.4039369 0.2807931 0.2358181
## MER 0.3806761 0.2968603 0.2257318
## DELL 0.3214074 0.1264958 0.1024763
## HPQ 0.5013769 0.2376912 0.0632573
## IBM 0.4113859 0.3263783 -0.0339072
## AA 1.0000000 0.6122479 0.0359865
## CAT 0.6122479 1.0000000 0.1150381
## PG 0.0359865 0.1150381 1.0000000
F.hat.o <- solve(crossprod(ind.dum))%*%(t(ind.dum)%*%t(R.rm))
E.hat.o <- R.rm-ind.dum%*%F.hat.o
diagD.hat.o <- diag(apply(E.hat.o,2,var))
Dinv.hat <- solve(diagD.hat.o)
H1 <- t(ind.dum)%*%Dinv.hat%*%ind.dum
Hmtx <- solve(H1)%*%t(ind.dum)%*%Dinv.hat
F.hat.g <- Hmtx%*%t(R.rm)
par(mfrow=c(3,1))
plot(F.hat.g[1,],type="l")
plot(F.hat.g[2,],type="l")
plot(F.hat.g[3,],type="l")
F.hat.gt <- t(F.hat.g)
E.hat.g <- R.rm-ind.dum%*%F.hat.g
diagD.hat.g <- apply(E.hat.g,2,var)
t(Hmtx)
## fin tech oth
## [1,] 0.2175398 0.0000000 0.0000000
## [2,] 0.2984352 0.0000000 0.0000000
## [3,] 0.2157723 0.0000000 0.0000000
## [4,] 0.2682527 0.0000000 0.0000000
## [5,] 0.0000000 0.2569000 0.0000000
## [6,] 0.0000000 0.3688383 0.0000000
## [7,] 0.0000000 0.3742617 0.0000000
## [8,] 0.0000000 0.0000000 0.1985438
## [9,] 0.0000000 0.0000000 0.3159870
## [10,] 0.0000000 0.0000000 0.4854691
cov.ind <- ind.dum%*%var(F.hat.gt)%*%t(ind.dum)+diag(diagD.hat.g)
sd.ind <- sqrt(diag(cov.ind))
corr.ind <- cov.ind/outer(sd.ind,sd.ind)
print(corr.ind)
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 1.00000000 0.3485291 0.31152914 0.3355360 0.13061681 0.15041073
## [2,] 0.34852911 1.0000000 0.34672390 0.3734429 0.14537315 0.16740327
## [3,] 0.31152914 0.3467239 1.00000000 0.3337980 0.12994028 0.14963168
## [4,] 0.33553596 0.3734429 0.33379804 1.0000000 0.13995364 0.16116248
## [5,] 0.13061681 0.1453731 0.12994028 0.1399536 1.00000000 0.26290444
## [6,] 0.15041073 0.1674033 0.14963168 0.1611625 0.26290444 1.00000000
## [7,] 0.15242249 0.1696423 0.15163301 0.1633180 0.26642080 0.30679471
## [8,] 0.09193178 0.1023177 0.09145561 0.0985033 0.05925829 0.06823840
## [9,] 0.11156333 0.1241671 0.11098549 0.1195382 0.07191259 0.08281036
## [10,] 0.13270135 0.1476932 0.13201402 0.1421872 0.08553794 0.09850052
## [,7] [,8] [,9] [,10]
## [1,] 0.15242249 0.09193178 0.11156333 0.13270135
## [2,] 0.16964230 0.10231770 0.12416711 0.14769318
## [3,] 0.15163301 0.09145561 0.11098549 0.13201402
## [4,] 0.16331804 0.09850330 0.11953817 0.14218718
## [5,] 0.26642080 0.05925829 0.07191259 0.08553794
## [6,] 0.30679471 0.06823840 0.08281036 0.09850052
## [7,] 1.00000000 0.06915109 0.08391795 0.09981797
## [8,] 0.06915109 1.00000000 0.14819746 0.17627657
## [9,] 0.08391795 0.14819746 1.00000000 0.21391953
## [10,] 0.09981797 0.17627657 0.21391953 1.00000000
glm(t(R.rm)[,1]~ind.dum[,1]+ind.dum[,2]+ind.dum[,3])
##
## Call: glm(formula = t(R.rm)[, 1] ~ ind.dum[, 1] + ind.dum[, 2] + ind.dum[,
## 3])
##
## Coefficients:
## (Intercept) ind.dum[, 1] ind.dum[, 2] ind.dum[, 3]
## -12.4908 -0.4257 4.8045 NA
##
## Degrees of Freedom: 9 Total (i.e. Null); 7 Residual
## Null Deviance: 375.3
## Residual Deviance: 321.5 AIC: 71.08