A class for handling Log-Normal joint PDFs. More...

#include <LogNormalJointPdf.h>

Inheritance diagram for QUESO::LogNormalJointPdf< V, M >:

Public Member Functions
Constructor/Destructor methods
	LogNormalJointPdf (const char *prefix, const VectorSet< V, M > &domainSet, const V &lawExpVector, const V &lawVarVector)
	Constructor. More...

	~LogNormalJointPdf ()
	Destructor. More...

Math methods
double	actualValue (const V &domainVector, const V domainDirection, V gradVector, M hessianMatrix, V hessianEffect) const
	Actual value of the Log-Normal PDF (scalar function). More...

double	lnValue (const V &domainVector, const V domainDirection, V gradVector, M hessianMatrix, V hessianEffect) const
	Logarithm of the value of the Log-Normal PDF (scalar function). More...

virtual void	distributionMean (V &meanVector) const
	Mean value of the underlying random variable. More...

virtual void	distributionVariance (M &covMatrix) const
	Covariance matrix of the underlying random variable. More...

double	computeLogOfNormalizationFactor (unsigned int numSamples, bool updateFactorInternally) const
	Computes the logarithm of the normalization factor. More...

const V &	lawExpVector () const
	Access to the vector of mean values and private attribute: m_lawExpVector. More...

const V &	lawVarVector () const
	Access to the vector of variance values and private attribute: m_lawVarVector. More...

Public Member Functions inherited from QUESO::BaseJointPdf< V, M >
	BaseJointPdf (const char *prefix, const VectorSet< V, M > &domainSet)
	Default constructor. More...

virtual	~BaseJointPdf ()
	Destructor. More...

virtual void	setNormalizationStyle (unsigned int value) const

void	setLogOfNormalizationFactor (double value) const
	Sets a logarithmic value to be used in the normalization factor (stored in the protected attribute m_normalizationStyle.) More...

virtual void	print (std::ostream &os) const
	Print method. Non-pure for backwards compatibility. More...

Public Member Functions inherited from QUESO::BaseScalarFunction< V, M >
void	setFiniteDifferenceStepSize (double fdStepSize)
	Sets the step size for finite differencing gradients. More...

void	setFiniteDifferenceStepSize (unsigned int i, double fdStepSize)

	BaseScalarFunction (const char *prefix, const VectorSet< V, M > &domainSet)
	Default constructor. More...

virtual	~BaseScalarFunction ()
	Destructor. More...

const VectorSet< V, M > &	domainSet () const
	Access to the protected attribute `m_domainSet:` domain set of the scalar function. More...

virtual double	lnValue (const V &domainVector) const
	Returns the logarithm of the function at `domainVector`. More...

virtual double	lnValue (const V &domainVector, V &gradVector) const
	Returns the logarithm of the function and its gradient at `domainVector`. More...

virtual double	lnValue (const V &domainVector, V &gradVector, const V &domainDirection, V &hessianEffect) const

Protected Attributes
V *	m_lawExpVector

V *	m_lawVarVector

bool	m_diagonalCovMatrix

Protected Attributes inherited from QUESO::BaseJointPdf< V, M >
unsigned int	m_normalizationStyle

double	m_logOfNormalizationFactor

Protected Attributes inherited from QUESO::BaseScalarFunction< V, M >
const BaseEnvironment &	m_env

std::string	m_prefix

const VectorSet< V, M > &	m_domainSet
	Domain set of the scalar function. More...

Additional Inherited Members
Protected Member Functions inherited from QUESO::BaseJointPdf< V, M >
double	commonComputeLogOfNormalizationFactor (unsigned int numSamples, bool updateFactorInternally) const
	Common method (to the derived classes) to compute the logarithm of the normalization factor. More...

Detailed Description

template<class V = GslVector, class M = GslMatrix>
class QUESO::LogNormalJointPdf< V, M >

A class for handling Log-Normal joint PDFs.

This class allows the mathematical definition of a Log-Normal Joint PDF.

Definition at line 50 of file LogNormalJointPdf.h.

Constructor & Destructor Documentation

template<class V , class M >

QUESO::LogNormalJointPdf< V, M >::LogNormalJointPdf	(	const char *	prefix,
		const VectorSet< V, M > &	domainSet,
		const V &	lawExpVector,
		const V &	lawVarVector
	)

Constructor.

Constructs a new object of the class, given a prefix and the domain of the PDF, a vector of mean values, lawExpVector, and a vector of covariance values lawVarVector (an alternative representation for a diagonal covariance matrix).

Definition at line 33 of file LogNormalJointPdf.C.

References QUESO::BaseEnvironment::displayVerbosity(), QUESO::LogNormalJointPdf< V, M >::lawExpVector(), QUESO::LogNormalJointPdf< V, M >::lawVarVector(), QUESO::BaseScalarFunction< V, M >::m_env, QUESO::BaseScalarFunction< V, M >::m_prefix, and QUESO::BaseEnvironment::subDisplayFile().

   :
   BaseJointPdf<V,M>(((std::string)(prefix)+"gau").c_str(),domainSet),
   m_lawExpVector     (new V(lawExpVector)),
   m_lawVarVector     (new V(lawVarVector)),
   m_diagonalCovMatrix(true)
 {
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 54)) {
     *m_env.subDisplayFile() << "Entering LogNormalJointPdf<V,M>::constructor() [1]"
                             << ": prefix = " << m_prefix
                             << std::endl;
   }
 
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 55)) {
     *m_env.subDisplayFile() << "In LogNormalJointPdf<V,M>::constructor()"
                           //<< ", prefix = "     << m_prefix
                             << ": meanVector = " << this->lawExpVector()
                       << ", Variances = "  << this->lawVarVector()
                             << std::endl;
   }
 
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 54)) {
     *m_env.subDisplayFile() << "Leaving LogNormalJointPdf<V,M>::constructor() [1]"
                             << ": prefix = " << m_prefix
                             << std::endl;
   }
 }

template<class V , class M >

QUESO::LogNormalJointPdf< V, M >::~LogNormalJointPdf ( )

Destructor.

Definition at line 66 of file LogNormalJointPdf.C.

 {
   delete m_lawVarVector;
   delete m_lawExpVector;
 }

Member Function Documentation

template<class V , class M >

double QUESO::LogNormalJointPdf< V, M >::actualValue	(	const V &	domainVector,
		const V *	domainDirection,
		V *	gradVector,
		M *	hessianMatrix,
		V *	hessianEffect
	)		const

virtual

Actual value of the Log-Normal PDF (scalar function).

This method calls lnValue() and applies the exponential to its result.

Implements QUESO::BaseJointPdf< V, M >.

Definition at line 88 of file LogNormalJointPdf.C.

References QUESO::queso_require_equal_to_msg.

 {
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 55)) {
     *m_env.subDisplayFile() << "Entering LogNormalJointPdf<V,M>::actualValue()"
                             << ", meanVector = "               << *m_lawExpVector
                             << ": domainVector = "             << domainVector
                             << ", domainVector.sizeLocal() = " << domainVector.sizeLocal()
                             << ", this->m_domainSet.vectorSpace().dimLocal() = " << this->m_domainSet.vectorSpace().dimLocal()
                             << std::endl;
   }
 
   queso_require_equal_to_msg(domainVector.sizeLocal(), this->m_domainSet.vectorSpace().dimLocal(), "invalid input");
 
   queso_require_msg(!(hessianMatrix || hessianEffect), "incomplete code for gradVector, hessianMatrix and hessianEffect calculations");
 
   double returnValue = 0.;
 
   V zeroVector(domainVector);
   zeroVector.cwSet(0.);
   if (domainVector.atLeastOneComponentSmallerOrEqualThan(zeroVector)) {
     // What should the gradient be here?
     returnValue = 0.;
   }
   else if (this->m_domainSet.contains(domainVector) == false) {
     // What should the gradient be here?
     returnValue = 0.;
   }
   else {
     // Already normalised
     returnValue = std::exp(this->lnValue(domainVector,domainDirection,gradVector,hessianMatrix,hessianEffect));
 
     if (gradVector) {
       (*gradVector) *= returnValue;
     }
   }
 
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 55)) {
     *m_env.subDisplayFile() << "Leaving LogNormalJointPdf<V,M>::actualValue()"
                             << ", meanVector = "   << *m_lawExpVector
                             << ": domainVector = " << domainVector
                             << ", returnValue = "  << returnValue
                             << std::endl;
   }
 
   return returnValue;
 }

template<class V , class M >

double QUESO::LogNormalJointPdf< V, M >::computeLogOfNormalizationFactor	(	unsigned int	numSamples,
		bool	updateFactorInternally
	)		const

virtual

Computes the logarithm of the normalization factor.

This routine calls BaseJointPdf::commonComputeLogOfNormalizationFactor().

Implements QUESO::BaseJointPdf< V, M >.

Definition at line 261 of file LogNormalJointPdf.C.

References QUESO::BaseJointPdf< V, M >::commonComputeLogOfNormalizationFactor().

 {
   double value = 0.;
 
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
     *m_env.subDisplayFile() << "Entering LogNormalJointPdf<V,M>::computeLogOfNormalizationFactor()"
                             << std::endl;
   }
   value = BaseJointPdf<V,M>::commonComputeLogOfNormalizationFactor(numSamples, updateFactorInternally);
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
     *m_env.subDisplayFile() << "Leaving LogNormalJointPdf<V,M>::computeLogOfNormalizationFactor()"
                             << ", m_logOfNormalizationFactor = " << m_logOfNormalizationFactor
                             << std::endl;
   }
 
   return value;
 }

template<class V , class M >

void QUESO::LogNormalJointPdf< V, M >::distributionMean ( V & meanVector ) const

virtual

Mean value of the underlying random variable.

Reimplemented from QUESO::BaseJointPdf< V, M >.

Definition at line 215 of file LogNormalJointPdf.C.

 {
   // FIXME - this is the mean of a non-truncated lognormal
   // distribution, and doesn't take into account a limited domainSet.
 
   if (m_diagonalCovMatrix) {
     unsigned int n_params = meanVector.sizeLocal();
     queso_assert_equal_to (n_params, this->lawExpVector().sizeLocal());
 
     for (unsigned int i = 0; i < n_params; ++i) {
       meanVector[i] = std::exp(this->lawExpVector()[i] + 0.5*this->lawVarVector()[i]);
     }
   }
   else {
     queso_error_msg("situation with a non-diagonal covariance matrix makes no sense");
   }
 }

template<class V , class M >

void QUESO::LogNormalJointPdf< V, M >::distributionVariance ( M & covMatrix ) const

virtual

Covariance matrix of the underlying random variable.

Reimplemented from QUESO::BaseJointPdf< V, M >.

Definition at line 235 of file LogNormalJointPdf.C.

 {
   // FIXME - this is the variance of a non-truncated lognormal
   // distribution, and doesn't take into account a limited domainSet.
 
   if (m_diagonalCovMatrix) {
     unsigned int n_params = this->lawExpVector().sizeLocal();
     queso_assert_equal_to (n_params, this->lawVarVector().sizeLocal());
     queso_assert_equal_to (n_params, covMatrix.numCols());
     queso_assert_equal_to (covMatrix.numCols(), covMatrix.numRowsGlobal());
 
     covMatrix.zeroLower();
     covMatrix.zeroUpper();
 
     for (unsigned int i = 0; i < n_params; ++i) {
       covMatrix(i,i) = (std::exp(this->lawVarVector()[i]) - 1) *
                        std::exp(2*this->lawExpVector()[i] + this->lawVarVector()[i]);
     }
   }
   else {
     queso_error_msg("situation with a non-diagonal covariance matrix makes no sense");
   }
 }

template<class V , class M >

const V & QUESO::LogNormalJointPdf< V, M >::lawExpVector ( ) const

Access to the vector of mean values and private attribute: m_lawExpVector.

Definition at line 74 of file LogNormalJointPdf.C.

Referenced by QUESO::LogNormalJointPdf< V, M >::LogNormalJointPdf().

 {
   return *m_lawExpVector;
 }

template<class V , class M >

const V & QUESO::LogNormalJointPdf< V, M >::lawVarVector ( ) const

Access to the vector of variance values and private attribute: m_lawVarVector.

Definition at line 81 of file LogNormalJointPdf.C.

Referenced by QUESO::LogNormalJointPdf< V, M >::LogNormalJointPdf().

 {
   return *m_lawVarVector;
 }

template<class V , class M >

double QUESO::LogNormalJointPdf< V, M >::lnValue	(	const V &	domainVector,
		const V *	domainDirection,
		V *	gradVector,
		M *	hessianMatrix,
		V *	hessianEffect
	)		const

virtual

Logarithm of the value of the Log-Normal PDF (scalar function).

The logarithm of the value of the Log-Normal density of diagonal covariance matrix (sigma^2) comes from the summation:

\[ lnValue =- \sum_i \frac{1}{domainVector_i * \sqrt{2 \pi * lawVarVector_i}} exp(-\frac{(\ln( domainVector_i) - lawExpVector_i)^2}{2 lawVarVector_i}) \]

as long as \( domainVector_i > 0 \) for all \( i \).

Reimplemented from QUESO::BaseScalarFunction< V, M >.

Definition at line 142 of file LogNormalJointPdf.C.

 {
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 55)) {
     *m_env.subDisplayFile() << "Entering LogNormalJointPdf<V,M>::lnValue()"
                             << ", meanVector = "   << *m_lawExpVector
                             << ": domainVector = " << domainVector
                             << std::endl;
   }
 
   queso_require_msg(!(domainDirection || hessianMatrix || hessianEffect), "incomplete code for gradVector, hessianMatrix and hessianEffect calculations");
 
   double returnValue = 0.;
 
   V zeroVector(domainVector);
   zeroVector.cwSet(0.);
   if (domainVector.atLeastOneComponentSmallerOrEqualThan(zeroVector)) {
     // What should the gradient be here?
     returnValue = -INFINITY;
   }
   else if (this->m_domainSet.contains(domainVector) == false) {
     // What should the gradient be here?
     returnValue = -INFINITY;
   }
   else {
     if (m_diagonalCovMatrix) {
       V diffVec(zeroVector);
       for (unsigned int i = 0; i < domainVector.sizeLocal(); ++i) {
         diffVec[i] = std::log(domainVector[i]) - this->lawExpVector()[i];
 
         // Compute the gradient of log of the PDF
         // The log of a log normal pdf is:
         // f(x) = -log(x \sigma sqrt(2 \pi)) - ((log(x) - \mu)^2 / (2 \sigma^2))
         // Therefore
         // \frac{df}{dx}(x) = -1/x - (log(x) - \mu) / (x \sigma^2)
         if (gradVector) {
           (*gradVector)[i] = -(1.0 / domainVector[i]) -
             diffVec[i] / (domainVector[i] * this->lawVarVector()[i]);
         }
       }
       returnValue = ((diffVec*diffVec)/this->lawVarVector()).sumOfComponents();
       returnValue *= -0.5;
 
       for (unsigned int i = 0; i < domainVector.sizeLocal(); ++i) {
         returnValue -= std::log(domainVector[i]);
         if (m_normalizationStyle == 0) {
           returnValue -= std::log(std::sqrt(2. * M_PI * this->lawVarVector()[i])); // Contribution of 1/(x\sqrt{2\pi\sigma^2})
         }
       }
     }
     else {
       queso_error_msg("situation with a non-diagonal covariance matrix makes no sense");
     }
     returnValue += m_logOfNormalizationFactor;
   }
 
   if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 55)) {
     *m_env.subDisplayFile() << "Leaving LogNormalJointPdf<V,M>::lnValue()"
                             << ", meanVector = "   << *m_lawExpVector
                             << ": domainVector = " << domainVector
                             << ", returnValue = "  << returnValue
                             << std::endl;
   }
 
   return returnValue;
 }

Member Data Documentation

template<class V = GslVector, class M = GslMatrix>

bool QUESO::LogNormalJointPdf< V, M >::m_diagonalCovMatrix

protected

Definition at line 103 of file LogNormalJointPdf.h.

template<class V = GslVector, class M = GslMatrix>

V* QUESO::LogNormalJointPdf< V, M >::m_lawExpVector

protected

Definition at line 101 of file LogNormalJointPdf.h.

template<class V = GslVector, class M = GslMatrix>

V* QUESO::LogNormalJointPdf< V, M >::m_lawVarVector

protected

Definition at line 102 of file LogNormalJointPdf.h.

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

src/stats/inc/LogNormalJointPdf.h
src/stats/src/LogNormalJointPdf.C

Public Member Functions

Protected Attributes

Additional Inherited Members

Detailed Description

template<class V = GslVector, class M = GslMatrix> class QUESO::LogNormalJointPdf< V, M >

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation

template<class V = GslVector, class M = GslMatrix>
class QUESO::LogNormalJointPdf< V, M >