--- title: "Model Misspecification" output: html_document: toc: true --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) ``` ```{r, include=FALSE} # loading the packages and datasets dataff <- read.table("ffdata.csv", stringsAsFactors = F, header = T, sep = ",") # packages used listofpackages <- c("plm", "foreign", "lmtest", "aod", "systemfit", "dygraphs") for (j in listofpackages){ if(sum(installed.packages()[, 1] == j) == 0) { install.packages(j) } library(j, character.only = T) } ``` ## Model misspecification ```{r} # 1962:1 - 2014:6 regdata <- subset(dataff, select = c(exret_mkt, MOM, SMB, HML, PR15), date >= 1962 & date <=2014 + 5/12 + 0.01) capm <- lm(PR15 ~ exret_mkt, data = regdata) summary(capm) ffreg <- lm(PR15 ~ exret_mkt + HML + SMB, data = regdata) summary(ffreg) wald.test(vcov(ffreg), coef(ffreg), H0 = c(0, 0), df = 630-4, L = matrix(c(0, 0, 0, 0, 1, 0, 0, 1), ncol = 4), verbose = T) # now the same by hand ssr1 <- sum(capm$residuals^2) ssr2 <- sum(ffreg$residuals^2) Fstat <- (dim(regdata)[1] - 4)/2 * (ssr1 - ssr2)/ssr2 cat("F-statistic value:", Fstat) # the semi-partial R2 R2_hml <- (coef(ffreg)[3]/vcov(ffreg)[3, 3]^0.5)^2 * (1 - summary(ffreg)$r.squared) / (630 - 4) R2_smb <- (coef(ffreg)[4]/vcov(ffreg)[4, 4]^0.5)^2 * (1 - summary(ffreg)$r.squared) / (630 - 4) # other F-tests wald.test(vcov(ffreg), coef(ffreg), H0 = c(-0.4), df = 630-4, L = matrix(c(0, 0, 1, -1), ncol = 4), verbose = T) # restricted regression ## we now impose the restriction the regression with changed regressors regdata$newY <- regdata$PR15 - 0.4 * regdata$SMB regdata$newX <- regdata$HML + regdata$SMB ffreg2 <- lm(newY ~ exret_mkt + newX, data = regdata) summary(ffreg2) ## ... check that the restriction has been imposed correctly ffreg3 <- lm(newY ~ exret_mkt + HML+SMB, data = regdata) summary(ffreg3) # system of equations eq_asset <- PR15 ~ exret_mkt + SMB + HML eq_hml <- HML ~ exret_mkt eq_smb <- SMB ~ exret_mkt eq_system <- list(asset = eq_asset, HML = eq_hml, SMB = eq_smb) sfit <- systemfit(eq_system, data = regdata) summary(sfit) # now we want to conduct the analysis of the impact of a shock to the market excess return with bootstrapping # 1. output containers and bootstrap specification nrep <- 1000 impact <- array(0, c(12, 3, nrep)) # 2. the loop for (i in 1:nrep){ # preparing the bt errors temp <- sample(1:dim(regdata)[1], replace = T, size = 12) errors1 <- resid(sfit)[temp, ] temp <- sample(1:dim(regdata)[1], replace = T, size = 12) errors2 <- resid(sfit)[temp, ] # setting the shock last <- nrow(regdata) regdata2 <- regdata[(last - 11):last, ] regdataold <- regdata[(last - 11):last, ] regdata2$exret_mkt[1:6] <- regdata2$exret_mkt[1:6] + 0.05 # output impact[, , i] <- as.matrix((predict(sfit, newdata = regdata2) + errors1) - (predict(sfit, newdata = regdataold) + errors2)) } # 3. preparing the data for plotting imp_mean <- apply(X = impact, 1:2, mean) imp_up <- apply(X = impact, 1:2, function(x) {quantile(x, prob = 0.95)}) imp_down <- apply(X = impact, 1:2, function(x) {quantile(x, prob = 0.05)}) # 4. plotting (changing the variable, be caregul for plotting parameters ) variable <- 1 plotdata <- cbind(imp_mean[, variable], imp_up[, variable], imp_down[, variable]) colnames(plotdata) <- c("imp_mean", "imp_up", "imp_down") dygraph(ts(plotdata, start = c(2014, 1), frequency = 12), main = "Impact of a market shock") %>% dyAxis("y", label = "return", valueRange = c(-0.04, 0.09)) %>% dySeries(c("imp_up", "imp_mean", "imp_down"), label = "Impact") ```