R : 두 변량 일반 샘플을 생성 할 때 'mu'또는 'Sigma'를 MASS :: mvrnorm()에서 벡터화 할 수 있습니까?

Question

는

n <- 3 
phi0 <- 1 
phi1 <- 0.1 
phi2 <- 0.2 
W2 <- runif(n,0,1) 
W3 <- runif(n,0,1) 
W <- cbind(W2,W3) 
phi <- rbind(phi0,phi1,phi2) 
rho <- 0.4 
sigma1 <- exp(as.numeric(model.matrix(~W) %*% phi)) 
sigma2<- 1 

library(MASS) 
#Sigma <- ??? 
mu <- rep(0,2) 
v <- mvrnorm(n, mu, Sigma) 

sigma1

이 분산 내 벡터, 아래의 데이터를 고려한다.

정상 이변 분포를 갖는 n 길이의 이진수 벡터 v=(v1,v2)을 생성하고 싶습니다. v의 i 번째 줄은 평균이 mu=(0,0)이고, 상관은 rho=0.4이고 한계 분산은 sigma=(sigma1, 1)입니다. 여기서 sigma1은 i 번째 줄 sigma1의 값을받습니다. 어떻게하면됩니까?

---

1D 케이스

편집자 설명은 rnorm는 musigma 및 벡터화를 수락

rnorm(3, 0, sqrt(sigma1)) 

가 N(0, sqrt(sigma1[i]))로부터 샘플들을 제공되도록. 기본적으로 OP는 mvrnorm에 대해 동일한 기능을 요구합니다.

Answer 1

염기성 솔루션

아니요, 벡터화 할 수 없습니다. for 루프를 작성하십시오.

v <- matrix(0, n, 2) 
for (i in 1:n) { 
    sig11 <- sigma1[i] 
    sig21 <- rho * sqrt(sig11) 
    Sigma <- matrix(c(sig11, sig21, sig21, 1), 2) 
    v[i, ] <- mvrnorm(1, c(0,0), Sigma) 
    } 

공분산에게 주어 고급 솔루션

Sigma

sig11^2    rho * sig11 * sig22 
rho * sig11 * sig22  sig22^2 

의 낮은 삼각 콜레 요인 L 다음

sig11   0 
rho * sig22 sqrt(1 - rho^2) * sig22 

하면 x <- rnorm(2), mu + L %*% x 것은 0,123,365의 샘플입니다.

이것은 당신의 데이터의 완전 벡터화 솔루션

# mu1, mu2, sig11, sig22, rho can be 
## length-n vectors 
## scalars 
## vectors than can be recycled to be length-n 
birnorm <- function (n, mu1, mu2, sig11, sig22, rho) { 
    x1 <- rnorm(n) 
    x2 <- rnorm(n) 
    z1 <- sig11 * x1 + mu1 
    z2 <- rho * sig22 * x1 + sqrt(1 - rho^2) * sig22 * x2 + mu2 
    cbind(z1, z2) 
    } 

을 제공합니다, 당신은

v <- birnorm(3, 0, 0, sqrt(sigma1), 1, 0.4) 

---

가이 기능을 사용할 수 있다는 이변 경우에

참고 mvrnorm을 다시 사용할 수 있습니다 bivariate case의 경우 mvrnorm을 다시 작성하십시오.

mvrnorm2 <- function (n, mu, Sigma) { 
    sig11 <- sqrt(Sigma[1]) 
    sig22 <- sqrt(Sigma[4]) 
    rho <- Sigma[2]/(sig11 * sig22) 
    birnorm(n, mu[1], mu[2], sig11, sig22, rho) 
    } 

그것은

mu <- c(0,0) 
Sigma <- matrix(c(1,0.5,0.5,1),2) 
library(microbenchmark) 
microbenchmark(mvrnorm(1000, mu, Sigma), mvrnorm2(1000, mu, Sigma))