Note that the directories used to store data are likely different on your computer, and such references will need to be changed before using any such code.
library(tidyverse)
df <- read.csv("../../Data/Session_2-1.csv", stringsAsFactors=FALSE)
df_full <- df
uol <- filter(df, isin == "SG1S83002349")
#clean_df <- subset(df,fyear==2017 & !is.na(revt) & !is.na(ni) & revt > 1)
# revt: Revenue, at: Assets
summary(uol[,c("revt", "at")])
revt at
Min. : 94.78 Min. : 1218
1st Qu.: 193.41 1st Qu.: 3044
Median : 427.44 Median : 3478
Mean : 666.38 Mean : 5534
3rd Qu.:1058.61 3rd Qu.: 7939
Max. :2103.15 Max. :19623
mod1 <- lm(revt ~ at, data = uol)
summary(mod1)
Call:
lm(formula = revt ~ at, data = uol)
Residuals:
Min 1Q Median 3Q Max
-295.01 -101.29 -41.09 47.17 926.29
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -13.831399 67.491305 -0.205 0.839
at 0.122914 0.009678 12.701 6.7e-13 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 221.2 on 27 degrees of freedom
Multiple R-squared: 0.8566, Adjusted R-squared: 0.8513
F-statistic: 161.3 on 1 and 27 DF, p-value: 6.699e-13
# Graph showing squared error
uolg <- uol[,c("at","revt")]
uolg$resid <- mod1$residuals
uolg$xleft <- ifelse(uolg$resid < 0,uolg$at,uolg$at - uolg$resid)
uolg$xright <- ifelse(uolg$resid < 0,uolg$at - uolg$resid, uol$at)
uolg$ytop <- ifelse(uolg$resid < 0,uolg$revt - uolg$resid,uol$revt)
uolg$ybottom <- ifelse(uolg$resid < 0,uolg$revt, uolg$revt - uolg$resid)
uolg$point <- TRUE
uolg2 <- uolg
uolg2$point <- FALSE
uolg2$at <- ifelse(uolg$resid < 0,uolg2$xright,uolg2$xleft)
uolg2$revt <- ifelse(uolg$resid < 0,uolg2$ytop,uolg2$ybottom)
uolg <- rbind(uolg, uolg2)
uolg %>% ggplot(aes(y=revt, x=at)) +
geom_point(aes(shape=point, group=point)) +
scale_shape_manual(values=c(NA,18)) +
geom_smooth(method="lm", se=FALSE) +
geom_errorbarh(aes(xmax=xright, xmin = xleft)) +
geom_errorbar(aes(ymax=ytop, ymin = ybottom)) +
theme(legend.position="none") + xlim(0, 20000) + ylim(0, 2500)
`geom_smooth()` using formula 'y ~ x'
```r
# tidyverse
uol <- uol %>%
mutate(revt_growth1 = revt / lag(revt) - 1)
# R way
uol$revt_growth2 = uol$revt / c(NA, uol$revt[-length(uol$revt)]) - 1
identical(uol$revt_growth1, uol$revt_growth2)
<!-- rnb-source-end -->
<!-- rnb-output-begin eyJkYXRhIjoiWzFdIFRSVUVcbiJ9 -->
[1] TRUE
<!-- rnb-output-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyBNYWtlIHRoZSBvdGhlciBuZWVkZWQgY2hhbmdlXG51b2wgPC0gdW9sICU+JVxuICBtdXRhdGUoYXRfZ3Jvd3RoID0gYXQgLyBsYWcoYXQpIC0gMSkgJT4lICAjIENhbGN1bGF0ZSBhc3NldCBncm93dGhcbiAgcmVuYW1lKHJldnRfZ3Jvd3RoID0gcmV2dF9ncm93dGgxKSAgICAgICAgIyBSZW5hbWUgZm9yIHJlYWRhYmlsaXR5XG4jIFJ1biB0aGUgT0xTIG1vZGVsXG5tb2QyIDwtIGxtKHJldnRfZ3Jvd3RoIH4gYXRfZ3Jvd3RoLCBkYXRhID0gdW9sKVxuc3VtbWFyeShtb2QyKVxuYGBgXG5gYGAifQ== -->
```r
```r
# Make the other needed change
uol <- uol %>%
mutate(at_growth = at / lag(at) - 1) %>% # Calculate asset growth
rename(revt_growth = revt_growth1) # Rename for readability
# Run the OLS model
mod2 <- lm(revt_growth ~ at_growth, data = uol)
summary(mod2)
<!-- rnb-source-end -->
<!-- rnb-output-begin eyJkYXRhIjoiXG5DYWxsOlxubG0oZm9ybXVsYSA9IHJldnRfZ3Jvd3RoIH4gYXRfZ3Jvd3RoLCBkYXRhID0gdW9sKVxuXG5SZXNpZHVhbHM6XG4gICAgIE1pbiAgICAgICAxUSAgIE1lZGlhbiAgICAgICAzUSAgICAgIE1heCBcbi0wLjU3NzM2IC0wLjEwNTM0IC0wLjAwOTUzICAwLjE1MTMyICAwLjQyMjg0IFxuXG5Db2VmZmljaWVudHM6XG4gICAgICAgICAgICBFc3RpbWF0ZSBTdGQuIEVycm9yIHQgdmFsdWUgUHIoPnx0fCkgIFxuKEludGVyY2VwdCkgIDAuMDkwMjQgICAgMC4wNTYyMCAgIDEuNjA2ICAgMC4xMjA0ICBcbmF0X2dyb3d0aCAgICAwLjUzODIxICAgIDAuMjc3MTcgICAxLjk0MiAgIDAuMDYzMSAuXG4tLS1cblNpZ25pZi4gY29kZXM6ICAwIMOi4oKsy5wqKirDouKCrOKEoiAwLjAwMSDDouKCrMucKirDouKCrOKEoiAwLjAxIMOi4oKsy5wqw6LigqzihKIgMC4wNSDDouKCrMucLsOi4oKs4oSiIDAuMSDDouKCrMucIMOi4oKs4oSiIDFcblxuUmVzaWR1YWwgc3RhbmRhcmQgZXJyb3I6IDAuMjQ0NCBvbiAyNiBkZWdyZWVzIG9mIGZyZWVkb21cbiAgKDEgb2JzZXJ2YXRpb24gZGVsZXRlZCBkdWUgdG8gbWlzc2luZ25lc3MpXG5NdWx0aXBsZSBSLXNxdWFyZWQ6ICAwLjEyNjcsXHRBZGp1c3RlZCBSLXNxdWFyZWQ6ICAwLjA5MzA3IFxuRi1zdGF0aXN0aWM6IDMuNzcxIG9uIDEgYW5kIDI2IERGLCAgcC12YWx1ZTogMC4wNjMwN1xuIn0= -->
Call: lm(formula = revt_growth ~ at_growth, data = uol)
Residuals: Min 1Q Median 3Q Max -0.57736 -0.10534 -0.00953 0.15132 0.42284
Coefficients: Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.09024 0.05620 1.606 0.1204
at_growth 0.53821 0.27717 1.942 0.0631 . — Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2444 on 26 degrees of freedom (1 observation deleted due to missingness) Multiple R-squared: 0.1267, Adjusted R-squared: 0.09307 F-statistic: 3.771 on 1 and 26 DF, p-value: 0.06307
<!-- rnb-output-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyBsY3Q6IHNob3J0IHRlcm0gbGlhYmlsaXRpZXMsIGNoZTogY2FzaCBhbmQgZXF1aXZhbGVudHMsIGViaXQ6IEVCSVRcbnVvbCA8LSB1b2wgJT4lXG4gIG11dGF0ZV9hdCh2YXJzKGxjdCwgY2hlLCBlYml0KSwgbGlzdChncm93dGggPSB+LiAvIGxhZyguKSAtIDEpKSAgIyBDYWxjdWxhdGUgMyBncm93dGhzXG5tb2QzIDwtIGxtKHJldnRfZ3Jvd3RoIH4gbGN0X2dyb3d0aCArIGNoZV9ncm93dGggKyBlYml0X2dyb3d0aCwgZGF0YT11b2wpXG5zdW1tYXJ5KG1vZDMpXG5gYGBcbmBgYCJ9 -->
```r
```r
# lct: short term liabilities, che: cash and equivalents, ebit: EBIT
uol <- uol %>%
mutate_at(vars(lct, che, ebit), list(growth = ~. / lag(.) - 1)) # Calculate 3 growths
mod3 <- lm(revt_growth ~ lct_growth + che_growth + ebit_growth, data=uol)
summary(mod3)
<!-- rnb-source-end -->
<!-- rnb-output-begin eyJkYXRhIjoiXG5DYWxsOlxubG0oZm9ybXVsYSA9IHJldnRfZ3Jvd3RoIH4gbGN0X2dyb3d0aCArIGNoZV9ncm93dGggKyBlYml0X2dyb3d0aCwgXG4gICAgZGF0YSA9IHVvbClcblxuUmVzaWR1YWxzOlxuICAgICBNaW4gICAgICAgMVEgICBNZWRpYW4gICAgICAgM1EgICAgICBNYXggXG4tMC40NjUzMSAtMC4xNTA5NyAgMC4wMDIwNSAgMC4xNzYwMSAgMC4zMTk5NyBcblxuQ29lZmZpY2llbnRzOlxuICAgICAgICAgICAgRXN0aW1hdGUgU3RkLiBFcnJvciB0IHZhbHVlIFByKD58dHwpICAgXG4oSW50ZXJjZXB0KSAgMC4wNzQ5OCAgICAwLjA0OTE1ICAgMS41MjYgIDAuMTQwMTggICBcbmxjdF9ncm93dGggICAwLjIzNDgyICAgIDAuMDczMTkgICAzLjIwOSAgMC4wMDM3NiAqKlxuY2hlX2dyb3d0aCAgLTAuMTE1NjEgICAgMC4wOTIyNyAgLTEuMjUzICAwLjIyMjMwICAgXG5lYml0X2dyb3d0aCAgMC4wMzgwOCAgICAwLjAyMjA4ICAgMS43MjQgIDAuMDk3NTEgLiBcbi0tLVxuU2lnbmlmLiBjb2RlczogIDAgw6LigqzLnCoqKsOi4oKs4oSiIDAuMDAxIMOi4oKsy5wqKsOi4oKs4oSiIDAuMDEgw6LigqzLnCrDouKCrOKEoiAwLjA1IMOi4oKsy5wuw6LigqzihKIgMC4xIMOi4oKsy5wgw6LigqzihKIgMVxuXG5SZXNpZHVhbCBzdGFuZGFyZCBlcnJvcjogMC4yMjI4IG9uIDI0IGRlZ3JlZXMgb2YgZnJlZWRvbVxuICAoMSBvYnNlcnZhdGlvbiBkZWxldGVkIGR1ZSB0byBtaXNzaW5nbmVzcylcbk11bHRpcGxlIFItc3F1YXJlZDogICAwLjMzLFx0QWRqdXN0ZWQgUi1zcXVhcmVkOiAgMC4yNDYyIFxuRi1zdGF0aXN0aWM6ICAzLjk0IG9uIDMgYW5kIDI0IERGLCAgcC12YWx1ZTogMC4wMjAzM1xuIn0= -->
Call: lm(formula = revt_growth ~ lct_growth + che_growth + ebit_growth, data = uol)
Residuals: Min 1Q Median 3Q Max -0.46531 -0.15097 0.00205 0.17601 0.31997
Coefficients: Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.07498 0.04915 1.526 0.14018
lct_growth 0.23482 0.07319 3.209 0.00376 che_growth -0.11561 0.09227 -1.253 0.22230
ebit_growth 0.03808 0.02208 1.724 0.09751 . — Signif. codes: 0 ‘*’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2228 on 24 degrees of freedom (1 observation deleted due to missingness) Multiple R-squared: 0.33, Adjusted R-squared: 0.2462 F-statistic: 3.94 on 3 and 24 DF, p-value: 0.02033
<!-- rnb-output-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuZGV0ZWN0b3IgPC0gZnVuY3Rpb24oKSB7XG4gIGRpY2UgPC0gc2FtcGxlKDE6Niwgc2l6ZT0yLCByZXBsYWNlPVRSVUUpXG4gIGlmIChzdW0oZGljZSkgPT0gMTIpIHtcbiAgICBcXGV4cGxvZGVkXFxcbiAgfSBlbHNlIHtcbiAgICBcXHN0aWxsIHRoZXJlXFxcbiAgfVxufVxuXG5leHBlcmltZW50IDwtIHJlcGxpY2F0ZSgxMDAwLGRldGVjdG9yKCkpXG4jIHAgdmFsdWVcbnAgPC0gc3VtKGV4cGVyaW1lbnQgPT0gXFxzdGlsbCB0aGVyZVxcKSAvIDEwMDBcbmlmIChwIDwgMC4wNSkge1xuICBwYXN0ZShcXHAtdmFsdWU6IFxcLCBwLCBcXC0tIEZhaWwgdG8gcmVqZWN0IEhfQVxuYGBgIn0= -->
```r
```r
detector <- function() {
dice <- sample(1:6, size=2, replace=TRUE)
if (sum(dice) == 12) {
\exploded\
} else {
\still there\
}
}
experiment <- replicate(1000,detector())
# p value
p <- sum(experiment == \still there\) / 1000
if (p < 0.05) {
paste(\p-value: \, p, \-- Fail to reject H_A
[1] \p-value: 0.974 -- Reject H_A that sun exploded\
```r
library(tidyverse)
read_csv('../../Data/Session_3-1.csv') %>%
ggplot(aes(y=revtq, x=atq)) +
geom_point() +
geom_smooth(method=\lm\) +
xlab(\Assets\) +
ylab(\Revenue\)
<!-- rnb-source-end -->
<!-- rnb-plot-begin -->
<img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAArwAAAGwCAMAAAB8TkaXAAAAxlBMVEUAAAAAADoAAGYAOpAAZrYzMzMzZv86AAA6ADo6AGY6kNtNTU1NTW5NTY5NbqtNjshmAABmADpmZmZmtv9uTU1uTW5uTY5ubo5ubqtuq+SOTU2OTW6OTY6Obk2OyP+QOgCQ2/+rbk2rbm6rbo6ryKur5OSr5P+2ZgC2kDq225C22/+2///Ijk3I///KysrW1tbbkDrbkJDb/7bb/9vb///kq27k///r6+v/tmb/yI7/25D/27b/5Kv//7b//8j//9v//+T///84JZNcAAAACXBIWXMAAA7DAAAOwwHHb6hkAAAdgElEQVR4nO2dDZvbtnKFmS/n1vJ13PomjZN0nY/bddusP9peW153bS///5+qSIkiQAIgQAKcM9SZJ0+0osjDM+S78AgkZ6uawVAalbQBBmNuEF6G2iC8DLVBeBlqg/Ay1AbhZaiNHPDuAxH8MDI2pAFiQ3kqhFdEA8SG8lQIr4gGiA3lqRBeEQ0QG8pTIbwiGiA2lKdCeEU0QGwoT4XwimiA2FCeCuEV0QCxoTwVwiuiAWJDeSqEV0QDxIbyVAiviAaIDeWpEF4RDRAbylMhvCIaIDaUp0J4RTRAbChPhfCKaIDYUJ4K4RXRALEBmEpVVfGbEV4JDRAbeKlUVQK9hFdEA8QGXiqEF18DxAZeKoQXXwPEBmAqrHnhNUBsKE+F8IpogNhQngrhFdEAsaE8FcIrogFiQ3kqhFdEA8SG8lQIr4gGiA3lqRBeEQ0QG8pTIbwiGiA24FJJmeUlvEIaIDbQUkm6vkZ4hTRAbCCk0uBKeDVpgNgASKXllfBq0gCxAZCKBS9rXg0aIDYAUqmGoy3vKkPXALGBkMqA3pTKgfCKaIDYgEiF8CrTALEBkQrhVaYBYgMjlVHNG4sv4RXRALEBlsoJ2mh6Ca+IBogNrFQ6aAkvtgaIDaxUCK8ODRAbWKmcoWXNC60BYgMslbTra4RXSAPEhvJUCK+IBogN5akQXhENEBvKUyG8IhogNpSnkhNeBkMkOPKuqgFiQ3kqhFdEA8SGSCqOCTGvRnjyjPCKaIDYkEjFdf2sNjBNuMWM8IpogNjAgde8utZ/THgRNUBsEF7CKyKBopGn5nXDy5oXUQPEBlAqzpp3ajPCK6EBYkN5KoRXRAPEhngqp2GW8GrSALEhnUr3vFptLknYNeGV0ACxIZ1KdaK3thbE75rwSmiA2JBOhfBq1ACxIZ5KAN4IigmviAaIDYBUfDVvzBhMeEU0QGzgphLVeoTwimiA2MBIxUFpRXhxNUBsQKQyxLR5F9fzifCKaIDYgEilsofZqqrixl3CK6QBYqNAKqm9FwxajbdxMoRXRAPERv5U0iZqjW2Sphm6XRNeCQ0QGwjwtvO8g7rB/tS/a8IroQFiAwDedoPaB2lYjvCKaIDYAKh5j/C2W/mmzPy7JrwSGiA2AFKp7HB96t814ZXQALEhmIp1D4N/eow1L54GiI3F5M2WsEkdwBtbehBeEQ0QGznIywFvPWSXndGBNUBsCMN7VqktUcKLrQFiQw5em976tMQhHwrCK6IBYkOu5j2pmPC65UNBeEU0QGzIpuKFN3rXhFdCA8SGcCpWzeuOV8FdE14JDRAb6Km8ekV48TRAbGCnckD3wYPgZoRXQgPEhlwq9ncyl0ZDbhOhXRNeCQ0QG2KpDGbDxhoduiwb8DRAbMikMr6PYahxQjdc8RJeIQ0QGyKpOG7CsTUi0SW8QhogNmI0pq4YzIbXo9GhG6FEeEU0QGxEaExeq50Hr08jAV3CK6QBYkMEXsdY3mkkoUt4hTRAbMjA69GImmCwNyO8EhogNpDgTUaX8AppgNjAgTd2gsHejPBKaIDYQIH39Rx0Ca+QBogNkamyURyL3VRy94RXSAPEBkIqJ3RnsEt4ZTRAbMinkjg5Ntg14ZXQALEhnMqrDt2ZGoRXRAPERvFUggXzq3OxS3g1aYDYyJaKB9LQVMUro9glvJo0QGzkSsUHqR/eV1axS3g1aYDYkIJ3eCWY8GrSALFRGl53OTG+iYHwatIAsZGk4eEzVPO6wnWrOeHVpAFiI0XDN7pOSdhbuZ+SILyaNEBsrACvtZnvkgTh1aQBYmNdeP1X0wivJg0QGzNq3sAjPL6tjluELgSXhPfuh5fDRZ9+2j1+V9dvd7vdX7sPg/uZZ2+rGiA20jXG429UzTvxlERBeG97Pru4/+2qfvtdXd9cceQVksCAN26aYeopiXLw3jz6ezPyHsbaluFPv7w8/u8wHt//fk14hSTk4XU8w+6M6ackipcNN8/asfYI792P7+pPP18fgN7t2sH3m0OEJBigcYBv5vrn5iHhDY5PSbx+Pc9enKfQhw28B1Qbbm93TTy7fdzCe/f9dd2PvsFfknm/W1vVALFh/93UxHB0vhlH3FMSpUfedpB9dG2PvO2n57q3gL2taoDYSIXXWjeK3ci+TdEO7M2i4f3l9KXNrHkJr5yEALyDlY9vAjain5JYo+ZtioVj3P/WVsDNgvs/OFUmIZFJI3HgHa/ttZHQh6E4vIe64VE/t9DP8/YLC9jbqgaIjUSNFHiTWojwCpsmDRAbqRquYdopkdiHgfBq0gCxUSiV5BYihFeTBoiNEqnMaSFCeDVpgNjIn8q8FiKEV5MGiI3cqcxtIUJ4NWmA2EifbQhJzO9+Q3g1aYDYiLgZ134XmG1Y0riJ8GrSALExoTG+ouaFd1HPMcKrSgPERi54T+gWsuHfjPBKaIDYSIPXU/NG/J3VZTb8mxFeCQ0QG2ON9OckXi9Gl/Dq0gCxMdJIvr13UanrtRG7GeGV0ACxsRTePOgSXl0aIDZi4Q21y3tdwEbsZoRXQgPExlTN2y90LD5dCBZMhfCKaIDYiNRwwXuuGAhvnNftaIDYmA2vcas54Y3zuh0NEBuxGgN2rackCG+c1+1ogNiYpTG41ZzwxnndjsZKNibnvdJ9jJ6SILxxXrejsY6N6VnbRB+upyQIb5zX7WiohNf9lAThjfO6HQ2F8Po6PRLeOK/b0VBX8/qblBLeOK/b0QCxEavxIPAwO+GN87odDRAbcRohdAlvrNftaIDYsDXcRcZUBxHCG+d1OxogNiwN5x0M0x1ECG+c1+1ogNiYgjemgwjhjfO6HQ0xGwM+g/DG3WpOeOO8bkdDysYQ0EDNG/uUBOGN87odDRx4PVPB8e3yCG+c1+1owMDruQiX0OmR8MZ53Y4GTM3rhDephQjhjfO6HQ0QG054E/swEN44r9vRmCcR+LI124f7AZ8kiRw25m1GeCU05v6Tb5JWIpUZfRgIb5zX7WhgwjurhQjhjfO6HQ1AeOf8LYkCNpI2I7wSGjg1bxcz0SW8sV63owFi46yxoL8u4Y3zuh0NEBsnjUVNSglvnNftaIDYaDUW9tclvHFet6Oxqo3Ak2z1vBmGWTbyaxBeEY01bQSeIc7RX5fwxnndjgYCvMv+gE+yjQIahFdEAwBe4a7mOTTc8L6pqqdvvvoH4S2lIV7znubGtKUy2MwF74uv/ufJ08/Pvya8pTSEbZyndXWn4oL345Onh//q91/+SXgLaYjaOBa7uXwQ3jiv29GQtPHKmNbVnYoL3vpNUzZ8fPIwjV3Cu6ZE4t9LO4d9IRgilcxf2N43X1BT2SW8K0ok/7W/YwyvpiGkMl/DDe+8KGBvqxoy8I4vBCOkMl+D8IpoFIF3AmfXvC5CKvM1XPB+fNIemopf2IppWBKzCgAnuwGhYp2hweA9IfzPv3LkLaVhSswsX0c2Qvcw+G413yq89fvkS2yMOdFCV1TndYPu6wy7AA0XvCwbimmUGHm9nZsC9+tuduR9kTryFrC3VY3lNW+sjfCt5tuD9/SF7QvWvMU01rIx9ZSEolRcm/lH3uQoYG+rGuvYmL7pUU0q7s0Ir4TGCjai7tfVkYp3Mxe8H77lPG9ZjdI2Yp+SUJBKaDMHvOl38hLe9SVCGtFPSeCnEtzMAW9zPyThLapR1EZCBxH0VCY2c8D7+TnhLazhl5icNzuv4NFIasOwPXjTL08Q3lwSk1cs+hWcGut3hgaDlzfmFNcoBG/yA8Hbg3duFLC3VY1UePuFAXhlOkMT3jiv29FIrHlNpL01r1BnaDh42behsEaihHM8tjXEOkOjwcu+DaU1csM7t9Hj9uDlo+/FNVIlXMVEryHaGZrwxnndjkZGG69ezUYXLZXkzRzwsm9Dfo3hX57MZmMJupuEl30bcmuM/2B1HhsL/pZELh9w8M6LAva2olEG3oUN+fP4ILxxXvVqlIB3ObpbhDe92iW8UxFR86Y9ypaBXI+P9SXyfmE7HMcZN5YVsLdVDYdE0kPEedDdJLyHeMEHMNM00sbNZfBm62q+UXgbfjnPG6+R2HxhCbwZu5pvFN4X6XdEEt7FEhFb9gUDyOFAg3dOzUB4l8IbEdYEA8jhAIOXz7DN0Ihnt1mzHryPisHcGMjhAIN3bhSwt0GNYz+8wfuI7Ua3msunkksiL7y8n7ecxjx4HU9JyKeSSyJvzcv7ectpzILXdau5fCq5JDLXvLwlspxGes3rfkoCIJVMEoRXlUaKhO9qmsJUMmu44OX9vMU14iX8F4LVpZJdwwkv7+ctrREpEXxKQlcqJTTc8M6LAva2qhEnEX5KQlUqRTRC8P43a95SGjESU09JKEqlkIYD3hdV1UySfXzCL2zFNKYlpm81V5NKMY0xvG++/PPz84dN2Zt6kbiAva1qTEnEPCWhJJWCGiN4G3Dr91/9R5X+V9gK2NuqRlgi7lZzFakU1RjB296V8+Hb9LkGwptJIvYpCQWpFNbwwTvnKbYC9raq4ZdYtav5NuH9y4zbeS8N3vS//ddv4bOxcldzwnuh8AZvqJlqVOq2sXpXc8JLeEfh/iwMb3IbBpDDgQVv1QXneUMaeeGd0UEE5HAgwbsgCthD1ghXDaMPK3PhwMasNgwgh4PwxnnVo+Fh1y0xs4MIyOEgvHFeNWv44Z3b/AbkcBDeOK+aNXzwzu/bBHI4CG+cV9UaViXRSSzpOQZyOAhvnNftaLQSC3uUIqUio0F4l2s4ph2mrr7VGdrrghwOwhvnFVPDMS82+Th7jdEYmvAS3mR4s7TXBTkchDfOK6aGQWr3wwS8eTpDgxwOwhvnVUwj/AclTHZ7er1iDbq6D0duCcJbUGMwkHo1Yjo3QXU1J7yE17OeI8C6mm8f3re73e6vLwcLP/20e/yu+eH+tyvCa6wY0oHrar59eG+uxssaZN9+d0R7+/BO1LyRAdjVfPPw3v9+3Y217fj76ZeXx//d/XD44e5f/+0C4M2gAdnVfPPwHqDdNaPrzbPjWNvCe/fju/rTz9cHsv/9WDZ8cwj/4H3pcbyaJu1is+GH9+77htHrBtUDt7e7Jp7dPj7C+/bZJdS8CzUcF4K1plJGovBsw81VOwA/urZH3sPLZcMb2dR8dEkCMBVBifLw/nKacDBr3rfHcfhi4Y2fGytqQ1QDGt6mQrj/42VT894eZ8fqZrbhVAFfxFSZV2MSXt+FYLxUJCXKzvO25cJP7cspLmqe16sx+x4GvFQkJXiFTURj6h4G8K7mhPeC4D2yGqMxcau5fCrZNAhvnFdpjVOVYJUNzhUnn5IQTyWfBuGN8yqtMYTXU/LqaQxNeAmvFVFPSYinkk+D8MZ5FdcY1LwueOOekpBPJZsG4Y3zCqcxYjf2AR+8VCQlCK+8RsLD7OiprCtBeKU1kvowYKeytgThLa6RszP0Bg5HRgnCW1qj/3Y20khvfqP/cOSUILylNbzwzukgov9w5JQgvBG31i7y4YN3VgcREPAILwa8MT0Tlvlw1rwzm9+AgEd4LwVeh8bsnmMg4BHeS4F3NPIuaJcHAh7hxYA3ruaNK4x9m9o175J2eSDgEV4QeGMicnie2Lb1sazTI8bhILwXCe9GupoT3guBty85FjfkRzkchFcTvEk1r2fdDOjCHA7CqwneBA33KJ2lIb/Gw1FSgvDm1vDdaZ6lM7S+w1FSgvDm1nDAm68ztL7DUVKC8GbXGLKbszO0wsNRUILwltUwv6bpPuO5NQhvnFfXwsS5r3k+7BkG3Wc8twbhjfPqWJY6czvDx2huTPcZz61BeOO8OpYVh7dMZ2gQ8JSnsgl4g/Qu+mMopTpDg4CnPJXNw+v4M1TxY3WxztAg4ClP5QLhja40CnaGBgFPeSqbgDe00Wx4/ReCdZ/x3BqEN86ra+EkiKOaNwrewp2hQcBTnop6eKOip7WOZ7eAj5wSKBqENxwn2mI0nGAavNYx7JbvDA0CnvJUNMDb8Rah4UYzDd7Jux51n/HcGoQ3GBa8HXpuBCvn7EMCvDG3mus+47k1CG8oKhPe7mcPg5WX3rOPELtxT0noPuO5NQhvKM40RsDrGXrjfMQ+JaH7jOfWILyhSIJ36qqF38eqnaFBwFOeigJ4z//QR9S8k7O4Xh/rdoYGAU95KhrgTdOYNfImPVap+4zn1iC8vhjdVDO94hhe871TY/225iDgKU8FG14bxICGseIIXmvBWEOkrTkIeMpTAYZ3NO0VB++oGg7CO6eFiO4znluD8Dqimgtv8DNbY14LEd1nPLcG4XXE+HrDaapssFK/8njx+DPTx9zuN7rPeG4NwuuI8WBauyA9vrU5t9fywDu/cZPuM55bg/C6YlQIeOBtFtg1RgS8S3qO6T7juTUIb5xXH7xVD28/DltrDXwsa5en+4zn1iC8cV73rpp3EN1iv8biHqW6z3huDcIb53W4wB5zh9/v3BoYnaFBwFOeCiS8JyinNEbYTrObpb2u7jOeW4PwGuEfRcPwmgtPPwylcTpDg4CnPBU4eAMlQAjewTJve12e8dwahLePKHjNWYXRzIIH3uOoyzOeW2Mj8GYJg93wOsbK482NdY7RVgylPDNkA2bk9T9D2WvYw6vnNhzzg77Y5XCVW2MjI28uew54j+Oo9c5X2TrC+J7GM55bg/DaYc0ddK8GvMGvZaOwphh4xnNrEF5PmDVAB2/VX5nYm1/d+g0sns8VQxXbuGQydJ/x3BqE1xNWfXDUcE5GGChbnxmX005LecZzaxBeT9gs1uNF5nr7IbzWlWDCW0iD8PrCB69rtb01LFeDy2mEt5AG4W3C89VrAK/7O9qw5h2gW51vluAZz61BePfhHjgnjWpiguE8/Lbo9usam/CM59YgvHvf7O5oDVfZUFn3oJ/QNSYlCG9BDcK7H8DrHGMtdu2V+4r3NDfWCxDeshqEt4kRjj54+zXMtft64cFgA1OaZzy3xsXDOywE3PDWU/COnpIY1xg847k1Lh3eEWKjAbP9qTZYHBcZ++5b2lCH8BbVILyOoXe8gtWSf7jF+VtaWJhnPLcG4R0ytrfpdMA7iP5bmldk2kdk6D7juTUuHV7X9YlhXdDe29AXvPbKo4tp3mcxecZza1wyvD7Kzt/Z+ne19VEfrgvB1irGzzzjuTUuGF5XJdB9WesZrM7wjmchTsWuJWhtOL6BcmHoPuO5NQivY1GPnQFvZUa79ul72ui+nMFvwGCrRaH7jOfWILyeRR1vx//XRjFxWulYMLhRdaKLCu8MY4RXGl73DMEEvN0nXbHrZNXDriC8xs6HEnOcEV4ReLsT5ThhJ+qse20MeI1NHxjF7qoD79Sh9u3E3D3hzaEhAW93phxn7LTIyWF/P+++e0rCUhhC6iS33KE2fiPdHBLe3BqC8LpGQt9o2YPXLugGXVvB/pVwj7kZD7Vt/rwLw5NnDZdEOruEVwJeaziNgbc7sx28o6cqx9ouqSWHabAPx74d8Frm99bYXO8dEolBeOXgHbHrBLcyxtcjvH2pO6bfVHWym+FQd2I+eM1/ABxZjumeG4QXB14fu/bq51LXRncAjCmY5TAZMYTXZtadjVOA8C7VEIN3dOZi2D0/JTHA4vSjk5W5h8mnc97DiGET5rMV29BgyRJ2Ca8AvMORy1jqib1B7oPByq6NQz6icfErVcMB1+HVNjcQrczOVfMDRIPwhtg9rv3gfK/5mFQ/vEMAJ+GOWbM2h/iwbZ8oCHiEdx68XnYdOBi3mrs+9EMzYigRXve1lNraVTS5e8KbWWN9ePt/O2MH3m5ubHhjjgfe/XBXoQVTta2LxfnwnncGAh7hXQqv7+Rb7J40PKT4oDnvylyyd6wyjKBsbS03fjZepw77xOcxAaJxUfB6zrx3EBtM6/qYdbE7+nbvvp1rjNpQdQDv3vBu7cPpIXw45geIxiXBa1IxNZD2lyQm1qusCYB+T+Ye7YV7y49zk37bYc3r2NRc32XCdzgWBIjGxcE7rCg98WA8rzsZe2tHzhF7bzHVwz9g1vrBlYvrwzG7TnpBwCO8M+DtIIpAt6qmIXfSa5E6WGPAlLVit8T6xJvLxOElvCU15OC1ikYPuinQ+uH1fiXc743qdW8M1n242Es41IS3pMbq8A6+DJlLcqBrMDlF7/mDiSph0aH26IGAR3gT4d1HFLwBdCdLiPMOTJjPo6p1dcyswM23U8cs/TAXkEDRuDB4p/ALjbqTW4/2YO21h9dcaFmLOGbRR7cyi5N5EjlsFNUgvAN0A+xObG0WsTaZHbw2SVG0BnIJhmV5nkQOG2U1CG8cuiaTwVWcTLbLVjzUZzOEt5CGALyL0B1S6rndIfdhmqFhmCG8ZTSQ4J0xxeBWy3+YZmiYXkaeQMAjvInwZkS3OpcP/YWyMLtrwytvo7gG4X2wYGI3SEmmwzRHI+QKBDzCmwZvEXTj6eUZz61xQfD60J1kk/CCalwOvDPRdV5Mti4CE14pjQuGN7Je2Nsl5ADWeHZ5xrNrXCy8kaPushRdPqQ1QGwoT2VFeJ2j7hS7GVIc+kDQALGhPJX14HWRO/09LUOKiBogNpSnshq8M9DNlCKiBogN5amIwJs+6C5JEVEDxIbyVATg9aN7rIu9kwcgZ0v5Gc+tcQHw2uzGjrVZUkTUALGhPJW14Q2MuqVSRNQAsaE8lQh4P/20e/zOs9D6LLAbo2AYsxt7kQHkbCk/47k1oOG9/+2qfvude6H9WWA3AXTVXR1TfsZza0DD++mXl/XdDy+bQfavL0/vu4XdZ1PwBufGSqeIqAFiQ3kq0/De/fiu/vTzdX3z7DjItvCeFnaf1fU3h/BK1K5a1782g5EUfpZuH7eANoweuL3dNfHstPD0clrT/ysyZ5x1/KIt2BZMA8SG8lSm4T2NroeqYbd7dO0becPw7peTu4c5W8rPeG4NaHiN8vb8PrXmhTnSKBogNpSnMg3v/W/HYvdQ896ep8VOC7vPCO/6Eiga0PAaU7qPrh0Lo+Z5UY40igaIDeWpRMAbHQXsbVUDxIbyVAiviAaIDeWpEF4RDRAbylMhvCIaIDaUp0J4RTRAbChPhfCKaIDYUJ4K4RXRALGhPBXCK6IBYkN5KoRXRAPEhvJUCK+IBogN5akQXhENEBvKUyG8IhogNpSnkhPeUAQes1g1QHyA2EDxsdAG4V01QGyg+CC8UQHiA8QGig/CGxUgPkBsoPgAh5fBKBaEl6E2CC9DbRBehtogvAy1URReZ5/JVeLtbtd0WLOfdXY9+Vw02tYWIQ8rWWl9SB+Su7/tdld5D0dJeJ19JteJmyvDgPNlhbhtYAl5WMlK60P6kDQNlu6+v856OErCa/fVWTPuf782DDhfypu4efT3pslmwMM6Vo4+pA/JbUPlzVXWw1ESXruj2ZrRNli7GnRWc/ZZKxrNmQh5WMvKqVOt+CGZOA7JNkrCa/eSXDMO/z41Q43d09LZ4bKsjQM0IQ9rWWl/ieQPSdMkLOvh2ObI28bNFUfe3kcboofk00/Phm1GgUdeuZq3jXCBtYaDO4ia14JXzsfd35pvjGpqXruX5JrR/AN0/8dLu6els8Nl0WjORMjDWla68kXykBzZzXs4tjvP2zbE5jzv2YfwIXnbtta/UjPPy2AUDcLLUBuEl6E2CC9DbRBehtogvAy1QXhXijfVw9DH//dfaxnZUBDedeLz83/58k//xx/+8ut6XjYThHedeP/lf3771P8x4Z0ThHedePHV/z7/+vD64duqqp72r5+fV9WXfzbvHp6XMSKD8K4SH588rN8c6oZ2hP3w7dPu9XND9Juv/tG875ZJe9UThHeVeN+Ae+Dywz8dC9/utVl+ILuFuVvGiA3Cu0q8+LqpEA7/e1FVX7cLjq9vqjYetqNu9xkjMgjvGvHxSctoO99w+Ll/PVQM7QqnL2zdZ4yoILxrxJsWya6ePZQJ3ev7L46zDOfZhu4zRkQQ3hWiLRjal7bGPZDavX5+fhh6DwQ3zHbLpN3qCcK7QnRD7psvfn1/qB6a0bZ7babKvmjr3a/PyxiRQXgZaoPwMtQG4WWoDcLLUBuEl6E2CC9DbRBehtogvAy1QXgZaoPwMtTG/wPri8rvb17s8QAAAABJRU5ErkJggg==" />
<!-- rnb-plot-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxucGxvdF9ub3JtIDwtIGZ1bmN0aW9uKGJvdW5kLCBvdXRlcikge1xuICB4IDwtIHNlcSgtb3V0ZXIsb3V0ZXIsbGVuZ3RoPTEwMClcbiAgaHggPC0gZG5vcm0oeClcbiAgcGxvdCh4LCBoeCwgdHlwZT1cXG5cXCwgeGxhYj1cXHogdmFsdWVzXFwsIHlsYWI9XFxOb3JtYWwgUERGXFwsIG1haW49XFxOb3JtYWwgRGlzdHJpYnV0aW9uXFwsIGF4ZXM9RkFMU0UpXG4gIGxpbmVzKHgsIGh4KVxuICBpIDwtIHggPCAtYm91bmRcbiAgcG9seWdvbihjKC1vdXRlcix4W2ldLC1ib3VuZCwtYm91bmQpLCBjKDAsaHhbaV0sbWF4KGh4W2ldKSwwKSwgY29sPVxccmVkXFwpXG4gIGkgPC0geCA+IGJvdW5kXG4gIHBvbHlnb24oYyhib3VuZCwgYm91bmQsIHhbaV0sb3V0ZXIpLCBjKDAsbWF4KGh4W2ldKSwgaHhbaV0sMCksIGNvbD1cXHJlZFxcKVxuICBheGlzKDEsIGF0PS1vdXRlcjpvdXRlciwgcG9zPTApXG59XG5wbG90X25vcm0oMS45NiwgNClcbmBgYFxuYGBgIn0= -->
```r
```r
plot_norm <- function(bound, outer) {
x <- seq(-outer,outer,length=100)
hx <- dnorm(x)
plot(x, hx, type=\n\, xlab=\z values\, ylab=\Normal PDF\, main=\Normal Distribution\, axes=FALSE)
lines(x, hx)
i <- x < -bound
polygon(c(-outer,x[i],-bound,-bound), c(0,hx[i],max(hx[i]),0), col=\red\)
i <- x > bound
polygon(c(bound, bound, x[i],outer), c(0,max(hx[i]), hx[i],0), col=\red\)
axis(1, at=-outer:outer, pos=0)
}
plot_norm(1.96, 4)
<!-- rnb-source-end -->
<!-- rnb-plot-begin -->
<img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAYAAAAEgCAMAAACKBVRjAAAAxlBMVEUAAAAAADoAAGYAOjoAOmYAOpAAZmYAZpAAZrY6AAA6ADo6AGY6OgA6OmY6ZpA6kLY6kNtmAABmADpmOgBmOjpmZgBmZmZmkLZmkNtmtrZmtttmtv+QOgCQZgCQZjqQkGaQkLaQtpCQttuQtv+Q27aQ29uQ2/+2ZgC2Zjq2ZpC2kGa225C227a229u22/+2/7a2/9u2///bkDrbkGbbtmbbtpDb27bb29vb/7bb/9vb////AAD/tmb/25D/27b//7b//9v////VTPQFAAAACXBIWXMAAA7DAAAOwwHHb6hkAAALN0lEQVR4nO2dDXvbthVGYc9a49ZtpWbttkbdunWzujXJZjbJNoeOpP//p4YPkqIkEoAkgu8L+p7niWNLInDvPSRAyiaktgIUhQ7guSMCwIgAMCIAjAgAIwLAiAAwIgCMCAAjAsCIADAiAIwIAAMUsFLqlf5vvVAvTt52s1SzR/dtoSxXX73df7x+5T8eDrZaqeuHbRfmpccNJAYrwOQ6kACt4FWHgE8vr2MF2Jc+MwHmEBhMgLq67+rkqNp9AnqPjKSABeiUawHvf6fUb35vqlqoq7+8VNev9RPv79TV99tfb9WNLe77l3qbL+8PBNi6fdBPvagf//WuGpNMH/r5qsU39RHwWvdlhyxXdB3CvHrpv6uG96L5+8+qCmB4wAJ04pWAakc22btvZ/9dVHv2navidls2e/qxANPO7NE9XrWlH28E2Bb/Vwlonu4VcBRN9/E1AFAB19/q/J2Ap9uqcnOb8uzt9j/6Cf3TO5379+4JXd3PHs0r510C7Df2cdPio9lwXle4arGeA8wP78wB0xbgvncvOYjm5m31fQKwAn7RxXQCqgHY/lc0p0e6Cu6rrbrhwz/vWkONoSXg6t4+rl9882bXiRNgWqwF2J155Y6QTgHH0eg2T5+pYsAKeNDF+NuiVVBTRPelnp2bA0SXaPMnNxoEBOgvqh7BawG2xb3q2se6BXREc9apQgxgAbqwvzWp1fmVZnfrE2AKe/Pjh0VIgD6hdKP2Tz0CdtV1Tzu9OwEd0UxVgJviIo8AV6d1jwD9UDUH6B8+/fnOzaFRR8ChgGd0BJjE1PEc0C2grEfjrkn40/JgaNr8sBvkDwTYCaHVVam8c8CkBdhD4OgsqPcImD0eFXrvQsw+ruv5tR2HbIumrUMB1/VZkO3QNFl323UWNG0B5hA4ug7onwO6JuGjtyJ+3p25l/WF2J6A5jpAt+uYVy/tug6YtgCT9+5K+A/mO+9Z0M1fi9ZYv+1+M+6dmQG+tA+Yi9i3h3NAcyW8fdLz9Vev7VmQfemb9pVwK5opChAMIgCMCAAjAsCIADAiAIwIACMCwIgAMCIAjAgAIwLAiAAwIgCMCAAjAsCIADCTE2B/PYYO4gRyijUGl09GCvKJNIomnWwU5BJnHKrne2IyCTMO5fmJlTyijEN5fyQljyijOEoli9yyCDKKjkxySC6HGOPoyiSD7DIIMY7uRPjT448wjp48+NPjjzCK3jTo86MPMA4RgMWTBXuC7PHF4cuCPEPy8OLwJkGeIXl4cfiT4E6RO7o4Ajlwp8gdXRyhHKhzpA4ujmAK1DlSBxdHOAXmJJljiyMiA+YkmWOLQwRgiUqAOEvi0OIQAWDiEuBNkzeyOCLj503TRbZevALHcS6xlaU10BKwWaZZFzAl0XUVAWmIryurAREAJm8Bp1SV1IAIAJO1gJNqSi6gWXoQsob7mUxHQJ6cGDtnqpxRxTEtAW4UGvfjOy5kUgJWdjnfbZFolfwUnFxQSgNVUGX9+RBPt9m8KzQlAZtlsyzsKs0CsQk4vZ6MBurT0GbkKXI5DT2jmswCmoGnnLAARgMiAEy2As6rJZ8BEQAm1/eCzqwkq4D8EAFYzg6bLt86IPNBINlcA09QQHF1n9O7EJMT4H4ZVmTzXugFZWQz0D4NzeUMVATAuaSKZAZEAJgsBVxWQy4DIgBMlm9FTE9AXlwaM1XOVMFEIgKwXBwyVc5UwcQhArAMEDFT0kyxxCECwExSQEbXAYMUj8gAUShxiAAwwwTMk/YuknoY4h6CBqoco4DV7LF4Qf97yaEqR2OgdYPGfFvOHtl/LzllAa+2T1882H+8DFc3FgNNHOb38utv7kXAyOziMDcGrObcQ9CAVeMTYO6N0WdC1CdB0xbAz6DBkmROEkYc0xaQwYXYpAVslsyzr2XgknEY2LsOIGfaAjJYrGboilEY2AXx9Pk9MI4IBq8Xm4Bb8kl4+HoxGGgNQeRrFCSoFpUA+kk4RbUIDGQzCSepFZMA9kl48gKaP4zgnIQTlQpvAB9BHNMXQD4HTF8A91lQskLBDewCoL5B6RkI4J6E09UJbQDdfxwJo0QXAN1/HCmjBFeg1b19N+6K8WosaY1oBJR20dyCcdGatDXCGji6DiD8u6DEFSIRUF8HEJ6Npq4Q1EAGR0Dy+nAI4J0D0tcHaYD/LGiE6pAIIGWMCIFVoBcwSoBwAcS3qY5TG5yBw55XZJ9hMlJlWAToI4HrJGi0wsAM7HVcso0/z0xAQfcpViOWBWVg1+9mSTb8j1sUuAC64X87clFABupu+Yb/sUuCFVAowj/NHbkkGAPEF2KjFwRigPetiPEjEwF7ACJDFINWwHPZHVkFIGfEqXcZBfisfMo9RgG/Lp1wj1GoZ2OAU4D6CBMwtnpKAeojTsDYXTMK0EVAChi3b0IBpgYiAIetP1TAqJ3TCXD1xwr4OOJMzCagqj9YwIj9cwlQdf5wAaMFQCWgyZ5AwFgRMAlo1Z9AwEgTAY8A1a4/g4BxgmARsF9+EgFjHAQkAg7KzyJgBAUUAg53fyIBySPBC1Ad5WcS8DHtG6TYG9RUd/W5BDgHqSwA33fvqz2hgIQWRn7zu8FTe1YBOwtDqkiWqHqWnFGn4Usf1bD/6Ys25m17oE0GaZi2SCLg4o152x5ok0Eapi2SCLh4Y962B9pkkIZpiyQCLt6Yt+2BNhmkYdoiiYCLN+Zte6BNBmmYtkiTESDEIQLAiAAwIgCMCAAjAsCIADAiAIwIACMCwIgAMCIAjAgAk1LAyrMKZuFbqH2zVMq/hOPTF33repWBFeD7twz36w3a4Mu4h4QCSs8ypMX1w7bsS2az1M94V7FbL/oWVjON9jbs3TLcrzdo2/kZC6+mE7Be9IezXsw9n+H9dGsW0Cz6K1X2rmznPgdk1VvD/i3D/fqD3voz7iWdgGL2gz+cwIeo9+9rpZr3fdJNoIaeLSP6NfiCDmbcRTIBT5/fB0bEIjCe+irVK8B+KrWvyEEB3n59QYcz7iKVADMWeMMpA9Ns6X26r4xu9/XtxCEB3n59QQcz7iaVAPNhTIFwNkvP86V/JdlkAgL9eoKOyLiLRALsUBAKx1Mn//6fbggK9dsfdFTGHQwvoDDn0kX19/LH64EXu1NtN2N2PVv01KHZ+MxJ2LNl1X5wAfmOoJvQujP2g7kQc1n01iL8YWZ9m4ZOQ/0C/P0GgjZQHAE7POGYp3o/xf7pNrgf9pYhdCHmK2CoX2/QzStOBPVWxMpztFZH83llDL1d4BEQ7NcXtHsBlwAhAhEARgSAEQFgRAAYEQBGBIARAWBEABgRAEYEgBEBYEQAGBEARgSAEQFgRAAYEQBGBIARAWBEABgRAEYEgBEBYEQAGBEARgSAEQFgRAAYEQBGBIARAWAmKGC9OPU+LSQiAIwIAJOrAHc/l7lhrrp58eq+voVLC7AO7JeiedUZt5COQa4CDO52YLd8xmr2aG5ONUs5tAWYn83Nj9ZSyWggYwH1LdPGw3oxX39z7x5rCbArzJgbI4NLdMDIV0CzaoOpdnVrcGmGmZYA97C24lunCUu+AnYLx+ixx9ygq4f761/2j4CyWT/ArkZGOALlK6B1W3t5/S9dcDsiHQxB+/fMr/wLFGHIVUB72Zj14juzmpsZZMp6CJpXI097p6c8P81VwN7KVSuzyI/b+fU5p66z0aMHHXcWZPZ8eyj4FyMDkakAO6Q3o7o7v7SrROhimx1dm1B/rK8D3KERWnISRKYCpoMIACMCwIgAMCIAjAgAIwLAiAAwIgCMCAAjAsCIADAiAIwIACMCwIgAMCIAjAgAIwLAiAAwIgCMCAAjAsCIADD/B7fTNZ2sWdnBAAAAAElFTkSuQmCC" />
<!-- rnb-plot-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuYW5vdmEobW9kMiwgbW9kMywgdGVzdD1cXENoaXNxXFwpXG5gYGBcbmBgYCJ9 -->
```r
```r
anova(mod2, mod3, test=\Chisq\)
<!-- rnb-source-end -->
<!-- rnb-output-begin eyJkYXRhIjoiQW5hbHlzaXMgb2YgVmFyaWFuY2UgVGFibGVcblxuTW9kZWwgMTogcmV2dF9ncm93dGggfiBhdF9ncm93dGhcbk1vZGVsIDI6IHJldnRfZ3Jvd3RoIH4gbGN0X2dyb3d0aCArIGNoZV9ncm93dGggKyBlYml0X2dyb3d0aFxuICBSZXMuRGYgICAgUlNTIERmIFN1bSBvZiBTcSBQcig+Q2hpKSAgXG4xICAgICAyNiAxLjU1MzQgICAgICAgICAgICAgICAgICAgICAgICBcbjIgICAgIDI0IDEuMTkxOCAgMiAgIDAuMzYxNjggICAwLjAyNjIgKlxuLS0tXG5TaWduaWYuIGNvZGVzOiAgMCDDouKCrMucKioqw6LigqzihKIgMC4wMDEgw6LigqzLnCoqw6LigqzihKIgMC4wMSDDouKCrMucKsOi4oKs4oSiIDAuMDUgw6LigqzLnC7DouKCrOKEoiAwLjEgw6LigqzLnCDDouKCrOKEoiAxXG4ifQ== -->
Analysis of Variance Table
Model 1: revt_growth ~ at_growth Model 2: revt_growth ~ lct_growth + che_growth + ebit_growth Res.Df RSS Df Sum of Sq Pr(>Chi)
1 26 1.5534
2 24 1.1918 2 0.36168 0.0262 — Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
<!-- rnb-output-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuYGBgclxuIyBFbnN1cmUgZmlybXMgaGF2ZSBhdCBsZWFzdCAkMU0gKGxvY2FsIGN1cnJlbmN5KSwgYW5kIGhhdmUgcmV2ZW51ZVxuIyBkZiBjb250YWlucyBhbGwgcmVhbCBlc3RhdGUgY29tcGFuaWVzIGV4Y2x1ZGluZyBOb3J0aCBBbWVyaWNhXG5kZl9jbGVhbiA8LSBmaWx0ZXIoZGYsIGRmJGF0PjEsIGRmJHJldnQ+MClcbiMgV2UgY2xlYW5lZCBvdXQgNTc4IG9ic2VydmF0aW9ucyFcbnByaW50KGMobnJvdyhkZiksIG5yb3coZGZfY2xlYW4pKSlcbmBgYFxuYGBgXG5gYGAifQ== -->
```r
```r
```r
# Ensure firms have at least $1M (local currency), and have revenue
# df contains all real estate companies excluding North America
df_clean <- filter(df, df$at>1, df$revt>0)
# We cleaned out 578 observations!
print(c(nrow(df), nrow(df_clean)))
```r
```r
forecast4 <-
lm(revt_lead ~ revt + act + che + lct + dp + ebit , data=df_clean)
tidy(forecast4)
```r
```r
forecast3.1 <-
lm(revt_lead ~ revt + act + che + lct + dp + ebit + factor(isin),
data=df_clean[df_clean$fic==\SGP\,])
# n=7 to prevent outputting every fixed effect
print(tidy(forecast3.1), n=15)
glance(forecast2)
anova(forecast1, forecast2, test="Chisq")
Analysis of Variance Table
Model 1: revt_lead ~ lct + che + ebit
Model 2: revt_lead ~ revt + act + che + lct + dp + ebit
Res.Df RSS Df Sum of Sq Pr(>Chi)
1 24 3059182
2 21 863005 3 2196177 1.477e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# Note the group_by -- without it, lead() will pull from the subsequent firm!
# ungroup() tells R that we finished grouping
df_clean <- df_clean %>%
group_by(isin) %>%
mutate(revt_lead = lead(revt)) %>%
ungroup()
forecast3 <-
lm(revt_lead ~ revt + act + che + lct + dp + ebit,
data=df_clean[df_clean$fic=="SGP",])
tidy(forecast3)
glance(forecast3)
forecast4 <-
lm(revt_lead ~ revt + act + che + lct + dp + ebit , data=df_clean)
tidy(forecast4)
glance(forecast4)
forecast3.1 <-
lm(revt_lead ~ revt + act + che + lct + dp + ebit + factor(isin),
data=df_clean[df_clean$fic=="SGP",])
# n=7 to prevent outputting every fixed effect
print(tidy(forecast3.1), n=15)
glance(forecast3.1)
anova(forecast3, forecast3.1, test="Chisq")
Analysis of Variance Table
Model 1: revt_lead ~ revt + act + che + lct + dp + ebit
Model 2: revt_lead ~ revt + act + che + lct + dp + ebit + factor(isin)
Res.Df RSS Df Sum of Sq Pr(>Chi)
1 324 14331633
2 304 13215145 20 1116488 0.1765
library(fixest)
forecast3.2 <-
feols(revt_lead ~ revt + act + che + lct + dp + ebit | isin,
data=df_clean[df_clean$fic=="SGP",])
NOTE: 29 observations removed because of NA values (LHS: 21, RHS: 8).
summary(forecast3.2)
OLS estimation, Dep. Var.: revt_lead
Observations: 331
Fixed-effects: isin: 21
Standard-errors: Clustered (isin)
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
RMSE: 199.8 Adj. R2: 0.843506
Within R2: 0.780586
df_clean %>%
filter(fic=="SGP") %>%
group_by(isin) %>%
mutate(mean_revt_lead=mean(revt_lead, na.rm=T)) %>%
slice(1) %>%
ungroup() %>%
ggplot(aes(x=mean_revt_lead)) + geom_histogram(aes(y = ..density..)) + geom_density(alpha=.4, fill="#FF6666")
`stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# Exports for the following week
save(df_clean, forecast2, uol, forecast4, file = "../../Data/Session_2_export.RData")