nag_ref_vec_multi_normal (g05eac) (PDF version)
g05 Chapter Contents
g05 Chapter Introduction
NAG C Library Manual

NAG Library Function Document

nag_ref_vec_multi_normal (g05eac)

+ Contents

    1  Purpose
    7  Accuracy

1  Purpose

nag_ref_vec_multi_normal (g05eac) sets up a reference vector for a multivariate Normal distribution with mean vector a  and variance-covariance matrix C , so that nag_return_multi_normal (g05ezc) may be used to generate pseudorandom vectors.

2  Specification

#include <nag.h>
#include <nagg05.h>
void  nag_ref_vec_multi_normal (const double a[], Integer n, const double c[], Integer tdc, double eps, double **r, NagError *fail)

3  Description

When the variance-covariance matrix is non-singular (i.e., strictly positive definite), the distribution has probability density function
f x = C -1 2π n exp - x-aT C -1 x-a
where n  is the number of dimensions, C  is the variance-covariance matrix, a  is the vector of means and x  is the vector of positions.
Variance-covariance matrices are symmetric and positive semidefinite. Given such a matrix C , there exists a lower triangular matrix L  such that LLT = C . L  is not unique, if C  is singular.
nag_ref_vec_multi_normal (g05eac) decomposes C  to find such an L . It then stores n , a  and L  in the reference vector r  for later use by nag_return_multi_normal (g05ezc). nag_return_multi_normal (g05ezc) generates a vector x  of independent standard Normal pseudorandom numbers. It then returns the vector a + Lx , which has the required multivariate Normal distribution.
It should be noted that this function will work with a singular variance-covariance matrix C , provided C  is positive semidefinite, despite the fact that the above formula for the probability density function is not valid in that case. Wilkinson (1965) should be consulted if further information is required.

4  References

Knuth D E (1981) The Art of Computer Programming (Volume 2) (2nd Edition) Addison–Wesley
Wilkinson J H (1965) The Algebraic Eigenvalue Problem Oxford University Press, Oxford

5  Arguments

1:     a[n]const doubleInput
On entry: the vector of means, a , of the distribution.
2:     nIntegerInput
On entry: the number of dimensions, n , of the distribution.
Constraint: n>0 .
3:     c[n×tdc]const doubleInput
Note: the i,jth element of the matrix C is stored in c[i-1×tdc+j-1].
On entry: the variance-covariance matrix of the distribution. Only the upper triangle need be set.
4:     tdcIntegerInput
On entry: the stride separating matrix column elements in the array c.
Constraint: tdcn .
5:     epsdoubleInput
On entry: the maximum error in any element of C , relative to the largest element of C .
Constraint: 0.0 eps 0.1 / n .
6:     rdouble **Output
On exit: reference vector for which memory will be allocated internally. This reference vector will subsequently be used by nag_return_multi_normal (g05ezc). If no memory is allocated to r (e.g., when an input error is detected) then r will be NULL on return, otherwise you should use the NAG macro NAG_FREE to free the storage allocated by r when it is no longer of use.
7:     failNagError *Input/Output
The NAG error argument (see Section 3.6 in the Essential Introduction).

6  Error Indicators and Warnings

On entry, tdc=value  while n=value . These arguments must satisfy tdcn .
On entry, eps=value  while 0.1/ n=value . These arguments must satisfy eps0.1 /n.
Dynamic memory allocation failed.
On entry, n=value.
Constraint: n1.
Matrix C is not positive semidefinite.
On entry, eps must not be less than 0.0: eps=value .

7  Accuracy

The maximum absolute error in LLT , and hence in the variance-covariance matrix of the resulting vectors, is less than n × maxeps,ε + n+3 ε / 2  times the maximum element of C , where ε  is the machine precision. Under normal circumstances, the above will be small compared to sampling error.

8  Further Comments

The time taken by nag_ref_vec_multi_normal (g05eac) is of order n 3 .
It is recommended that the diagonal elements of C  should not differ too widely in order of magnitude. This may be achieved by scaling the variables if necessary. The actual matrix decomposed is C + E = LLT , where E  is a diagonal matrix with small positive diagonal elements. This ensures that, even when C  is singular, or nearly singular, the Cholesky factor L  corresponds to a positive definite variance-covariance matrix that agrees with C  within a tolerance determined by eps.

9  Example

The example program prints five pseudorandom observations from a bivariate Normal distribution with means vector
1.0 2.0
and variance-covariance matrix
2.0 1.0 1.0 3.0 ,
generated by nag_ref_vec_multi_normal (g05eac) and nag_return_multi_normal (g05ezc) after initialization by nag_random_init_repeatable (g05cbc).

9.1  Program Text

Program Text (g05eace.c)

9.2  Program Data


9.3  Program Results

Program Results (g05eace.r)

nag_ref_vec_multi_normal (g05eac) (PDF version)
g05 Chapter Contents
g05 Chapter Introduction
NAG C Library Manual

© The Numerical Algorithms Group Ltd, Oxford, UK. 2012