Deploying to gh-pages from @ 2c2e862 🚀

Author: gregorkastner

articles/Chapter09.html

@@ -30,6 +30,8 @@ abline(h = probs, lty = 3)
 mtext(probs, side = 2, at = probs, adj = c(0, 1), cex = .8, col = "dimgrey")
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-3-1.png)
+
 ### Example 9.2: CHF exchange rate data - single predictions
 
 We load the data and then plot the pdf and cdf for the predictive
@@ -64,6 +66,8 @@ abline(h = probs, lty = 3)
 mtext(probs, side = 2, at = probs, adj = c(0, 1), cex = .8, col = "dimgrey")
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-4-1.png)
+
 We inspect the parameters of the Student-t distribution.
 
 ``` r
@@ -98,6 +102,8 @@ for (i in seq_along(Ns)) {
 }
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-6-1.png)
+
 ### Example 9.4: Road safety data - posterior predictive credible interval
 
 We verify that a 95% posterior predictive interval is given by \[1, 9\]
@@ -232,6 +238,8 @@ abline(h = q_t, col = 2, lty = 2, lwd = 1.5)
 legend("topleft", c("Normal", "Student t"), lty = 1:2, col = c(4, 2), lwd = 1.5)
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-14-1.png)
+
 ### Example 9.9: Road safety data - sampling-based prediction
 
 We use a sampling-based approach to obtain draws from the posterior
@@ -265,6 +273,8 @@ abline(h = probs, lty = 3)
 mtext(probs, side = 2, at = probs, adj = c(0, 1), cex = .8, col = "dimgrey")
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-16-1.png)
+
 ### Example 9.11: Predicting the probability of future successes
 
 We illustrate the variance of the purely sampling-based estimator versus
@@ -333,6 +343,8 @@ boxplot(pk2, xlab = "k", range = 0, main = "Rao-Blackwellized",
 points(pk3, col = 3, cex = 1.5, pch = 16)
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-20-1.png)
+
 ## Section 9.4 Posterior Predictive Distributions in Regression Analysis
 
 ### Example 9.12: Road Safety Data; potential outcome analysis
@@ -506,6 +518,8 @@ matplot(as.vector(time(e.pred)), pred.int, col = "blue", type = "l",
 abline(v = 1994.75, col = "red")
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-25-1.png)
+
 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
@@ -581,6 +595,8 @@ for (p in 2:4) {
 }
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-28-1.png)
+
 ### Example 9.15: US GDP data - multi-step forecasting
 
 Now, we want to “sample the future” up to 12 steps ahead for $`p = 2`$.
@@ -644,6 +660,8 @@ for (i in 1:4) {
 }
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-30-1.png)
+
 And now we plot the predictions.
 
 ``` r
@@ -665,6 +683,8 @@ polygon(pxs, c(quants["5%", 1], quants["95%", ], rev(quants["5%", ])),
 abline(h = 0, lty = 3)
 ```
 
+![](Chapter09_files/figure-html/unnamed-chunk-31-1.png)
+
 ### Example 9.16: US GDP data - forecasting non-linear functionals
 
 Let us compute the probability of seeing negative growth rates a least
@@ -709,3 +729,107 @@ polygon(pxs, c(quants["5%", 1], quants["95%", ], rev(quants["5%", ])),
         col = rgb(1, 0, 0, .2), border = NA)
 abline(h = 0, lty = 3)
 ```
+
+![](Chapter09_files/figure-html/unnamed-chunk-34-1.png)
+
+### Section 2.5.3: Forecasting Volatility via GARCH(1,1) Models
+
+We start by defining the parameter transformation.
+
+``` r
+
+unrestricted2restricted <- function(theta) {
+  alpha0 <- exp(theta[1])
+  exptheta2 <- exp(theta[2])
+  exptheta3 <- exp(theta[3])
+  denom <- 1 + exptheta2 + exptheta3
+  alpha1 <- exptheta2 / denom
+  gamma1  <- exptheta3 / denom
+  c(alpha0, alpha1, gamma1)
+}
+```
+
+Next, we define the log likelihood function and the log posterior.
+
+``` r
+
+loglik <- function(y, theta) {
+  pars <- unrestricted2restricted(theta)
+  sigma2 <- rep(NA_real_, length(y))
+
+  # initialize at unconditional variance
+  sigma2[1] <- pars[1] / (1 - pars[2] - pars[3])
+
+  for (i in 2:length(y)) {
+    sigma2[i] <- pars[1] + pars[2] * y[i - 1]^2 + pars[3] * sigma2[i - 1]
+  }
+
+  -0.5 * log(2 * pi) * length(y) - sum(log(sigma2) + y^2 / sigma2)
+}
+
+logpost <- function(y, theta, priormeans = rep(0, 3), priorsds = rep(10, 3)) {
+  loglik(y, theta) +
+    dnorm(theta[1], priormeans[1], priorsds[1], log = TRUE) +
+    dnorm(theta[2], priormeans[2], priorsds[2], log = TRUE) +
+    dnorm(theta[3], priormeans[3], priorsds[3], log = TRUE)
+}
+```
+
+Finally, we define our random walk MH posterior sampler.
+
+``` r
+
+sampler <- function(y, M = 20000L, nburn = 5000L, propsds = c(0.1, 0.1, 0.1)) {
+  # allocate storage
+  thetas <- matrix(NA_real_, M, 3)
+
+  # starting values
+  theta <- rep(0, 3)
+  currentlogpost <- logpost(y, theta)
+
+  naccepts <- 0L
+  for (m in seq_len(nburn + M)) {
+    prop <- rnorm(3, theta, propsds)
+    proplogpost <- logpost(y, prop)
+    logR <- proplogpost - currentlogpost
+
+    if (log(runif(1)) < logR) {
+      theta <- prop
+      currentlogpost <- proplogpost
+      if (m > nburn) naccepts <- naccepts + 1L
+    }
+
+    if (m > nburn) thetas[m - nburn, ] <- theta
+  }
+  thetas
+}
+```
+
+We are now ready to apply the sampler to our EUR-CHF exchange rate
+percentage log returns. After running the sampler, we need to transform
+back to our original parameters.
+
+``` r
+
+data("exrates", package = "stochvol")
+y <- 100 * diff(log(exrates$USD / exrates$CHF))
+
+# Run the sampler
+thetas <- sampler(y)
+
+# Transform back to the usual parameterization
+draws <- t(apply(thetas, 1, unrestricted2restricted))
+colnames(draws) <- c("alpha0", "alpha1", "gamma1")
+
+# Visually check convergence
+plot.ts(thetas)
+```
+
+![](Chapter09_files/figure-html/unnamed-chunk-38-1.png)
+
+``` r
+
+plot.ts(draws)
+```
+
+![](Chapter09_files/figure-html/unnamed-chunk-38-2.png)

articles/Chapter09.md

@@ -32,7 +32,7 @@ knitr::kable(cbind(BF, logBF, PrM1, PrM2))
 |  403.4287935 |     6 | 0.9975274 | 0.0024726 |
 | 1096.6331584 |     7 | 0.9990889 | 0.0009111 |
 
-### Example 10.5: CHF exchange rate data - Testing the null hypothesis $`\mu = 0`$
+### Example 10.5: CHF exchange rate data - Testing for zero mean
 
 Before producing Figure 10.1, we begin with a concrete example. We load
 the data, compute the required parameters, and evaluate the two marginal
@@ -801,7 +801,7 @@ argmax <- which(probs[5:8, ] == max(probs[5:8, ]), arr.ind = TRUE)
 plot(xis, post_xis[argmax[1, 1], argmax[1, 2], ], type = "l", lwd = 2,
      xlim = c(0.27, 0.43),
      ylim = c(0, 1.17 * max(post_xis[probs[5:8, ] > 0.0001])),
-     lty = 2, main = "Most probable and BMA posterior", xlab = expression(xi))
+     lty = 2, main = "MAP and BMA posterior", xlab = expression(xi))
 lines(xis, bma_post_xis, lwd = 2)
 
 for (i in seq_len(nrow(weighted_post_xis))) {
@@ -812,7 +812,7 @@ for (i in seq_len(nrow(weighted_post_xis))) {
   }
 }
 
-legend("topright", c("Most probable posterior", "BMA posterior",
+legend("topright", c("MAP posterior", "BMA posterior",
                      "Weighted components"),
        col = 1, lty = c(2, 1, 1), lwd = c(2, 2, 1))
 abline(h = 0, lwd = 1.5)

articles/Chapter09_files/figure-html/unnamed-chunk-25-1.png

@@ -13,7 +13,7 @@ articles:
   Chapter09: Chapter09.html
   Chapter10: Chapter10.html
   Chapter11: Chapter11.html
-last_built: 2026-05-20T18:12Z
+last_built: 2026-05-20T21:07Z
 urls:
   reference: https://gregorkastner.github.io/BayesianLearningCode/reference
   article: https://gregorkastner.github.io/BayesianLearningCode/articles