rp_measure.log.R (11322B)
1 2 R version 4.1.3 (2022-03-10) -- "One Push-Up" 3 Copyright (C) 2022 The R Foundation for Statistical Computing 4 Platform: x86_64-apple-darwin17.0 (64-bit) 5 6 R is free software and comes with ABSOLUTELY NO WARRANTY. 7 You are welcome to redistribute it under certain conditions. 8 Type 'license()' or 'licence()' for distribution details. 9 10 Natural language support but running in an English locale 11 12 R is a collaborative project with many contributors. 13 Type 'contributors()' for more information and 14 'citation()' on how to cite R or R packages in publications. 15 16 Type 'demo()' for some demos, 'help()' for on-line help, or 17 'help.start()' for an HTML browser interface to help. 18 Type 'q()' to quit R. 19 20 > #!/usr/bin/env Rscript 21 > # 22 > # rp_measure.R 23 > # 24 > # This code runs the predictive regression 25 > # We save in output the expected excess return estimate 26 > # 27 > # 28 > # (c) Valentin Haddad, Erik Loualiche & Matthew Plosser 29 > # 30 > # Last updated on June 4th 2019 31 > # 32 > ################################################################################## 33 > 34 > ################################################################################## 35 > message("Log file for code executed at\n") 36 Log file for code executed at 37 38 > message(format(Sys.time(), "%a %b %d %X %Y")) 39 Mon Jun 06 12:38:03 2022 40 > ################################################################################## 41 > 42 > 43 > ################################################################################## 44 > # APPEND REQUIRED PACKAGES 45 > 46 > # See this https://stackoverflow.com/questions/4090169/elegant-way-to-check-for-missing-packages-and-install-them 47 > using<-function(...) { 48 + libs<-unlist(list(...)) 49 + req<-unlist(lapply(libs,require,character.only=TRUE)) 50 + need<-libs[req==FALSE] 51 + if(length(need)>0){ 52 + install.packages(need) 53 + lapply(need,require,character.only=TRUE) 54 + } 55 + } 56 > 57 > package_to_load <- c("crayon", "devtools", "wesanderson", "ggplot2", "statar", 58 + "stringr", "lubridate", "lmtest", "sandwich", "stargazer", "data.table") 59 > 60 > using(package_to_load) 61 Loading required package: crayon 62 Loading required package: devtools 63 Loading required package: usethis 64 Loading required package: wesanderson 65 Loading required package: ggplot2 66 67 Attaching package: ‘ggplot2’ 68 69 The following object is masked from ‘package:crayon’: 70 71 %+% 72 73 Loading required package: statar 74 Loading required package: stringr 75 Loading required package: lubridate 76 77 Attaching package: ‘lubridate’ 78 79 The following objects are masked from ‘package:base’: 80 81 date, intersect, setdiff, union 82 83 Loading required package: lmtest 84 Loading required package: zoo 85 86 Attaching package: ‘zoo’ 87 88 The following objects are masked from ‘package:base’: 89 90 as.Date, as.Date.numeric 91 92 93 Attaching package: ‘lmtest’ 94 95 The following object is masked from ‘package:crayon’: 96 97 reset 98 99 Loading required package: sandwich 100 Loading required package: stargazer 101 102 Please cite as: 103 104 Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables. 105 R package version 5.2.3. https://CRAN.R-project.org/package=stargazer 106 107 Loading required package: data.table 108 109 Attaching package: ‘data.table’ 110 111 The following objects are masked from ‘package:lubridate’: 112 113 hour, isoweek, mday, minute, month, quarter, second, wday, week, 114 yday, year 115 116 > 117 > 118 > check_file = file.exists("log/R-session-info.log.R") 119 > sink("log/R-session-info.log.R", append=check_file) 120 > cat(bold("\n\n# -----\n# Session info for rp_measure.csv\n\n")) 121 > session_info() 122 > sink() 123 > ################################################################################## 124 > 125 > 126 > ################################################################################## 127 > dt_predict <- fread("./tmp/predict.csv") 128 > dt_predict[, datem := as.monthly(ISOdate(str_sub(dateym, 1, 4), str_sub(dateym, 5, 6), 1)) ] 129 > ################################################################################## 130 > 131 > 132 > ################################################################################## 133 > r1_3 = lm(rmrf_y3 ~ dp + cay + rf, data = dt_predict[ year(datem) < 2011]) 134 > nw_r1_3 = coeftest(r1_3, df = Inf, vcov = NeweyWest(r1_3, lag = 12, prewhite = FALSE) ) 135 > r2_3 = lm(rmrf_y3 ~ dp + cay + rf, data = dt_predict) 136 > nw_r2_3 = coeftest(r2_3, df = Inf, vcov = NeweyWest(r2_3, lag = 12, prewhite = FALSE) ) 137 > 138 > stargazer(r1_3, nw_r1_3, r2_3, nw_r2_3, type="text") 139 140 =========================================================================================== 141 Dependent variable: 142 ----------------------------------------------------------------------- 143 rmrf_y3 rmrf_y3 144 OLS coefficient OLS coefficient 145 test test 146 (1) (2) (3) (4) 147 ------------------------------------------------------------------------------------------- 148 dp 3.632*** 3.632*** 3.370*** 3.370*** 149 (0.440) (1.074) (0.436) (1.135) 150 151 cay 2.156*** 2.156*** 1.814*** 1.814*** 152 (0.261) (0.604) (0.246) (0.630) 153 154 rf -0.938*** -0.938** -1.246*** -1.246*** 155 (0.177) (0.364) (0.163) (0.378) 156 157 Constant -0.019 -0.019 0.011 0.011 158 (0.015) (0.037) (0.013) (0.036) 159 160 ------------------------------------------------------------------------------------------- 161 Observations 236 264 162 R2 0.384 0.344 163 Adjusted R2 0.376 0.336 164 Residual Std. Error 0.074 (df = 232) 0.074 (df = 260) 165 F Statistic 48.181*** (df = 3; 232) 45.459*** (df = 3; 260) 166 =========================================================================================== 167 Note: *p<0.1; **p<0.05; ***p<0.01 168 > 169 > # --- output for the readme 170 > star = stargazer(r2_3, type="text", style = "aer", 171 + covariate.labels = c("D/P ratio", "cay", "T-bill (three-month)"), 172 + dep.var.labels = "Future Excess Returns", 173 + omit.stat = c("ser", "adj.rsq") ) 174 175 =========================================================== 176 Future Excess Returns 177 ----------------------------------------------------------- 178 D/P ratio 3.370*** 179 (0.436) 180 181 cay 1.814*** 182 (0.246) 183 184 T-bill (three-month) -1.246*** 185 (0.163) 186 187 Constant 0.011 188 (0.013) 189 190 Observations 264 191 R2 0.344 192 F Statistic 45.459*** (df = 3; 260) 193 ----------------------------------------------------------- 194 Notes: ***Significant at the 1 percent level. 195 **Significant at the 5 percent level. 196 *Significant at the 10 percent level. 197 > 198 > star[1] = "~~~R" 199 > star[length(star)+1] = "~~~" 200 > cat(star, sep = '\n', file = './tmp/reg_update.txt') 201 > ################################################################################## 202 > 203 > 204 > ################################################################################## 205 > # OUTPUT PREDICTED VALUE 206 > dt_exp_rmrf <- cbind(dt_predict[!is.na(rmrf_y3), -c("datem")], exp_rmrf = predict(r2_3)) 207 > 208 > fwrite(dt_exp_rmrf, "./output/predict.csv") 209 > ################################################################################## 210 > 211 > 212 > ################################################################################## 213 > # PLOT 214 > dt_plot <- dt_exp_rmrf[, .( 215 + date=as.Date(ISOdate(str_sub(dateym,1, 4), as.integer(str_sub(dateym, 5, 6)), 1)), 216 + dp, cay, rf, rmrf_y3, exp_rmrf)] 217 > dt_plot[] 218 date dp cay rf rmrf_y3 exp_rmrf 219 1: 1952-03-01 0.05817138 0.01510493 0.0159 0.18092953 0.21495611 220 2: 1952-06-01 0.05739649 0.02483727 0.0170 0.21642173 0.22862789 221 3: 1952-09-01 0.05709103 0.01484008 0.0171 0.23193277 0.20933933 222 4: 1952-12-01 0.05522191 0.02216598 0.0209 0.22202729 0.21159354 223 5: 1953-03-01 0.05455042 0.02152118 0.0201 0.26058722 0.20915790 224 --- 225 260: 2016-12-01 0.02284466 -0.01902802 0.0051 0.11899389 0.04744519 226 261: 2017-03-01 0.02185142 -0.02168661 0.0074 0.01294375 0.03640904 227 262: 2017-06-01 0.02147595 -0.02432134 0.0098 0.07365837 0.02737347 228 263: 2017-09-01 0.02132813 -0.02799587 0.0103 0.08875791 0.01958676 229 264: 2017-12-01 0.02075782 -0.02490222 0.0132 0.12023222 0.01966250 230 > 231 > 232 > p0 <- dt_plot[, .(date, dp, cay, rf, rmrf_y3) ] %>% 233 + melt(id.vars="date") %>% 234 + ggplot(aes(date, value, colour = variable)) + 235 + geom_line(alpha=0.75, size=0.25) + geom_point(shape=1, size = 1, alpha=0.5) + 236 + theme_bw() 237 > # p0 238 > 239 > p1 <- dt_plot[, .(date, exp_rmrf, rmrf_y3) ] %>% 240 + melt(id.vars="date") %>% 241 + ggplot(aes(date, 100*value, colour = variable)) + 242 + geom_line(alpha=0.75, size=0.25) + geom_point(shape=1, size = 1, alpha=0.5) + 243 + xlab("") + ylab("Returns (percent)") + 244 + theme_bw() + 245 + theme(legend.position = c(0.3, 0.9)) + 246 + scale_colour_manual(name = "", 247 + breaks = c("exp_rmrf", "rmrf_y3"), 248 + values = c(wes_palette("Zissou1")[1], wes_palette("Zissou1")[5]), 249 + labels=c("Expected", "Realized")) + 250 + guides(colour = guide_legend(nrow = 1)) 251 > ggsave("./output/predict.png", p1, width = 8, height=6) 252 > ################################################################################## 253 > 254 > 255 > ################################################################################## 256 > 257 > proc.time() 258 user system elapsed 259 2.456 0.224 2.887