stylistic changes

Author: bettinagruen

vignettes/Chapter08.Rmd

@@ -607,7 +607,7 @@ gen.proposal.poisson <- function(y, X, e, b0 = 0, B0 = 100, t.max = 20) {
   betas <- matrix(NA_real_, ncol = t.max, nrow = d)
   beta.new <- matrix(c(log(mean(y/e)), rep(0, d - 1)), nrow = d) 
   
-  B0.inv=solve(B0)
+  B0.inv <- solve(B0)
   for (t in seq_len(t.max)) {
     beta.old <- beta.new
     

vignettes/Chapter09.Rmd

@@ -17,7 +17,7 @@ knitr::opts_chunk$set(
   fig.height = 6
 )
 knitr::opts_knit$set(global.par = TRUE)
-pdfplots <- FALSE # default: FALSE; set this to TRUE only if you like pdf figures
+pdfplots <- TRUE # default: FALSE; set this to TRUE only if you like pdf figures
 ```
 
 ```{r, include = FALSE}
@@ -349,8 +349,12 @@ points(pk3, col = 3, cex = 1.5, pch = 16)
 # Section 9.4 Posterior Predictive Distributions in Regression Analysis
 
 ## Example 9.12: Road Safety Data; potential outcome analysis
-We are now interested in predicting the number of children who would have been killed or seriously injured without the legal intervention on October 1, 1994.
-To do so we reuse functions defined in Example 8.8.
+
+We are now interested in predicting the number of children who would
+have been killed or seriously injured without the legal intervention
+on October 1, 1994. To do so we reuse functions defined in Example
+8.8.
+
 ```{r}
 gen.proposal.poisson <- function(y, X, e, b0 = 0, B0 = 100, t.max = 20) {
   N <- length(y)
@@ -359,7 +363,7 @@ gen.proposal.poisson <- function(y, X, e, b0 = 0, B0 = 100, t.max = 20) {
   betas <- matrix(NA_real_, ncol = t.max, nrow = d)
   beta.new <- matrix(c(log(mean(y/e)), rep(0, d - 1)), nrow = d)
 
-  B0.inv=solve(B0)
+  B0.inv <- solve(B0)
   for (t in seq_len(t.max)) {
     beta.old <- beta.new
 
@@ -417,7 +421,7 @@ sample_beta<- function(y,X,e, b0, B0, qmean, qvar, beta.old){
     beta <- beta.old
     acc  <- 0
   }
-  return(res = list(beta=beta, acc=acc))
+  return(res = list(beta = beta, acc = acc))
 }
 
 poisson <- function(y, X, e, b0 = 0, B0 = 100, burnin = 1000L, M = 10000L) {
@@ -449,36 +453,41 @@ poisson <- function(y, X, e, b0 = 0, B0 = 100, burnin = 1000L, M = 10000L) {
   return(list(beta.post = beta.post, accept = mean(acc)))
 }
 ```
-Our goal is to predict the number of killed or seriously injured children
-from October 1994  (i.e when the  legal intervention giving priority to pedestrians  became effective) using only data before that time point. Hence we estimate the model in Example 8.8. with an intercept and a holiday effect using only the information up to September 1994.
+
+Our goal is to predict the number of killed or seriously injured
+children from October 1994 onward (i.e., when the legal intervention
+giving priority to pedestrians became effective) using only data
+before that time point. Hence we estimate the model in Example
+8.8. with an intercept and a holiday effect using only the information
+up to September 1994.
 
 ```{r}
 data("accidents", package = "BayesianLearningCode")
-y <- accidents[1:93, "children_accidents"]
-t=length(y)
-e <- accidents[1:93, "children_exposure"]
+y <- window(accidents[, "children_accidents"], end = c(1994, 9))
+e <- window(accidents[, "children_exposure"],  end = c(1994, 9))
+t <- length(y)
 
 X <- cbind(intercept = rep(1, length(y)),
-           holiday = c(rep(rep(c(0, 1, 0), c(6, 2, 4)),7), rep(0,6), rep(1,2),0))
-M=10000
+           holiday = rep(rep(c(0, 1, 0), c(6, 2, 4)), length.out = t))
+M <- 10000
 set.seed(1)
-res <- poisson(y, X, e, b0 = 0, B0 = 100,M=M)
+res <- poisson(y, X, e, b0 = 0, B0 = 100, M = M)
 ```
 
 We define the covariates for the time points from October 1994 on and use the 
 draws from the posterior distribution to predict the values of the time series.
 
 ```{r}
-X.pred <- cbind(intercept = rep(1, 99),
-                holiday = c( rep(0,3),rep(rep(c(0, 1, 0), c(6, 2, 4)),8)))
-e.pred <- accidents[94:192, "children_exposure"]
-t.pred=length(e.pred)
+e.pred <- window(accidents[, "children_exposure"], start = c(1994, 10))
+t.pred <- length(e.pred)
+X.pred <- cbind(intercept = rep(1, t.pred),
+                holiday = c(rep(0, 3), rep(rep(c(0, 1, 0), c(6, 2, 4)),
+                                           length.out = t.pred - 3)))
 
-lambda = matrix(e.pred, ncol=M, nrow=t.pred)*exp(X.pred%*%t(res$beta.post))
+lambda <- matrix(e.pred, ncol = M, nrow = t.pred) * exp(X.pred %*% t(res$beta.post))
 
 set.seed(2)
-pred<-matrix(rpois(M*t.pred, lambda), ncol=M, nrow=t.pred)
-
+pred <- matrix(rpois(M * t.pred, lambda), ncol = M, nrow = t.pred)
 ```
 
 We reuse also the function to determine mean and quantiles of the draws.
@@ -489,35 +498,39 @@ res.mcmc <- function(x, lower = 0.025, upper = 0.975) {
                   paste0(upper * 100, "%"))
   res
 }
-pred.int<-t(round(apply(pred, 1, res.mcmc), 3))
+pred.int <- t(round(apply(pred, 1, res.mcmc), 3))
 ```
 
-Then we plot the predictive mean together with the (equal-tailed) 95\%  prediction intervals.
+Then we plot the predictive mean together with the (equal-tailed) 95%
+prediction intervals.
 
 ```{r, echo = -c(1:2)}
 if (pdfplots) {
   pdf("9-4_1.pdf", width = 8, height = 4)
   par(mar = c(2.7, 2, 1.5, .5), mgp = c(1.6, .6, 0))
 }
 par(mfrow = c(1, 1))
-plot(1:192, accidents[,"children_accidents"], type="p" ,ylim=c(0,7),
-     xaxt="n",xlab="", ylab="",
-     main="Number of children killed or seriously injured")
-axis(1, at = seq(1, 192, by = 12), las=2, labels=1987:2002)
-matplot(94:192, pred.int, col="blue", type="l",lty=c(2,1,2),lwd=2,add=T)
-abline(v=94, col="red")
-
+plot(time(accidents), accidents[, "children_accidents"], type = "p", ylim = c(0, 7),
+     xlab = "", ylab = "",
+     main = "Number of children killed or seriously injured")
+matplot(as.vector(time(e.pred)), pred.int, col = "blue", type = "l",
+        lty = c(2, 1, 2), lwd = 2, add = TRUE)
+abline(v = 1994.75, col = "red")
 ```
 
-We see that the prediction intervals after the intervention are  much too wide which again indicates that there is an intervention effect. Note that whereas the risk 
-for a child to be killed or seriously injured in this model is constant the predicted mean number of killed and seriously injured children decreases from 1994 due to the decreasing number of exposed in that time period.  
+We see that the prediction intervals after the intervention are much
+too wide which again indicates that there is an intervention
+effect. Note that whereas the risk for a child to be killed or
+seriously injured in this model is constant the predicted mean number
+of killed and seriously injured children decreases from 1994 due to
+the decreasing number of exposures in that time period.
 
 # Section 9.5: Bayesian Forecasting of Time Series
 
 ## Example 9.13: US GDP data - one-step-ahead forecasting
 
-For creating the design matrix for an AR model, we re-use the function from
-Chapter 7.
+For creating the design matrix for an AR model, we re-use the function
+from Chapter 7.
 
 ```{r}
 ARdesignmatrix <- function(dat, p = 1) {