Class to generate normal distributed, correlated random numbers given the mean values and the covariance matrix of all dimensions. More...

#include <MultivariateNormalGenerator.h>

Public Member Functions
	MultivariateNormalGenerator ()
	default constructor to allow later initialization

	MultivariateNormalGenerator (int n, const double mean, const double cov)
	constructor with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

	MultivariateNormalGenerator (const Eigen::VectorXd &mean, const Eigen::MatrixXd &cov)
	constructor with Eigen matrix interface.

Eigen::VectorXd	generate () const
	Generate a set of correlated random numbers with the previouly set mean and covariance.

void	reset ()
	reset the generator setting the size to 0.

size_t	size () const
	Return the number of elements to be generated on generate()

void	generate (double *output) const
	Generate a set of correlated random numbers with the previouly set mean and covariance and store them in buffer output.

ROOT::Math::XYZVector	generateVec3 () const
	Generate a set of correlated random numbers with the previously set mean and covariance and return a ROOT::Math::XYZVector.

TVectorD	generateVecT () const
	Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVectorT<double>

bool	setMeanCov (int n, const double mean, const double cov)
	set the mean and covariance for the distribution with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

bool	setMeanCov (const Eigen::VectorXd &mean, const Eigen::MatrixXd &cov)
	set the mean and covariance for the distribution.

template<class value_type >
bool	setMeanCov (const TVectorT< value_type > &mean, const TMatrixTBase< value_type > &cov)
	set mean and covariance matrix from ROOT vector/matrix objects, e.g.

template<class value_type >
bool	setMeanCov (const ROOT::Math::XYZVector &mean, const TMatrixTBase< value_type > &cov)
	set the mean and covariance for the distribution.

Private Attributes
Eigen::VectorXd	m_mean
	Member to store the mean values of the distribution.

Eigen::MatrixXd	m_transform
	Member to store the transformation matrix for standard normal distributed random values.

Detailed Description

Class to generate normal distributed, correlated random numbers given the mean values and the covariance matrix of all dimensions.

This class can be used to generate normal distributed random values according to a given covariance matrix (assuming the covariance matrix is positive semi-definite).

To use it first set the desired mean values and covariance matrix using setMeanCov() and then call generate() to generate one set of values.

Warning: : setMeanCov() will not work for all matrices. Please check the return value when setting the covariance matrix.

To generate normal distributed random values according to a covariance matrix we need to decompose the covariance matrix $M$ into $M= A A^T$ . Given the vector of mean values as $\mu$ and a vector of standard normal distributed random values $(\mu=0, \sigma=1)$ as n we can obtain a set of correlated random values $x = \mu + A * n$ .

Usually the Cholesky decomposition is chosen as it computes $M = L L^T$ . However the Cholesky composition only works for positive definite matrices so it does not work if for example one of the values is fixed and has no error.

To ease this restriction a little we use the LDLT decomposition given as $M = P^T L D L^T P$ where $P$ is a permutation matrix and D a diagonal matrix. We then can use $A = P^T L \sqrt(D)$ to caluclate the correlated values also for positive semi-definite covariance matrices if the elements of D are positive.

Definition at line 54 of file MultivariateNormalGenerator.h.

Constructor & Destructor Documentation

◆ MultivariateNormalGenerator() [1/3]

MultivariateNormalGenerator ( )

inline

default constructor to allow later initialization

Definition at line 57 of file MultivariateNormalGenerator.h.

57{}

◆ MultivariateNormalGenerator() [2/3]

MultivariateNormalGenerator	(	int	n,
		const double *	mean,
		const double *	cov
	)

inline

constructor with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

Parameters

n	dimensionality
mean	pointer to the n mean values of the distribution
cov	pointer to the n*n covariance values in row major layout

Definition at line 65 of file MultivariateNormalGenerator.h.

    {
      setMeanCov(n, mean, cov);
    }

◆ MultivariateNormalGenerator() [3/3]

MultivariateNormalGenerator	(	const Eigen::VectorXd &	mean,
		const Eigen::MatrixXd &	cov
	)

inline

constructor with Eigen matrix interface.

Parameters

mean	Vector of mean values
cov	Matrix containing the covariance values

Definition at line 73 of file MultivariateNormalGenerator.h.

    {
      setMeanCov(mean, cov);
    }

Member Function Documentation

◆ generate() [1/2]

Eigen::VectorXd generate ( ) const

inline

Generate a set of correlated random numbers with the previouly set mean and covariance.

Returns: Vector containing the generated random numbers

Definition at line 81 of file MultivariateNormalGenerator.h.

    {
      //To get the correlated multivariate normal distribution, we multiply
      //standard normal distributed values (mean=0, sigma=1) with a
      //transformation matrix we obtained from an LDLT decomposition and add
      //the mean values.
      Eigen::VectorXd x(m_mean.rows());
      for (int i = 0; i < m_mean.rows(); ++i) x(i) = gRandom->Gaus();
      return m_mean + (m_transform * x);
    }

◆ generate() [2/2]

void generate ( double * output ) const

inline

Generate a set of correlated random numbers with the previouly set mean and covariance and store them in buffer output.

Parameters

output pointer to array where generated values will be stored.

Definition at line 101 of file MultivariateNormalGenerator.h.

    {
      Eigen::VectorXd x = generate();
      for (int i = 0; i < x.rows(); ++i) { output[i] = x(i); }
    }

◆ generateVec3()

ROOT::Math::XYZVector generateVec3 ( ) const

inline

Generate a set of correlated random numbers with the previously set mean and covariance and return a ROOT::Math::XYZVector.

Optimally, the set mean and covariance matrix should be of dimension three, otherwise just the first size() elements of the ROOT::Math::XYZVector are set and the remaining elements are zero. If size() is bigger than 3 the remaining values will be discarded.

Definition at line 113 of file MultivariateNormalGenerator.h.

    {
      Eigen::VectorXd x = generate();
      std::array<double, 3> tmp = {0.};
      for (unsigned int i = 0; i < std::min(3u, (unsigned int)size()); ++i) {
        tmp[i] = x(i);
      }
      return ROOT::Math::XYZVector(tmp[0], tmp[1], tmp[2]);
    }

◆ generateVecT()

TVectorD generateVecT ( ) const

inline

Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVectorT<double>

Definition at line 126 of file MultivariateNormalGenerator.h.

    {
      Eigen::VectorXd x = generate();
      TVectorD output(x.rows());
      output.SetElements(x.data());
      return output;
    }

◆ reset()

void reset ( )

reset the generator setting the size to 0.

Subsequent calls to generate will return 0-sized results until the generator is reinitialized using setMeanCov()

Definition at line 14 of file MultivariateNormalGenerator.cc.

{
  m_mean.resize(0);
  m_transform.resize(0, 0);
}

◆ setMeanCov() [1/2]

bool setMeanCov	(	const Eigen::VectorXd &	mean,
		const Eigen::MatrixXd &	cov
	)

set the mean and covariance for the distribution.

Parameters

mean	Vector of mean values
cov	Matrix containing the covariance values

Returns: true if covariance could be decomposited, false otherwise

Definition at line 36 of file MultivariateNormalGenerator.cc.

{
  reset();
  if (mean.rows() != cov.rows()) {
    B2ERROR("Mean values and covariance matrix need to be of the same dimension");
    return false;
  }
  if (cov.rows() != cov.cols()) {
    B2ERROR("Covariance matrix needs to be a square matrix");
    return false;
  }
 
  // The usual way to calulate multivariate normal distributed values is to use
  // the cholesky decomposition C = ll' and then calculate y = mean + l * x
  // where x is a vector of standard normal distributed values. However this
  // only works for positive definite matrices so we use the LDLT composition
  // which gives us C = P'LDL'P and we compute L=P'Lsqrt(D) to get the same
  // result but in a more robust way which works also for semi-definite
  // matrices
  auto ldlt = cov.ldlt();
  if (ldlt.info() != Eigen::Success) {
    B2ERROR("Cannot compute LDLT decomposition of covariance "
            "matrix, maybe not positive semi-definite?");
    return false;
  }
  Eigen::MatrixXd L = ldlt.matrixL();
  Eigen::MatrixXd D = ldlt.vectorD().asDiagonal();
  if (D.minCoeff() < 0) {
    B2ERROR("MultivariateNormalGenerator: Negative values when computing LDL^T "
            "decomposition, cannot compute M=AA^T, resulting random numbers "
            "will not be correct");
    return false;
  }
  auto P = ldlt.transpositionsP().transpose();
  m_transform = P * L * D.cwiseSqrt();
  //And save the mean values
  m_mean = mean;
  return true;
}

◆ setMeanCov() [2/2]

bool setMeanCov	(	int	n,
		const double *	mean,
		const double *	cov
	)

set the mean and covariance for the distribution with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

Parameters

n	dimensionality
mean	pointer to the n mean values of the distribution
cov	pointer to the n*n covariance values in row major layout

Returns: true if covariance could be decomposited, false otherwise

Definition at line 20 of file MultivariateNormalGenerator.cc.

{
  Eigen::VectorXd emean(n);
  Eigen::MatrixXd ecov(n, n);
  //Fill the vector and matrix from the buffers in row major layout
  for (int i = 0; i < n; ++i) {
    emean[i] = mean[i];
    for (int j = 0; j < n; ++j) {
      ecov(i, j) = cov[i * n + j];
    }
  }
  //And delegate to the correct function
  return setMeanCov(emean, ecov);
}

◆ size()

size_t size ( ) const

inline

Return the number of elements to be generated on generate()

Definition at line 96 of file MultivariateNormalGenerator.h.

96{ return m_mean.rows(); }

Member Data Documentation

◆ m_mean

Eigen::VectorXd m_mean

private

Member to store the mean values of the distribution.

Definition at line 171 of file MultivariateNormalGenerator.h.

◆ m_transform

Eigen::MatrixXd m_transform

private

Member to store the transformation matrix for standard normal distributed random values.

Definition at line 174 of file MultivariateNormalGenerator.h.

The documentation for this class was generated from the following files:

framework/utilities/include/MultivariateNormalGenerator.h
framework/utilities/src/MultivariateNormalGenerator.cc

Public Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ MultivariateNormalGenerator() [1/3]

◆ MultivariateNormalGenerator() [2/3]

◆ MultivariateNormalGenerator() [3/3]

Member Function Documentation

◆ generate() [1/2]

◆ generate() [2/2]

◆ generateVec3()

◆ generateVecT()

◆ reset()

◆ setMeanCov() [1/2]

◆ setMeanCov() [2/2]

◆ size()

Member Data Documentation

◆ m_mean

◆ m_transform