Generating Correlated Random Numbers

by M.A. (Thijs) van den Berg


This article describes common methods that are used in generating correlated random numbers.

Contents

[edit] Generating two sequences of correlated random numbers

Generating two sequences of random numbers with a given correlation is done in two simple steps:

  1. Generate two sequences of uncorrelated normal distributed random numbers X_1, X_2 \,
  2. Define a new sequence Y_1 = \rho X_1 + \sqrt{1-\rho^2} X_2

This new Y_1 \, sequence will have a correlation of \rho\, with the X_1 \, sequence.

[edit] Before and after correlating

Image:RandNormal2C0.png Image:RandNormal2C0.80.png

[edit] Generating multiple (more than two) sequences of correlated random variables

A general way to generate correlated (normal distributed) random numbers -with a given correlation matrix C \,- is done by finding a matrix U such that

U^T U = C \,

Using this matrix, one can generate correlated random numbers R_c\, from uncorrelated numbers R\, by multiplying them with this matrix.

R_c = R U \,

There are multiple way to find such a matrix. The two most common methods are

  1. A Cholesky decomposition of the correlation matrix
  2. An Eigenvector decomposition of the correlation matrix (also known as a spectral decomposition)

[edit] Example using a Cholesky decomposition

We want to generate random numbers for three variables with the following correlations matrix:

C=\begin{bmatrix} 1 & 0.6 & 0.3 \\ 0.6 & 1 & 0.5 \\ 0.3 & 0.5& 1 \end{bmatrix}\,

Doing a Cholesky decomposition on the correlation matrix give the following matrix:

U=\begin{bmatrix} 1 & 0.6 & 0.3 \\ 0 & 0.8 & 0.4 \\ 0 & 0 & 0.866 \end{bmatrix}\,

Depending on the algorithm used for the Cholesky decomposition, you can also get the transposed of this matrix. This does doesn't matter.

We can now transform three uncorrelated random numbers

R=\begin{bmatrix} -0.3999 & -1.6041 & -1.0106 \end{bmatrix}\,

into correlated numbers Rc by multiplying them with the U matrix.

R_c = R U \,

= \begin{bmatrix} -0.3999 & -1.6041 & -1.0106 \end{bmatrix} \begin{bmatrix} 1 & 0.6 & 0.3 \\ 0 & 0.8 & 0.4 \\ 0 & 0 & 0.866 \end{bmatrix}

= \begin{bmatrix} -0.3999 & -1.5232 & -1.6368 \end{bmatrix}

[edit] Example using a Eigenvector decomposition

The correlation matrix C\, has the following eigenvectors

E_1=\begin{bmatrix}0.5700\\ 0.6372\\ 0.5188\end{bmatrix}, E_2=\begin{bmatrix}0.6091\\ 0.0960\\ -0.7872\end{bmatrix}, E_3=\begin{bmatrix}0.5515\\ -0.7647\\ 0.3334\end{bmatrix}

and eigenvalues

\lambda_1=1.9439, \lambda_2=0.7069, \lambda_3=0.3493 \,

If we define

V= E_i Diag( \sqrt{ \lambda_i})
= \begin{bmatrix} 0.5700 & 0.6091 & 0.5515 \\ 0.6372 & 0.0960 & -0.7647 \\ 0.5188 & -0.7872 & 0.3334 \end{bmatrix} \begin{bmatrix} 1.3942 & 0 & 0 \\ 0 & 0.8402 & 0 \\ 0 & & 0.5911 \end{bmatrix}
= \begin{bmatrix} 0.7947 & 0.5121 & 0.3259 \\ 0.8884 & 0.0808 & -0.4520 \\ 0.7232 & -0.6619 & 0.1970 \end{bmatrix}

we get a matrix V\, that has the property

V V^T = C \,

which gives us the final result for the decomposition

U^T U = C, U=V^T \,


We can now transform three uncorrelated random numbers

R=\begin{bmatrix} -0.3999 & -1.6041 & -1.0106 \end{bmatrix}\,

into correlated numbers Rc by multiplying them with the U matrix.

R_c = R U \,

= \begin{bmatrix} -0.3999 & -1.6041 & -1.0106 \end{bmatrix} \begin{bmatrix} 0.7947 & 0.8884 & 0.7232 \\ 0.5121 & 0.0808 & -0.6619 \\ 0.7232 & -0.4520 & 0.1970 \end{bmatrix}
= \begin{bmatrix} -1.4687 & -0.0280 & 0.5734 \end{bmatrix}

[edit] Example Matlab code

The matlab code below decomposes the correlation matrix C into an upper matrix U using a Cholesky decomposition. Next 10 vectors with 3 random normal numbers are multiplied with the upper matrix to generate 10 correlated draws. 


>> C= [ 1.0 0.6 0.3;   0.6 1.0 0.5;   0.3 0.5 1.0]

C =

    1.0000    0.6000    0.3000
    0.6000    1.0000    0.5000
    0.3000    0.5000    1.0000

>> U = chol(C)

U =

    1.0000    0.6000    0.3000
         0    0.8000    0.4000
         0         0    0.8660

>> randn(10,3)*U

ans =

   -0.4326   -0.4089    0.0505
   -1.6656   -0.4187   -1.3665
    0.1253   -0.3955    0.4209
    0.2877    1.9192    2.3656
   -1.1465   -0.7970   -0.9976
    1.1909    0.8057    1.1459
    1.1892    1.5669    1.8695
   -0.0376    0.0248   -1.3678
    0.3273    0.1199   -1.1880
    0.1746   -0.5611    0.2141

[edit] Comments

Name (required):

Website:

Comment:


[edit] correlman said ...

Thanks, excellent!

--correlman 11:19, 4 May 2007 (CEST)

[edit] Jos said ...

Great! Would be nice to have an additional page on fixing invalid correlation matrices.

--Jos 11:36, 4 May 2007 (CEST)

[edit] maquiavelo said ...

Superb! Finally a clear explanation after googling in vain a hundred other articles.

--maquiavelo 18:13, 20 May 2007 (CEST)

[edit] Denisa said ...

I really hope this is the answer for all my work till now ,because i needed such an algorithm to simulate fast fading :)

--Denisa 18:35, 29 May 2007 (CEST)

[edit] Pepe said ...

Clear, concise, has examples - what to need more? Thumbs up!

--Pepe 12:29, 30 May 2007 (CEST)

[edit] harry said ...

Amazing n cool . I can't believe its soo.. easy Finally a relief from useless notation found elsewhere..........

--harry 03:34, 16 June 2007 (CEST)

[edit] Borat said ...

Very Nice

--Borat 16:07, 14 August 2007 (CEST)

[edit] Zhang said ...

Good job!

--Zhang 19:10, 17 August 2007 (CEST)

[edit] Rama said ...

thanks..it is great help for me :)

--Rama 10:47, 19 August 2007 (CEST)

[edit] fnord said ...

How do you extend the method to the case where one or more of the cross-correlations are negative?

--fnord 17:22, 7 September 2007 (CEST)

[edit] sitmo said ...

Negative correlations is not different. The Cholesky decomposition of C to get U will also work with negative correlations (as long as the correlation matrix is 'valid/positive definite')

--sitmo 09:45, 12 September 2007 (CEST)

[edit] sw3quant said ...

I used this for Monte Carlo, but needed to modify it. I used the covar matrix (not the corellMatrix). Multiply by the transpose of Q, not by Q. we checked through this at work using R and generating a vector of 1M gaussian variables. We used teh MAPCK implementation in C# which is pretty good.

--sw3quant 11:18, 12 October 2007 (CEST)

[edit] Tolga said ...

my understanding is that, in both methods, we could really have used ANY orthogonal matrix U for the right-multiplications. After all, the correlations of the correlated random variables that are produced are actually given in the components of the orthonormal eigenvectors, not by the elements of the initial symmetric correlation matrix that yield those orthonormal eigenvectors.

Am i right??

But if so, then what is the point of starting out with a correlation matrix?

I am also wondering why/whether it is necessary for the squares of the correlations to add up to 1 when combining the uncorrelated random sequences to produce the correlated ones.


--Tolga 02:44, 24 December 2007 (CET)

[edit] Ricardo said ...

First of all, thank you for this straight foward solution.

I need something quite similar, but with log-normal random variables.

My idea is to generate N(0,1) corralated variables and then "transform" them to LN (miu,sigma) correlated variables.

The theory presented remains valid?

Thank you!

--Ricardo 14:39, 10 January 2008 (CET)

[edit] admin said ...

Ricardo, Yes is does. If you take Y=exp( N(0,1) ) then Y is lognormal distributes and will take values between 0 and +infinity


--Admin 15:49, 12 March 2008 (CET)

[edit] chakazoulou said ...

this correlation don't give rho

--chakazoulou 11:49, 29 March 2008 (CET)

[edit] chakazoulou said ...

except if var(X1) = 1 = Var(x2)

--chakazoulou 13:16, 29 March 2008 (CET)

[edit] Sofia said ...

Thank you admin for your very useful solution. I have a question. Does this formula apply to exponentially distributed variables as well? That is, if X1 and X2 are exp. distr. variables, then Y1 is correlated with X1 with coefficient rho? Thanks again!

--Sofia 11:01, 6 May 2008 (CEST)

[edit] Ricardo said ...

Can you add a section on how to do this when the correlation matrix has complex values instead of all real? Thanks.


--Ricardo 22:38, 8 May 2008 (CEST)

[edit] Rob said ...

Great. Please your help on how to generate multiple sequences of correlated random vaiables not normally distributed for instance Gamma distribed?

--Rob 05:56, 15 May 2008 (CEST)

Retrieved from "http://sitmo.com/doc/Generating_Correlated_Random_Numbers"