# Back of the Envelope

Observations on the Theory and Empirics of Mathematical Finance

## Simulating Correlated Stochastic Differential Equations (or How to Simulate Heston Stochastic Volatility Model)

I notice that students new to computational finance often make mistakes in simulating correlated Brownian motion paths. Here is a ready reckoner.

Let’s take the example of generating paths for asset prices using the Heston stochastic volatility model:

\begin{aligned} dS &= rSdt + \sqrt{\nu} S dW_1 \\ d \nu &= \kappa (\theta - \nu) dt + \sigma \sqrt{\nu} dW_2 \end{aligned}

where $\nu$ is the instantaneous variance of asset returns, and the increment in Brownian motions $dW_1$ and $dW_2$ are correlated with correlation coefficient $\rho$, i.e. $\langle dW_1, dW_2 \rangle = \rho dt$.

The simplest way to generate paths $S$ and $\nu$ is to use the Euler discretization (there are better methods available of course, for Heston in particular) as:

\begin{aligned} ln S(t + \Delta t) &= ln S(t) + r \Delta t - \frac{1}{2} \nu (t) \Delta t + \sqrt{\nu (t)} \epsilon_1 \sqrt{\Delta t} \\ \nu (t + \Delta t) &=\nu (t) + \kappa (\theta - \nu (t) ) \Delta t + \sigma \sqrt{\nu (t)} \epsilon_2 \sqrt{\Delta t} \end{aligned}

where $\epsilon_1$ and $\epsilon_2$ are standard Gaussian random variables with correlation $\rho$

To generate correlated standard Gaussian random variables, i.e. $\epsilon_1$ and $\epsilon_2$, the most popular method is to use what is called the Cholesky decomposition. Given two uncorrelated standard Gaussian random variables $Z_1$ and $Z_2$ (easily done both in Excel and in R), Cholesky decomposition can be used to generate $\epsilon_1$ and $\epsilon_2$ as:

\begin{aligned} \epsilon_1 &= Z_1 \\ \epsilon_2 &= \rho Z_1 + \sqrt{1 - \rho^2} Z_2 \end{aligned}

If, God forbid, your job requires simulating three correlated stochastic differential equations, say when you are using a Double Heston or a Double Lognormal model, then you would need to simulate three jointly correlated Gaussian random variables.

In that case if the correlation structure is $\langle dW_1, dW_2 \rangle = \rho_{12} dt$$\langle dW_2, dW_3 \rangle = \rho_{23} dt$ and $\langle dW_1, dW_3 \rangle = \rho_{13} dt$, then given three uncorrelated Gaussian random variables $Z_1$, $Z_2$ and $Z_3$, one could use Cholesky decomposition to generate generate $\epsilon_1$, $\epsilon_2$ and $\epsilon_3$ as:

\begin{aligned} \epsilon_1 &= Z_1 \\ \epsilon_2 &= \rho_{12} Z_1 + \sqrt{1 - \rho^2_{12}} Z_2 \\ \epsilon_3 &=\rho_{13} Z_1 + \rho^*_{23} Z_2 + \rho^*_{33} Z_3\end{aligned}

where

$\displaystyle \rho^*_{23} = \frac{\rho_{23} - \rho_{12} \rho_{13}}{\sqrt{1 - \rho^2_{12}}}$ and $\displaystyle \rho^*_{33} = \sqrt{\frac{1 - \rho^2_{12} - \rho^2_{23} - \rho^2_{13} + 2 \rho_{12} \rho_{23} \rho_{13}}{1 - \rho^2_{12}}}$.