Forecast bayesian GARCH model

Question

I am using this package in R to do Bayesian estimation of GARCH models. I want to forecast $y_t$ (i.e. the mean equation), but it seems that the package has no built-in function for this. The model setup in the package is as follows:

Now I was thinking how I should implement a forecasting procedure. I thought of doing this:

1. Estimate model to obtain parameter posteriors
2. Take mean of each posterior to get a point estimate
3. Obtain $h_t$
4. Simulate from normal and gamma to obtain $y_t$.

My questions/problems:

1. Is this a proper way of doing it?
2. This might take very long, because every forecast means that I have to simulate three times, first to get the posteriors, than for $\epsilon_t$ and then for $\bar{\omega}_t$
3. Are there any better options (preferably already coded)?

Answer 1

I am sure you already know $E[y_t]=0$ as it follows directly from the first equation you have listed. So when you say "I want to forecast $y_t$", I interpret it as "I would like to find the forecast distribution of $y_t$". Note that because $\varepsilon \stackrel{iid}{\sim} \mathcal{N}(0,1)$ we know that $$y_t \stackrel{iid}{\sim} \mathcal{N}(0,\frac{\nu-2}{\nu}\varpi_t h_t)$$

where $\frac{\nu-2}{\nu}\varpi_t h_t$ is the variance of $y_t$. For notation purpose define $$f(y_t\,|\,\nu,\varpi_t,h_t):=\mathcal{N}(0,\frac{\nu-2}{\nu}\varpi_t h_t)$$

Since the variance of $y_t$ is a random variable, the marginal distribution of $y_t$, $f(y_t)$, is not directly apparent. In this case, we can use a Gibbs sampler approach to draw from $f(y_t\,|\,\nu,\varpi_t,h_t)$ many times, each for different values of $h_t$ and $\varpi_t$ of whom were previously drawn from their own distributions. The resulting sample will in effect be drawn from $f(y_t)$ .

I am not familiar with this R function. Posting the output of the function would be helpful. I am assuming it returns a posterior sample of the $\alpha_0,\alpha_1,\beta$ and $\varpi$. Suppose the sample has $G$ draws. Then

1. For each posterior draw $g=1,2,\cdots,G$, compute $h^{(g)}_t$
2. simulate $y^{(g)}_t \sim \mathcal{N}(0,\frac{\nu-2}{\nu}\varpi^{(g)} h^{(g)}_t)$. Equivalently, you could simulate $G$ $\varepsilon_t$'s and multiply each one by a respective $(\frac{\nu-2}{\nu}\varpi^{(g)}_t h^{(g)}_t)^{1/2}$

The second step can be done multiple times each draw.

I would not take the sample mean of the posterior and use it to generate an estimate for $h_t$. This is because there is a non-linearity in the recursion between $\beta$ and the other variables. So $h_t$ evaluated at the mean of the posterior for the other variables may not be equal to $E[h_t]$ which could inturn lead to bias.

As far as something already coded goes...I do not know of anything off hand. If you are unsatisfied or would like help with coding, then editing your post to include the code you already have and the resulting output would help users answer the question more effectively.