Paella1
Paella1

Reputation: 15

Generating random variables that are correlated with one vector but not between each other

I have a vector s1 containing normally distributed random variables. I want to generate 4 more normally distributed random vectors, each of which has its own correlation with s1 and its own variance. Let's call them s2 to s5.

If I use a mvrnorm() with a covariance matrix, I have to designate the covariances between s1 and each of the other vectors, which is fine. But I also have to designate covariances between each of the other vectors (E.g. between s2 & s3), which is not fine. I will end up with a correlation between s2 & s3, and there is no reason why there should be one.

How can I do generate s2 to s5 with designated (and different) standard deviations and designated covariances with s1, WITHOUT forcing correlations between s2 to s5?

edit: Here's the covariance matrix AFTER setting rho(3,2) to zero

           [,1]      [,2]       [,3]
[1,] 0.00022500 0.0002625 0.00044625
[2,] 0.00026250 0.0006250 0.00000000
[3,] 0.00044625 0.0000000 0.00122500

Upvotes: 0

Views: 324

Answers (1)

thothal
thothal

Reputation: 20329

Just set the corresponding elements in your covariance matrix to 0:

library(MASS)
set.seed(1)
(sig <- matrix(c(5, .5, .8, .5, 1, 0, .8, 0, .5), 3))
#      [,1] [,2] [,3]
# [1,]  5.0  0.5  0.8
# [2,]  0.5  1.0  0.0  ## <- 0 = covariance between s2 and s3
# [3,]  0.8  0.0  0.5
x <- mvrnorm(1e5, rep(0, 3), sig)
cov(x)
#           [,1]         [,2]         [,3]
# [1,] 5.0356870 0.5100643820 0.8004814044
# [2,] 0.5100644 1.0042540190 0.0008037978
# [3,] 0.8004814 0.0008037978 0.4972328657

## with empirical = TRUE you can force the cov matrix to match exactly sig
cov(mvrnorm(1e5, rep(0, 3), sig, empirical = TRUE))
#      [,1]          [,2]          [,3]
# [1,]  5.0  5.000000e-01  8.000000e-01
# [2,]  0.5  1.000000e+00 -2.267044e-15
# [3,]  0.8 -2.267044e-15  5.000000e-01

Update based on comments

If the problem is to find a positive definite correlation matrix, you can use Matrix::nearPD to find the nearest positive definite matrix:

set.seed(1)
sig <- structure(c(0.000225,  0.0002625, 0.00044625, 
                   0.0002625 , 0.000625, 0, 
                   0.00044625, 0       , 0.001225), 
                 .Dim = c(3L, 3L))
cov(mvrnorm(1e5, rep(0, 3), Matrix::nearPD(sig, TRUE, TREU)$mat, empirical = TRUE))

#            V1           V2           V3
# V1 1.00000000 2.625000e-04 4.462500e-04
# V2 0.00026250 1.000000e+00 3.614917e-15
# V3 0.00044625 3.614917e-15 1.000000e+00

Upvotes: 1

Related Questions