I want to simulate some exchange rates with a DCC GARCH. I know the package rmgarch but I want to code the simulation my self. The following are the main equations of the model:

$r_t = a_t$

$a_t = H^{1/2}_{t} Z_t$

$H_t = D_t R_t D_t$

$R_t = Q^{*-1}_t Q_t Q^{*-1}_t$

$Q_t = (1-a-b)\bar{Q} + a\epsilon_{t-1}\epsilon'_{t-1} + b Q_{t-1}$

My question is: if I all ready have the univariate GARCH parameters for $D_t$, $a$ and $b$ and $\bar{Q}$ how do I get $\epsilon_{t-1}$ to perform the simulations? and what are the steps to perform it.

  • 1
    $\begingroup$ you set the first errors term to zero (their mean) and then derive the other ones by fitting data to your model (recursively for each $t$ - by going backward from the origin) $\endgroup$
    – Malick
    Commented Jul 14, 2018 at 14:39
  • $\begingroup$ @Malick I don't think is like that because with what data will i compare it if it's a simulation? $\endgroup$ Commented Jul 14, 2018 at 16:46
  • $\begingroup$ @Malick but isn't $Z_t$ the variable that follows a N(0,1)? $\endgroup$ Commented Jul 16, 2018 at 16:25
  • $\begingroup$ Sorry Alejandro, I have deleted my answer as it was wrong, I'll try to come up with a proper answer. $\endgroup$
    – Malick
    Commented Jul 18, 2018 at 10:58
  • $\begingroup$ @Malick thanks anyway if you manage to come up with the right answer i would appreciated $\endgroup$ Commented Jul 18, 2018 at 15:17

1 Answer 1


Simulating a DCC-GARCH(1,1) model

Given that you already have a given set of proper defined parameters for the DCC-GARCH model, the standardized residuals $\varepsilon_{t-1}$ are recovered from the univariate GARCH models and fed into the DCC structure to yield your simulated correlation dynamics. The simulation of a DCC-GARCH model can be defined in 3 overall steps:

  1. Simulate univariate GARCH dynamics using a filter algorithm for pre-defined parameters and recover the standardized residuals. After that, construct a general function that simulates $d$ GARCH models with different sets of parameters (or the same set of parameters).

  2. Use the standardized residuals from (1) to simulate the DCC dynamics and recover the correlation structure. From construction, we know that $\varepsilon_{t-1} \overset{iid}{\sim} N(0,1)$ and thus you can simulate the DCC dynamics and back out the correlation without the GARCH component.

  3. Recover the covariance matrix via the univariate variances and the correlation structure using the following formula, $\Sigma_t=D_t^{1/2}R_t D_t^{1/2}$, with $D_t = diag(\sigma_{i,t}^2)$

Throughout the answer, the code will be written in R. Moreover, the code will be non-optimized.

1. Simulating multiple univariate GARCH models

Formally let us define the (demeaned) returns in a univariate setting,

\begin{equation} r_t = \sigma_t \epsilon_t \iff \frac{r_t}{\sigma_t} = \epsilon_t, \end{equation}

where $\epsilon_t \overset{iid}{\sim}N(0,1)$ are univariate standardized residuals. This further implies that $r_t \sim N(0, \sigma_t^2)$, which we will be using in the simulation. The variance dynamics of the univariate returns follows a GARCH(1,1) on the form:

\begin{equation} \sigma^2_t = \omega + \alpha r_{t-1}^2 + \beta \sigma^2_{t-1}. \end{equation}

The way to simulate a GARCH model we initialize the simulation by sampling $r_1 \sim N(0,\sigma^2_1)$, where $\sigma^2_1 = \frac{\omega}{1-\alpha - \beta}$ is the unconditional variance. After that, we recusively update the variance dynamics and simulate the returns using the rnorm() function.

The function is presented below under the name GARCHSim and has been extensively commented. For ease of implementation, we Construct a function (GARCHmultisim) that performs multiple simulations in one function call. Here, the parameters variable (params) will be a $d \times 3$ matrix with $d$ denoting the amount of GARCH models to simulate. The function outputs a list of two matrices with dimensions $sims \times d$, where one contains the standardized residuals and the other contains the variance.

2 + 3. Simulating the DCC structure using standardized residuals

Constructing the filter algorithm for the DCC structure is, more or less, the same procedure as with the GARCH model. We need to update $Q_t$ for each time-point and then reconstruct the correlation matrix following $R_t = \widetilde{Q}_t^{-1/2} Q_t \widetilde{Q}_t^{-1/2}$ with $\widetilde{Q}_t$ containing the diagonal elements of $Q_t$. Then we can reconstruct the covariance matrix, $\Sigma_t$, using the relation in point 3.

We further model $Q_t$ following a GARCH(1,1)-like structure where we impose weak-stationarity (ie. $a+b<1$) to ensure that $Q_t$ is positive definite and implying that $R_t$ is positive definite: \begin{equation*} Q_t = \bar{Q} (1-a-b) + a (\eta_{t-1}\eta_{t-1}^\intercal) + b Q_{t-1}, \end{equation*}

where $\bar{Q}$ is the unconditional covariance matrix of the multivariate standardized residuals $\eta_t$. We can recover both parameters following

\begin{align*} \bar{Q} &= \frac{1}{T} \sum_{t=1}^{T} \eta_t \eta_t^\intercal\\ \eta_t &= D_t^{-\frac{1}{2}} r_t \sim N(0, R_t), \end{align*} which we will utilize in our simulation procedure. The function has been named DCCsim and can be found in the end of the code snippet.

Graphical illustration

I have provided a graphical illustration of a simulated bivariate DCC-GARCH(1,1) model with arbitrary parameters given by:

$$\omega_{1}=0.02, \quad \alpha_1 = 0.08, \quad \beta_1 = 0.89, \quad \omega_{2}=0.02, \quad \alpha_2 = 0.05, \quad \beta_2 = 0.94, \quad a=0.1, \quad b=0.89,$$

with $\omega_i, \: \alpha_i, \: \beta_i$ for $i=1,2$ being the parameters for the first and second GARCH process, which are also depicted below as the red and green one respectively.


Some comments

  • The correlation structure in the DCC model is dependent on the residuals of the GARCH models. Essentially, the univariate GARCH models do not capture the dependence structure and thus it will be contained within the residuals. However, the residuals follows a univariate standard Gaussian distribution (and are iid) thus having no correlation structure whatsoever. This implies that the unconditional correlation tends toward zero as $sims \rightarrow \infty$.

  • You can also simulate multivariate normal variates with dynamic correlations in the DCC model, where the correlations are estimated from real world stock prices (that is, you can use the mvrnorm function in R with a changing covariance matrix, based on different correlations), or you can try to use alternative distributions for the residuals.

  • The Constant Conditional Correlation (CCC) model is recovered by reconstructing the covariance matrix (DCC_COVAR) where you replace acor with the unconditional correlation.

  • Any extensions in the univariate GARCH models is valid due to the separation of the correlation- and volatility dynamics.

The above implementation provides a basic simulation framework. Consequently, there might be more efficient ways to construct realistic simulated correlations within the DCC structure or completely alternative ways to simulate this model. I hope this atleast provides some insight for you and future readers.

Code appendix

The code is presented below.

# Univariate GARCH simulation
GARCHSim <- function(omega, alpha, beta, sims) {
  #sims (integer): Number of simulations performed. 
  ## initialize the vector of simulated returns and variances
  ret = numeric(sims)
  sigma2 = numeric(sims)
  ## initialize the variance at time t = 1 with the unconditional variance.
  sigma2[1] = omega/(1.0 - alpha - beta)

  ## sample the first observations
  ret[1] = rnorm(1, mean = 0, sd = sqrt(sigma2[1]))
  ## loop over sims.
  for (t in 2:sims) {
    #update volatility
    sigma2[t] = omega + alpha * ret[t - 1]^2 + beta * sigma2[t - 1]

    # sample new observarion
    ret[t] = rnorm(1, mean = 0, sd = sqrt(sigma2[t]))
  ## we return a list with three components: the sampled returns, the variance, 
  ## and the standardized residuals. 
  lOut = list()
  lOut[["ret"]] = ret
  lOut[["sigma2"]] = sigma2

  # gives standardized residuals for demeaned returns. 
  lOut[["stdres"]] = ret/sqrt(sigma2)

# Multiple simulations of the GARCH model.
GARCHmultisim <- function(params, sims){

  multisim_res <- matrix(0L, nrow=sims, ncol = nrow(params))
  multisim_sigma <- matrix(0L, nrow=sims, ncol = nrow(params))

  # rows in parameter matrix determines amount of garch simulations.
  for(i in 1:nrow(params)){

    # matrix() transforms it into a vector. 
    multisim_res[, i] <- matrix(GARCHSim(params[i,1], params[i,2], params[i,3], sims)[["stdres"]])
    multisim_sigma[, i] <- matrix(GARCHSim(params[i,1], params[i,2], params[i,3], sims)[["sigma2"]])

  lOut = list()

  lOut[["res"]] = multisim_res
  lOut[["sigma"]] = multisim_sigma


# Simulation for the Dynamic Conditional Correlations (DCC).
DCCsim <- function(GARCH_params, a, b, sims){

  stdres <- GARCHmultisim(GARCH_params, sims)[["res"]]

  Gsigs <- GARCHmultisim(GARCH_params, sims)[["sigma"]]

  mQ = cor(stdres)
  iN = ncol(stdres)
  iT = nrow(stdres)
  # initialize the array for the correlations
  aCor = array(0, dim = c(iN, iN, iT))

  # initialize the array for the Q matrices
  aQ = array(0, dim = c(iN, iN, iT))
  ## initialization at the unconditional cor
  aCor[,, 1] = mQ
  aQ[,,1] = mQ

  # main loop
  for (t in 2:iT) {
    # update the Q matrix
    aQ[,, t] = mQ * (1 - a - b) + a * t(stdres[t - 1, , drop = FALSE]) %*%
    stdres[t - 1, , drop = FALSE] + b * aQ[,,t - 1]
    ## Compute the correlation as Q_tilde^{-1/2} Q Q_tilde^{-1/2} 
    aCor[,, t] = diag(sqrt(1/diag(aQ[,, t]))) %*% aQ[,, t] %*% 
    diag(sqrt(1/diag(aQ[,, t]))) 

    #Getting univariate variance matrix D_t:

    msigma <- array(0L, dim = c(ncol(Gsigs), ncol(Gsigs), nrow(Gsigs)))

    for(i in 1:ncol(Gsigs)){

      msigma[i, i, ] <- Gsigs[, i]


    # constructing the covariance matrix
    DCC_COVAR <- array(0, dim = c(ncol(Gsigs), ncol(Gsigs), nrow(Gsigs)))

    DCC_returns <- matrix(0, nrow = nrow(Gsigs), ncol = ncol(Gsigs))
    for(i in 1:sims){

      DCC_COVAR[,,i] <- msigma[,,i]^0.5 %*% aCor[,,i] %*% msigma[,,i]^0.5

      # Getting returns:
      DCC_returns[i,] <- mvrnorm(1, rep(0,ncol(Gsigs)), DCC_COVAR[,,i])


  lOut = list()

  # correlations
  lOut[["aCor"]] = aCor
  # univariate variances
  lOut[["Gsigs"]] = Gsigs
  # covariances
  # returns
  lOut[["DCC_returns"]] = DCC_returns

  • $\begingroup$ How is it done for the CCC-GARCH model? @Pleb $\endgroup$
    – Adrian
    Commented Jun 21, 2023 at 16:56
  • $\begingroup$ @Adrian As specified under my "some comments" section in my answer, you recover the CCC model by setting acor equal to the unconditional correlation. This can be done by setting $a=0$ and $b=0$ in the DCC model, so the unconditional correlation is recovered from the factorization $R = \widetilde{Q}^{-1/2} Q \widetilde{Q}^{-1/2}$. $\endgroup$
    – Pleb
    Commented Jun 21, 2023 at 18:15
  • 1
    $\begingroup$ @THAT'SMYQUANTMYQUANTITATIV Yes, only the dependency structure is captured (used) in the DCC-part. The DCC model is regarded as the $Q_t$ model dynamics. To estimate $a$ and $b$, you only need your standardised residuals (not the original $\Sigma$). When you have estimated $a$ and $b$ you can now forecast $Q_{t+1|t}$, which you can use to recover $R_t$ (as seen above). Since you have estimated all your $N$ univariate GARCH models you can recover $D_{t+1|t}$ and hence, you can recover your 1-step ahead covariance matrix, $\Sigma_{t+1|t}$, following the fomula observed in step 3 of my answer. $\endgroup$
    – Pleb
    Commented May 5 at 13:52
  • 1
    $\begingroup$ @THAT'SMYQUANTMYQUANTITATIV Yes, you are right, there will be no dependence when simulating iid GARCH models, as also highlighted in my "some comments" section of my answer. This is also the reason I advocate for simulating multivariate correlated returns using the mvrnorm (second bullet-point in "some comments"), that you can either feed into the GARCH simulations (instead of using rnorm) and further into the DCC structure (through the GARCH residuals) or you can back out the residuals by standardizing the simulated mult. returns directly and feed them into the DCC model. $\endgroup$
    – Pleb
    Commented May 6 at 6:34
  • 1
    $\begingroup$ @THAT'SMYQUANTMYQUANTITATIV No, not really. However, one risk-based application for testing mult. GARCH models, is to estimate risk-spillover effects (from markets to certain companies or from markets to markets). A paper titled CoVaR (2011) provides an introduction to their risk-measure and they also give ways to estimate the it using GARCH models (see the appendix of the paper). I thought that would be a good applied way to test your models (how do they fare in terms of risk-spillover estimation), but I also provided it as inspiration. Good luck :-) $\endgroup$
    – Pleb
    Commented May 7 at 10:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.