Class for the multivariate Normal distribution. More...

#include <STK_MultiLaw_Normal.h>

Inheritance diagram for STK::MultiLaw::Normal< RowVector >:

Public Types
typedef MultiLaw::IMultiLaw< RowVector >	Base

Public Member Functions
	Normal (RowVector const &mu, ArraySquareX const &sigma)
	Constructor.

virtual	~Normal ()
	destructor.

RowVector const &	mu () const
	@return the location parameter

ArraySquareX const &	sigma () const
	@return the variance-covariance matrix

ArraySquareX const &	squareroot () const
	@return the square root of the variance-covariance matrix

SymEigen< ArraySquareX > const &	decomp () const
	@return the eigenvalue decomposition

void	setParameters (RowVector const &mu, ArraySquareX const &sigma)
	update the parameters specific to the law.

virtual Real	pdf (RowVector const &x) const
	compute the probability distribution function (density) of the multivariate normal law

Real	lpdf (RowVector const &x) const
	compute the log probability distribution function.

template<class Array >
Real	lnLikelihood (Array const &data) const
	compute the log likehood of a data set.

virtual void	rand (RowVector &x) const
	simulate a realization of the Multivariate Law and store the result in x.

template<class Array >
void	rand (ArrayBase< Array > &x) const
	simulate a realization of the Multivariate Law and store the result in x (using a reference vector).

Public Member Functions inherited from STK::MultiLaw::IMultiLaw< RowVector >
virtual	~IMultiLaw ()
	destructor.

Public Member Functions inherited from STK::Law::ILawBase
String const &	name () const

Protected Attributes
RowVector	mu_
	The position parameter.

ArraySquareX	sigma_
	The covariance parameter.

SymEigen< ArraySquareX >	decomp_
	the decomposition in eigenvalues of the covariance matrix

ArrayDiagonalX	invEigenvalues_
	inverse of the eigenvalues of sigma_

ArraySquareX	squareroot_
	The square root of the matrix `Sigma_`.

Protected Attributes inherited from STK::Law::ILawBase
String	name_
	Name of the Law.

Additional Inherited Members
Protected Member Functions inherited from STK::MultiLaw::IMultiLaw< RowVector >
	IMultiLaw (String const &name)
	Constructor.

Protected Member Functions inherited from STK::Law::ILawBase
	ILawBase (String const &name)
	Constructor.

	~ILawBase ()
	destructor.

Detailed Description

template<class RowVector>
class STK::MultiLaw::Normal< RowVector >

Class for the multivariate Normal distribution.

In probability theory and statistics, the "multivariate normal distribution" or "multivariate Gaussian distribution", is a generalization of the one-dimensional (univariate) normal distribution to higher dimensions. A random vector is said to be multivariate normally distributed if every linear combination of its components has a univariate normal distribution. The multivariate normal distribution is often used to describe, at least approximately, any set of (possibly) correlated real-valued random variables each of which clusters around a mean value.

The multivariate normal distribution of a $ p$ -dimensional random vector

$\mathbf{X} = \left( X_1, X_2, \ldots, X_p \right)'$

can be written in the following notation

$\mathbf{X}\ \sim\ \mathcal{N}(\mu,\ \Sigma).$

with $ p $ -dimensional mean vector

$\mu = \left( \mathrm{E}[X_1], \mathrm{E}[X_2], \ldots, \mathrm{E}[X_k] \right)'$

and $p \times p$ covariance matrix

$\Sigma = [\mathrm{Cov}(X_i, X_j)]_{i=1,2,\ldots,p;\ j=1,2,\ldots,p}$

Definition at line 100 of file STK_MultiLaw_Normal.h.

Member Typedef Documentation

◆ Base

template<class RowVector >

typedef MultiLaw::IMultiLaw<RowVector> STK::MultiLaw::Normal< RowVector >::Base

Definition at line 103 of file STK_MultiLaw_Normal.h.

Constructor & Destructor Documentation

◆ Normal()

template<class RowVector >

STK::MultiLaw::Normal< RowVector >::Normal	(	RowVector const &	mu,
		ArraySquareX const &	sigma
	)

inline

Constructor.

Parameters

mu	mean of the Normal distribution
sigma	covariance matrix of the Normal distribution

Definition at line 108 of file STK_MultiLaw_Normal.h.

         : Base(_T("MultiLaw::Normal"))
               , mu_()
               , sigma_()
               , decomp_(sigma)
    { setParameters(mu, sigma);}

References STK::MultiLaw::Normal< RowVector >::mu(), STK::MultiLaw::Normal< RowVector >::setParameters(), and STK::MultiLaw::Normal< RowVector >::sigma().

◆ ~Normal()

template<class RowVector >

virtual STK::MultiLaw::Normal< RowVector >::~Normal ( )

inlinevirtual

destructor.

Definition at line 116 of file STK_MultiLaw_Normal.h.

116{}

Member Function Documentation

◆ decomp()

template<class RowVector >

SymEigen< ArraySquareX > const & STK::MultiLaw::Normal< RowVector >::decomp ( ) const

inline

@return the eigenvalue decomposition

Definition at line 124 of file STK_MultiLaw_Normal.h.

124{ return decomp_;}

References STK::MultiLaw::Normal< RowVector >::decomp_.

◆ lnLikelihood()

template<class RowVector >

template<class Array >

Real STK::MultiLaw::Normal< RowVector >::lnLikelihood ( Array const & data ) const

inline

compute the log likehood of a data set.

sum the values of the log-pdf at the points stored in x.

Parameters

data	the multivariate values to compute the lpdf.

Returns: the value of the log-pdf

Definition at line 197 of file STK_MultiLaw_Normal.h.

    {
      // check ranges
      if (data.cols() != mu_.range() )
      { STKRUNTIME_ERROR_NO_ARG(MultiLaw::Normal::lnLikelihood(x),data.cols() != mu_.range());}
      // get dimensions of the samples and sum over all ln-likelihood values
      const int first = data.beginRows(), last = data.lastIdxRows();
      Real sum = 0.0;
      for (int i=first; i<= last; i++)
      {
        // compute lpdf using \sum_i ||(x_i-mu)'P||_{D^{-1}}^2
        Real res = 0.5 * ((data.row(i) - mu_) * decomp_.rotation()).wnorm2(invEigenvalues_) ;
        sum += res;
      }
      sum += data.sizeRows()*( invEigenvalues_.size() * Const::_LNSQRT2PI_
                             + 0.5 * log((double)decomp_.det())
                             );
      return -sum;
    }

References STK::MultiLaw::Normal< RowVector >::decomp_, STK::MultiLaw::Normal< RowVector >::invEigenvalues_, STK::MultiLaw::Normal< RowVector >::lnLikelihood(), STK::MultiLaw::Normal< RowVector >::mu_, STKRUNTIME_ERROR_NO_ARG, and STK::sum().

Referenced by STK::GaussianAAModel< Array >::computeProjectedLnLikelihood(), and STK::MultiLaw::Normal< RowVector >::lnLikelihood().

◆ lpdf()

template<class RowVector >

Real STK::MultiLaw::Normal< RowVector >::lpdf ( RowVector const & x ) const

inlinevirtual

compute the log probability distribution function.

Give the value of the log-pdf at the point x.

Parameters

x	the multivariate value to compute the lpdf.

Returns: the value of the log-pdf

Implements STK::MultiLaw::IMultiLaw< RowVector >.

Definition at line 179 of file STK_MultiLaw_Normal.h.

    {
      // check ranges
      if (x.range() != mu_.range() )
      { STKRUNTIME_ERROR_NO_ARG(MultiLaw::Normal::lpdf(x),x.range() != mu_.range());}
      // compute pdf using ||(x-mu)'P||_{D^{-1}}^2
      Real res = 0.5 * ((x - mu_) * decomp_.rotation()).wnorm2(invEigenvalues_)
               + invEigenvalues_.size() * Const::_LNSQRT2PI_
               + 0.5 * log((double)decomp_.det());
      return -res;
    }

References STK::MultiLaw::Normal< RowVector >::decomp_, STK::MultiLaw::Normal< RowVector >::invEigenvalues_, STK::MultiLaw::Normal< RowVector >::lpdf(), STK::MultiLaw::Normal< RowVector >::mu_, and STKRUNTIME_ERROR_NO_ARG.

Referenced by STK::MultiLaw::Normal< RowVector >::lpdf(), and STK::MultiLaw::Normal< RowVector >::pdf().

◆ mu()

template<class RowVector >

RowVector const & STK::MultiLaw::Normal< RowVector >::mu ( ) const

inline

@return the location parameter

Definition at line 118 of file STK_MultiLaw_Normal.h.

118{ return mu_;}

References STK::MultiLaw::Normal< RowVector >::mu_.

Referenced by STK::MultiLaw::Normal< RowVector >::Normal(), and STK::MultiLaw::Normal< RowVector >::setParameters().

◆ pdf()

template<class RowVector >

virtual Real STK::MultiLaw::Normal< RowVector >::pdf ( RowVector const & x ) const

inlinevirtual

compute the probability distribution function (density) of the multivariate normal law

$f(x) = \frac{1}{ (2\pi)^{k/2}|\Sigma|^{1/2} } \exp\!\left( {-\tfrac{1}{2}}(x-\mu)'\Sigma^{-1}(x-\mu) \right),$

Give the value of the pdf at the point x.

Parameters

x	the multivariate value to compute the pdf.

Returns: the value of the pdf

Implements STK::MultiLaw::IMultiLaw< RowVector >.

Definition at line 162 of file STK_MultiLaw_Normal.h.

    {
      // check determinant is not 0
      if (decomp_.det() == 0.)
      { STKRUNTIME_ERROR_NO_ARG(MultiLaw::Normal::pdf(x),|sigma| == 0.);}
      // check ranges
      if (x.range() != mu_.range() )
      { STKRUNTIME_ERROR_NO_ARG(MultiLaw::Normal::pdf(x),x.range() != mu_.range());}
      // compute pdf in log-space
      return std::exp((double)lpdf(x));
    }

References STK::MultiLaw::Normal< RowVector >::decomp_, STK::MultiLaw::Normal< RowVector >::lpdf(), STK::MultiLaw::Normal< RowVector >::mu_, STK::MultiLaw::Normal< RowVector >::pdf(), STK::MultiLaw::Normal< RowVector >::sigma(), and STKRUNTIME_ERROR_NO_ARG.

Referenced by STK::MultiLaw::Normal< RowVector >::pdf().

◆ rand() [1/2]

template<class RowVector >

template<class Array >

void STK::MultiLaw::Normal< RowVector >::rand ( ArrayBase< Array > & x ) const

inline

simulate a realization of the Multivariate Law and store the result in x (using a reference vector).

The class RowVector have to derive from IContainerRef.

Parameters

[out] x the simulated value.

Definition at line 234 of file STK_MultiLaw_Normal.h.

    {
      // fill it with iid N(0,1) variates
      x.randGauss();
      // rotate with squareroot_ and translate with mu_
      x.asDerived() = x.asDerived() * squareroot_ + Const::Vector<Real>(x.rows()) * mu_;
    }

References STK::MultiLaw::Normal< RowVector >::mu_, and STK::MultiLaw::Normal< RowVector >::squareroot_.

◆ rand() [2/2]

template<class RowVector >

virtual void STK::MultiLaw::Normal< RowVector >::rand ( RowVector & x ) const

inlinevirtual

simulate a realization of the Multivariate Law and store the result in x.

Parameters

[out] x the simulated value.

Implements STK::MultiLaw::IMultiLaw< RowVector >.

Definition at line 221 of file STK_MultiLaw_Normal.h.

    {
      // fill it with iid N(0,1) variates
      x.randGauss();
      // rotate with squareroot_ and translate with mu_
      x = x * squareroot_ + mu_;
    }

References STK::MultiLaw::Normal< RowVector >::mu_, and STK::MultiLaw::Normal< RowVector >::squareroot_.

◆ setParameters()

template<class RowVector >

void STK::MultiLaw::Normal< RowVector >::setParameters	(	RowVector const &	mu,
		ArraySquareX const &	sigma
	)

inline

update the parameters specific to the law.

Definition at line 126 of file STK_MultiLaw_Normal.h.

    {
      // check dimensions
      if (mu.range() != sigma.range())
      { STKRUNTIME_ERROR_NO_ARG(MultiLaw::Normal::setParameters(mu, sigma),mu.range() != sigma.range());}
      mu_    = mu;
      sigma_ = sigma;
      // decomposition of the covariance matrix
      decomp_.setData(sigma_);
      if (!decomp_.run())
      {
 
        throw runtime_error(decomp_.error());
      }
      // compute the inverse of the eigenvalues of sigma_ and the squareroot_
      // matrix needed by rand
      invEigenvalues_.resize(mu_.range());
      squareroot_.resize(mu_.range());
      // get dimension
      int rank = decomp_.rank(), end = mu_.begin() + rank;
      for (int j=mu_.begin(); j< end; j++)
      { invEigenvalues_[j] = 1./decomp_.eigenValues()[j];}
      for (int j=end; j< mu_.end(); j++) { invEigenvalues_[j] = 0.;}
 
      squareroot_ = decomp_.rotation() * decomp_.eigenValues().sqrt().diagonalize() * decomp_.rotation().transpose();
    }