rmnorm: Random number generator for (conditional) multivariate normal distribution

Description

This function generates random numbers (i.e. variates) from (conditional) multivariate normal distribution.

Usage

rmnorm(
  n,
  mean,
  sigma,
  given_ind = numeric(),
  given_x = numeric(),
  dependent_ind = numeric(),
  is_validation = TRUE,
  n_cores = 1L
)

Value

This function returns a numeric matrix which rows a random variates from (conditional) multivariate normal distribution with mean equal to mean and covariance equal to sigma. If given_x and given_ind are also provided then random variates will be from conditional multivariate normal distribution. Please, see details section of cmnorm to get additional information on the conditioning procedure.

Arguments

n: positive integer representing the number of random variates to be generated from (conditional) multivariate normal distribution. If given_ind is not empty vector then n should be be equal to nrow(given_x).
mean: numeric vector representing expectation of multivariate normal vector (distribution).
sigma: positively defined numeric matrix representing covariance matrix of multivariate normal vector (distribution).
given_ind: numeric vector representing indexes of multivariate normal vector which are conditioned at values given by given_x argument.
given_x: numeric vector which i-th element corresponds to the given value of the given_ind[i]-th element (component) of multivariate normal vector. If given_x is numeric matrix then it's rows are such vectors of given values.
dependent_ind: numeric vector representing indexes of unconditional elements (components) of multivariate normal vector.
is_validation: logical value indicating whether input arguments should be validated. Set it to FALSE to get performance boost (default value is TRUE).
n_cores: positive integer representing the number of CPU cores used for parallel computing. Currently it is not recommended to set n_cores > 1 if vectorized arguments include less then 100000 elements.

Details

This function uses Cholesky decomposition to generate multivariate normal variates from independent standard normal variates.

Examples

Run this code

# Consider multivariate normal vector:
# X = (X1, X2, X3, X4, X5) ~ N(mean, sigma)

# Prepare multivariate normal vector parameters
  # expected value
mean <- c(-2, -1, 0, 1, 2)
n_dim <- length(mean)
  # correlation matrix
cor <- c(   1,  0.1,  0.2,   0.3,  0.4,
          0.1,    1, -0.1,  -0.2, -0.3,
          0.2, -0.1,    1,   0.3,  0.2,
          0.3, -0.2,  0.3,     1, -0.05,
          0.4, -0.3,  0.2, -0.05,     1)
cor <- matrix(cor, ncol = n_dim, nrow = n_dim, byrow = TRUE)
  # covariance matrix
sd_mat <- diag(c(1, 1.5, 2, 2.5, 3))
sigma <- sd_mat %*% cor %*% t(sd_mat)

# Simulate random variates from this distribution
rmnorm(n = 3, mean = mean, sigma = sigma)

# Simulate random variate from (X1, X3, X5 | X1 = -1, X4 = 1)
given_x <- c(-1, 1)
given_ind = c(1, 4)
rmnorm(n = 1, mean = mean, sigma = sigma, 
       given_x = given_x, given_ind = given_ind)

# Simulate random variate from (X1, X3, X5 | X1 = -1, X4 = 1)
# and (X1, X3, X5 | X1 = 2, X4 = 3)
given_x = rbind(c(-1, 1), c(2, 3))
rmnorm(n = nrow(given_x), mean = mean, sigma = sigma, 
       given_x = given_x, given_ind = given_ind)