Is it true (see here, footnote 2, p.22 / p.14, without proof) that we can obtain two discretized brownian motions $W_t^1, W_t^2$ with correlation $\rho$ by doing

$$d W_t^1 \sim \mathcal N(0,\sqrt{dt})$$

$$d W_t^2 = \rho\, d W_t^1 + \sqrt{1-\rho^2} dZ_t$$ with $dZ_t \sim \mathcal N(0,\sqrt{dt})$ ?

If this is true, we could easily simulate them in Python by doing:

enter image description here

import numpy as np
import matplotlib.pyplot as plt

n = 1000000; dt = 0.01; rho = 0.8
dW1 = np.random.normal(0, np.sqrt(dt), n)
dW2 = rho * dW1 + np.sqrt(1 - rho **2) * np.random.normal(0, np.sqrt(dt), n)
W1 = np.cumsum(dW1)
W2 = np.cumsum(dW2)

Is this correct?


First you need to correct the formula to: $$ W_t^2 = \rho W_t^1 + \sqrt{1-\rho^2} Z_t, $$ where $Z_t$ is a BM independent of $W_t^1$ If you calculate the variance and the covariance, then you see that it is true: $$ V[W_t^1] = t $$ and $$ V[W_t^2] = \rho^2 V[W_t^1] + (1-\rho^2) V[Z_t] = \rho^2 t + (1-\rho^2) t = t, $$ which is the desired variance.

For the covariance you get $$ Cov[W_t^1,W_t^2] = \rho Cov[W_t^1, W_t^1] + \sqrt{1-\rho^2} Cov[W_t^1, Z_t] $$ which gives (because $Cov[W_t^1, W_t^1] = V[W_t^1]$ and by indpenence of $W_t^1$ and $Z_t$: $$ Cov[W_t^1,W_t^2] = \rho t + 0, $$ and noting that $\sqrt{V[W_t^1]} = \sqrt{V[W_t^2]} = \sqrt{t}$ we get $$ Cor[W_t^1,W_t^2] = \frac{Cov[W_t^1,W_t^2]}{\sqrt{V[W_t^1]} \sqrt{V[W_t^2]} } = \frac{\rho t}{t} = \rho. $$

  • $\begingroup$ Thanks for this proof @Richard. Do you think my simulation code is ok? It seems ok. Another thing : I think your first equation is a "correlation" link between brownian motions $W_t^1$ and $W_t^2$, whereas my equation was a link between brownian motion increments $d W_t^1$ and $d W_t^2$. By summing the increments (cumulative sum), it seems that it would leed to the same formula as yours. (of course when dealing with discretized brownian motion) $\endgroup$
    – Basj
    Feb 23 '16 at 20:50
  • $\begingroup$ Hi, and no, I don't think so. You need the square-root because constant multiplicators enter variance with their square. and no matter whether we speak of increments or the processes it is variance and covariance that matters. Writing $W_t = W_t - W_0$ with $W_0 =0$ you can see that it does not matter whether we speak of increments or the process. In general we can write $W_t = W_s + W_t-W_s$ and $W_t-W_s$ is independent of $W_s$ for $t>s$. $\endgroup$
    – Ric
    Feb 24 '16 at 7:04
  • $\begingroup$ not that in the answer below the second parameter in $N(x,y)$ is the variance and not the volatility. $\endgroup$
    – Ric
    Feb 24 '16 at 7:04
  • $\begingroup$ Oh yes sure, it was a mistake, I forgot the square root! Edited in my question. Now our 2 equations are equivalent (yours about brownian motion, mine about increments $d W_t$) $\endgroup$
    – Basj
    Feb 24 '16 at 8:46
  • $\begingroup$ Yes, but in your second equastion you have a $Z$ which is just one random variable. MAybe you should put a $dZ_t$ there. Then if you integrate you have the same as I have. $\endgroup$
    – Ric
    Feb 24 '16 at 8:56

Here is the general approach you can follow to generate two correlated random variables. Let's suppose, X and Y are two random variable, such that: $$X \sim N(\mu_1, \sigma_1^2)$$ $$Y \sim N(\mu_2, \sigma_2^2)$$ and $$cor(X,Y)=\rho$$ Now consider: $y=bx + e_i$, where $x$ $(=\frac{X-\mu_1}{\sigma_1}$) and $y$ $(=\frac{Y-\mu_2}{\sigma_2}$) both follow standard normal distribution , such that $cor(x,y)=\rho. $ For standard normal variate, $b= \rho$. So we have:

$$y=\rho x + e_i$$

Now, here is the algorithm, you can follow:

1) Generate $n$ standard normal variate for $x$.

2) Since, $e_i \sim N(0, 1-\rho^2)$. So generates $n$ normal variate as $e_i$ from normal distribution with mean 0 and variance $1-\rho^2$.

3) Get $y=\rho x + e_i$.

4) Convert your standard normal numbers back to Normal (remember correlation is independent of change of origin and scale)

R code for simulating correlated GBM:

corGBM <- function(n, r, t=1/365, plot=TRUE) {
#n is number of samples 
#r is correlation
#t is tick step
se <- sqrt(1 - r^2) #standard deviation of error
e <- rnorm(n, mean=0, sd=se)
y <- r*x + e

X <- cumsum(x* sqrt(t))
Y <- cumsum(y* sqrt(t))
Max <- max(c(X,Y))
Min <- min(c(X,Y))

if(plot) {
plot(X, type="l", ylim=c(Min, Max))
lines(Y, col="blue")

#sample result
[1] 0.8523341

enter image description here

  • $\begingroup$ There are some problems in your R code I think : a) you aren't generating brownian motion but only increments. You have to cumsum them to get brownian motion. b) you define r2 but you don't use it c) even if both notations work, why writing r ** 2 and then r^2? d) you don't call the function correlatedvalue. Can you include code to plot the two correlated brownian motions? (More generally I think your method is more or less the same as the one I used in Python in my original post, but it's interesting to have it in R for learning purpose, thanks!) $\endgroup$
    – Basj
    Feb 23 '16 at 20:44
  • $\begingroup$ @Basj Thanks for your valuable inputs. I have added R code for simulating two geometric Brownian motion. I hope it would help you. $\endgroup$
    – Neeraj
    Feb 24 '16 at 15:24

Your Answer

