In the **rugarch** package there are two main ways to simulate a GARCH process:

- The ugarchsim method which takes an already estimated object of class
*uGARCHfit*. - The ugarchpath method which takes a specification of class
*uGARCHspec*with fixed parameters.

Before proceeding to a demonstration to show how to obtain equivalence between the 2, it is important to say a few words about the initialization of the simulation which has led to some confusion. Consider the typical GARCH(1,1) process:

\[ \begin{gathered}

{y_{t + 1}} = {\mu _{t + 1}} + {\varepsilon _{t + 1}},{\text{ }}{\varepsilon _{t + 1}} = {z_{t + 1}}\sqrt {{h_{t + 1}}}\\

{h_{t + 1}} = \omega + {\alpha _1}\varepsilon _t^2 + {\beta _1}{h_t}\\

\end{gathered} \]

Now assume that we wish to simulate \( t+1 \) given estimated information upto time \( t \). As can be seen from the equation, we require a value for \( \varepsilon_{t}^2 \). In the **rugarch** package, the **ugarchsim** method has an option (‘startMethod’) of using either the estimated object’s last known values at time \( t \) to initialize the simulation (“sample”) or their unconditional long run values (“unconditional”). In the former case, \( \varepsilon_{t}^2 \) is based on the last value from the estimated object, whilst in the latter case \( z_{t} \) is set to zero which means that \( \varepsilon_{t}^2 \) is also zero. In either case, it is clear that there is no distributional ‘uncertainty’ in the value of \( h_{t+1} \). At \( t+2 \) however, we need a value for \( \varepsilon_{t+1}^2 \) which we do not have from the estimated object which stops at time \( t \) so that we need to use simulated \( z_{t+1} \) from the standardized distribution used in conjunction with the value of \( h_{t+1} \) calculated in the previous step, and hence the distributional uncertainty now enters the simulation.

The following code shows how to obtain equivalence between the 2 simulation methods.

data(sp500ret) spec1 = ugarchspec() fit = ugarchfit(spec1, sp500ret) sim1 = ugarchsim(fit, n.sim = 1, m.sim = 5, startMethod = 'sample', rseed = 10) spec2 = spec1 setfixed(spec2) < - as.list(coef(fit)) sim2 = ugarchpath(spec2, presigma = as.numeric(tail(sigma(fit), 1)), preresiduals = as.numeric(tail(residuals(fit), 1)), prereturns = tail(sp500ret[, 1], 1), n.sim = 1, m.sim = 5, rseed = 10) sigma(sim1) ## [,1] [,2] [,3] [,4] [,5] ## T+1 0.0248 0.0248 0.0248 0.0248 0.0248 # simulated sigma sigma(sim1) == sigma(sim2) ## [,1] [,2] [,3] [,4] [,5] ## T+1 TRUE TRUE TRUE TRUE TRUE # simulated returns fitted(sim1) == fitted(sim2) ## [,1] [,2] [,3] [,4] [,5] ## T+1 TRUE TRUE TRUE TRUE TRUE

This short note described the simulation methods used in the **rugarch** package and its rational based on distributional uncertainty which does not affect \( h_{t+1} \). There is also a method called ugarchboot which can take into account parameter uncertainty (via simulation of the parameter distribution) which then leads to uncertainty in the \( h_{t+1} \) value.