I know that for the $2$-dimensional case: given a correlation $\rho$ you can generate the first and second values, $ X_1 $ and $X_2$, from the standard normal distribution. Then from there make $X_3$ a linear combination of the two $X_3 = \rho X_1 + \sqrt{1-\rho^2}\,X_2$ then take $$ Y_1 = \mu_1 + \sigma_1 X_1, \quad Y_2 = \mu_2 + \sigma_2 X_3$$

So that now $Y_1$ and $Y_2$ have correlation $\rho$.

How would this be scaled to $n$ variables? With the condition that the end variables satisfy a given correlation matrix? I'm guessing at least n variables will need to be generated then a reassignment through a linear combination of them all will be required... but I'm not sure how to approach it.

  • $\begingroup$ If you happen to be using Matlab, see my function here and/or if you have the Statistics Toolbox, see the mvnrnd function. $\endgroup$
    – horchler
    Commented Jul 17, 2013 at 21:09
  • $\begingroup$ @horchler Actually using R, almost done writing the script... but I'll check yours out seeing as I'm getting a close but not correct result. $\endgroup$ Commented Jul 17, 2013 at 21:24
  • $\begingroup$ R doesn't appear to have a builtin cholcov function (just chol) so you'll just need to make sure that you actually use correlation matrices (ones on the diagonal) rather than covariance matrices to meet the positive semi-definite criterion required for Cholesky decomposition. You can use R's cov2cor to convert if needed. $\endgroup$
    – horchler
    Commented Jul 17, 2013 at 21:38
  • $\begingroup$ I couldn't find the appropriate library so I started making my own algorithm based on the link provided by @JosephK in his answer $\endgroup$ Commented Jul 17, 2013 at 21:40

2 Answers 2


If you need to generate $n$ correlated Gaussian distributed random variables $$ \bf Y \sim \mathcal N(\bf \mu, \Sigma) $$ where $\textbf{Y} = (Y_1,\dots,Y_n)$ is the vector you want to simulate, $\mu =(\mu_1,\dots, \mu_n)$ the vector of means and $\Sigma$ the given covariance matrix,

  1. you first need to simulate a vector of uncorrelated Gaussian random variables, $\bf Z $
  2. then find a square root of $\Sigma$, i.e. a matrix $\bf C$ such that $\bf C \bf C^\intercal = \Sigma$.

Your target vector is given by $$ \bf Y = \bf \mu + \bf C \bf Z. $$

A popular choice to calculate $\bf C$ is the Cholesky decomposition.

  • 2
    $\begingroup$ Exactly what I was looking for. Thanks a lot $\endgroup$ Commented Jul 17, 2013 at 20:39
  • 1
    $\begingroup$ It's not altogether clear that this is correct: The question says "the correlation between the observed outcomes will be the same as in the matrix". If you choose from a multivariate normal with a certain correlation, generally the sample correlation will not equal the population correlation. If the idea is to make the sample correlation equal to the specified value, then one is sampling from the conditional distribution given that value. That can be done; I've written an algorithm for it in case $n=2$. $\endgroup$ Commented Jul 17, 2013 at 21:45
  • 2
    $\begingroup$ @MichaelHardy Where is this algorithm? Also, I'm primarily interested in cases for larger $n$, but most specifically I care about the case where $n=3$. $\endgroup$ Commented Jul 18, 2013 at 14:54
  • 2
    $\begingroup$ All clear, but the statement "then find a square root of Σ, i.e. a matrix C such that CC⊺=Σ" may be confusing as C is not strictly the square root of Σ $\endgroup$
    – dayum
    Commented Sep 10, 2018 at 2:05


I wanted to add that comment to the previous answer :

Why is it true ?

Define $Y = \mu + C Z$, and recall that each component of $Z$ is i.i.d. normal :

$$ \Cov( Y ) = \Cov(CZ + \mu ) = \Cov(CZ ) = C \cdot \Cov (Z) \cdot C^T = C I C^T = \Sigma $$ which is the desired result : $$ \implies Y \sim N (\mu, \Sigma) $$

  • 3
    $\begingroup$ Wow, thank you for this wonderful explanation. This is really helped me! $\endgroup$ Commented Oct 16, 2021 at 22:13

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .