4
$\begingroup$

Given a distribution A with a mean of $\mu_1$ and standard deviation of $\sigma_1$, how can I generate:

  • Distribution B with a mean of $\mu_2$ and standard deviation of $\sigma_2$ and a correlation of $X_1$ with distribution A
  • Distribution C with a mean of $\mu_3$ and standard deviation of $\sigma_3$ and a correlation of $X_2$ with distribution B and $X_3$ with distribution A

Can someone please tell me if this even makes sense? My naive approach was the following:

  1. Generate A with the given parameters
  2. Generate B with the given parameters and then see if the generated values have the specified correlation with A. If not, regenerate B until this correlation is achieved.
  3. Generate C using the approach in Step 2.

However, I am not quite sure if this approach will terminate. Is there a better way to achieve this? I'd love to see an example in R.

$\endgroup$
  • 2
    $\begingroup$ There are restrictions on the correlations that are achievable. The restrictions depend on the correlations and on the actual distributions. If you have additional requirements, such as that distributions A, B, and C all be in the same location-scale family, that will considerably simplify matters, so please provide any such details you can. It is possible your question is already answered here, including the R code: see stats.stackexchange.com/q/24257 for instance. $\endgroup$ – whuber Jan 14 '14 at 16:32
  • $\begingroup$ @whuber: Thank you. In my case, there is evidence pointing to A following Gaussian but not always. Therefore, while my original intent was to understand this in the case of Gaussians, I am interested in understanding when it is not Gaussians. I am not an expert in this so I'm not quite sure I understand the second half of your question concerning location-scale families. I read this on Wiki and based on my understanding, it looks my answer is "yes" but I'm not quite sure how this will affect the solution. Can you kindly shed more light on how this will affect the solution? $\endgroup$ – Legend Jan 14 '14 at 18:34
  • $\begingroup$ For an example of the limitations on correlations between lognormal distributions, please see the nice answer at stats.stackexchange.com/questions/41734. A partial answer to your general question appears at stats.stackexchange.com/questions/62146/…. $\endgroup$ – whuber Jan 14 '14 at 19:17
  • $\begingroup$ @whuber: Thank you very much for the references. I'll look at them in detail. $\endgroup$ – Legend Jan 14 '14 at 19:19
8
$\begingroup$

If you mean the individual distributions are Gaussian, then sampling from a multivariate normal with mean vector $\mathbf{\mu}$ and covariance matrix $\mathbf{\Sigma}$ will generate such data.

Here is an R example using the function mvrnorm() in package MASS (which comes with R):

## means of individual distributions
mu1 <- 5
mu2 <- 10
mu3 <- 0
## variance
sigma1 <- 5
sigma2 <- 1
sigma3 <- 0.5
## Correlations
X1 <- 0.5
X2 <- 0.1
X3 <- 0.8

## load package
require("MASS")

We need to supply n, the number of values from each distribution, mu the mean vector, and Sigma the covariance matrix. In the code below I form these from the scalars entered above.

set.seed(1)
dat <- mvrnorm(100, mu = c(mu1, mu2, mu3),
               Sigma = matrix(c(sigma1, X1    ,     X3,
                                X1    , sigma2,     X2,
                                X3    , X2    , sigma3),
                              ncol = 3, byrow = TRUE),
               empirical = TRUE)

I used empirical = TRUE to specify empirical not population parameters for $\mathbf{\mu}$ and $\mathbf{\Sigma}$. This results in the covariance matrix of dat having exactly the values we specified:

R> cov(dat)
     [,1] [,2] [,3]
[1,]  5.0  0.5  0.8
[2,]  0.5  1.0  0.1
[3,]  0.8  0.1  0.5

as do the column means:

R> colMeans(dat)
[1]  5.000e+00  1.000e+01 -8.882e-18

If you use the default, empirical = FALSE, then you get random samples from a population which will have different sample mean vector and sample covariance matrix from the specified one as you have only seen n examples from that larger population:

set.seed(1)
dat2 <- mvrnorm(100, mu = c(mu1, mu2, mu3),
                Sigma = matrix(c(sigma1, X1    ,     X3,
                                 X1    , sigma2,     X2,
                                 X3    , X2    , sigma3),
                               ncol = 3, byrow = TRUE))

R> cov(dat2)
       [,1]    [,2]    [,3]
[1,] 4.0441 0.39858 0.61120
[2,] 0.3986 0.91110 0.04842
[3,] 0.6112 0.04842 0.48782
R> colMeans(dat2)
[1]  5.24138 10.06668  0.02448
$\endgroup$
  • 5
    $\begingroup$ This is a well-written and clear answer. Please note, though, that if the question being posed truly is about Gaussian distributions only, then it is a duplicate of stats.stackexchange.com/questions/24257/… and eventually will be closed as such (and your answer can be merged with that thread). If it is not just about Gaussians, then the issue becomes much more complicated and needs a different kind of answer. $\endgroup$ – whuber Jan 14 '14 at 16:49
  • $\begingroup$ @whuber Thanks for pointing out the duplicate; I only noticed that in your comment above after I was mostly done with this answer. If they do want Gaussians then please do merge this one with the existing Q&A. $\endgroup$ – Gavin Simpson Jan 14 '14 at 16:51
  • $\begingroup$ +1 Thank you Gavin. This example is just perfect. I just asked @whuber another question and am waiting for his response as well for the case when it is non-Gaussian. It would be greatly benefit me if you can append your answers with any other details regarding this. $\endgroup$ – Legend Jan 14 '14 at 18:35
  • 1
    $\begingroup$ @Legend Thanks. For non-Gaussian's, the only approach i am aware of is that if Iman and Conover which works to achieved rank-order-correlated distributions and which works for any distribution IIRC. There is an R implementation but I'll have to look that up. $\endgroup$ – Gavin Simpson Jan 14 '14 at 18:55
  • $\begingroup$ @Legend I found it, it is the cornode() function in package mc2d $\endgroup$ – Gavin Simpson Jan 14 '14 at 18:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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