Estimating GARCH Models"
In bayesforecast: Bayesian Time Series Modeling with Stan

knitr::opts_chunk$set(collapse = TRUE,comment = "#>")
knitr::opts_chunk$set(echo = TRUE)

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  dev = "png",
  dpi = 150,
  fig.asp = 0.8,
  fig.width = 10,
  out.width = "60%",
  fig.align = "center")

library(bayesforecast)
library(ggplot2)

Introduction

The foreign exchange consists of converting one currency into another at a specific rate known as the foreign exchange rate. The conversion rates are constantly changing as they are driven by the market forces of supply and demand. In this vignette, we will model the volatility of the series of daily observations of the foreign exchange between Germany and the United Kingdom proposed by Ardia y Hoogerheide (2010) using the stan_garch function of the bayesforecast package. The records correspond from January 3, 1984, to December 31, 1991. However, for this work, we will use the first 350 observations.

The steps for estimating the model are:

Plot the data and identify any unusual observations.
Create de GARCH Model through the stan_garch function of the **bayesforecast** package.
Plot and observe the residuals of the model. If the residuals look like white noise, we proceed to make the prediction. Otherwise, we will choose another model. Plot the data and identify any unusual observations.

Plotting the data:

autoplot(demgbp,main ="DEM/GBP Foreign exchange",ylab ="log-returns",x = "Days" )

The model proposed by Ardia y Hoogerheide (2010) in their bayesGARCH package is a GARCH(1,1) with t-student innovations of the form:

$$ y_t = \mu + \left(\dfrac{v-2}{v}\lambda_t \sigma_t \right)^{1/2} \epsilon_t,\ \sigma^2_t = \sigma_0 + \alpha_1\epsilon_{t-1}^2 + \beta_1\sigma^2_{t-1}\ \epsilon_t \sim N(0,1)\ \lambda_t \sim IG(v/2,v/2) $$

where $\sigma_0 > 0$,$\alpha_1,\beta_1$ and $v > 2$, moreover $N(0,1)$ denotes the standard normal distribution; $IG$ denotes the inverted gamma distribution. The restriction on the degrees of freedom parameter $v$ ensures the conditional variance to be finite and the restrictions on the GARCH parameters $\sigma_0, \alpha_1$ and $\beta$ guarantee its positivity.

This model can be run in bayesforecast like:

set.seed(6551)
sf1 = stan_garch(ts = demgbp,order = c(1,1,0),genT = TRUE,chains = 1)

sf1 = stan_garch(ts = demgbp,order = c(1,1,0),genT = TRUE,chains = 1)
summary(sf1)

summary(sf1)

In the model proposed by Ardia y Hoogerheide (2010) they use truncated normal priors on the GARCH parameters $\beta$ and $\alpha=(\sigma_0,\alpha_1)'$

$$ p(\alpha)\propto\phi_{N_2}(\alpha|\mu_\alpha,\Sigma_\alpha) 1 {\alpha\in R_+^2}\ p(\beta)\propto\phi_{N_1}(\beta|\mu_\beta,\Sigma_\beta) 1 {\beta\in R_+}\ $$ where $\mu_\bullet$ and $\Sigma_\bullet$ are the hyperparameters, $1{\bullet}$ is the indicator function and $\phi_{N_d}$ is the d-dimensional normal density.

The prior distribution of vector $\lambda=(\lambda_1,...\lambda_T)'$ conditional on $v$ is found by noting that the components $\lambda_t$ are independent and identically distributed from the inverted gamma, which yields

$$ p(\lambda|v)=\left(\frac{v}{2} \right)^{\frac{Tv}{2}}\left[ \Gamma\left(\frac{v}{2}\right)\right]^{-T}\left(\prod_{t=1}^T\lambda_t\right)^{-\frac{v}{2}-1}exp\left[-\frac{1}{2}\sum_{t=1}^T\frac{v}{\lambda_t}\right] $$ furthermore the prior distribution on the degrees of freedom parameters is a translated exponencial with parameters $\lambda^*>0$ y $\delta\geq2$

$$ p(v)=\lambda^exp[-\lambda^(v-\delta)]1{v>\delta} $$ The priors used by Ardia y Hoogerheide (2010) are slightly similar to the prior by default in bayesforecast, these were defined by default following the recommendations proposed by Stan.

prior_summary(sf1)

Based on these previously defined parameters, the obtained posteriors are:

mcmc_plot(sf1)

It is observed that the posteriors converge. Additionally, when following the model's fit, it does not entirely fit the data because the model does not capture the dependency structure in the location parameter. Then we proceed to review the model residuals. The check_residuals function estimates the posterior mean of the residuals and plots them. Note that the plot is not sufficient to corroborate the normality and stationarity assumptions, but they are an initial indicator of the adjustment.

check_residuals(sf1)

The residual series (Upper part) seems to be stationary. However, the histogram and quantile graph (middle part) show that the model has heavy tails due to the series's high volatility.

Forecast

Based on the previous results, we finally predict the model for the next few days:

autoplot(object = forecast(sf1,h = 12),ylab="log-returns")

References

Ardia, D., & Hoogerheide, L. F. (2010). Bayesian estimation of the GARCH(1,1) model with student-t innovations. The R Journal, 2(2), 41-47.
Ardia, D. (2008). Financial Risk Management with Bayesian Estimation of GARCH Models: Theory and Applications. volume 612 series Lecture Notes in Economics and Mathematical Systems. Springer-Verlag, Berlin, Germany. https://doi.org/10.1007/978-3-540-78657-3
Vrontos, I., Dellaportas, P., & Politis, D. (2000). Full Bayesian Inference for GARCH and EGARCH Models. Journal of Business & Economic Statistics, 18(2), 187-198.
Engle, R. (1982). Autoregressive Conditional Heteroscedasticity with Estimates of the Variance of United Kingdom Inflation. Econometrica, 50(4), 987-1007. url: http://www.jstor.org/stable/1912773.
Bollerslev, T. (1986). Generalized autoregressive conditional heteroskedasticity. Journal of Econometrics. 31(3), 307-327. doi: https://doi.org/10.1016/0304-4076(86)90063-1