мврнорм (от масс) против рмвнорм (от мвтнорм)

Question

Я генерирую большой объем данных из многомерного нормального распределения для моделирования. Интересно, кто-нибудь знает, какая команда для этого наиболее эффективна? Если это mvrnorm (из пакета «MASS») или rmvnorm (из пакета «mvtnorm»).

Answer 1

На такие вопросы легко ответить, выбрав разные подходы. Пусть

library(microbenchmark)
library(MASS)
library(mvtnorm)

n <- 10000
k <- 50
mu <- rep(0, k)
rho <- 0.2
Sigma <- diag(k) * (1 - rho) + rho 

Таким образом, у нас есть 50 переменных с единичной дисперсией и корреляцией 0,2. Генерируя 10000 наблюдений, мы получаем

microbenchmark(mvrnorm(n, mu = mu, Sigma = Sigma),
               rmvnorm(n, mean = mu, sigma = Sigma, method = "eigen"),
               rmvnorm(n, mean = mu, sigma = Sigma, method = "svd"),
               rmvnorm(n, mean = mu, sigma = Sigma, method = "chol"),
               times = 100)
# Unit: milliseconds
#                                                    expr      min       lq     mean   median        uq      max neval cld
#                      mvrnorm(n, mu = mu, Sigma = Sigma) 65.04667 73.02912 85.30384 81.70611  92.69137 148.6959   100  a 
#  rmvnorm(n, mean = mu, sigma = Sigma, method = "eigen") 71.14170 81.08311 95.12891 88.84669 100.62174 237.0012   100   b
#    rmvnorm(n, mean = mu, sigma = Sigma, method = "svd") 71.32999 81.30640 93.40939 88.54804 104.00281 208.3690   100   b
#   rmvnorm(n, mean = mu, sigma = Sigma, method = "chol") 71.22712 78.59898 94.13958 89.04653 108.27363 158.7890   100   b

Таким образом, возможно, mvrnorm работает немного лучше. Поскольку у вас есть конкретное приложение, вы должны установить значения n, k и Sigma, более подходящие для этого приложения.

Поскольку вы, похоже, не ограничены этими двумя подходами, вы можете изучить Rcpp альтернативы; см., например, 1 , 2 , 3 .