queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
QUESO::GaussianLikelihoodFullCovarianceRandomCoefficient< V, M > Class Template Reference

A class that represents a Gaussian likelihood with full covariance and random coefficient. More...

#include <GaussianLikelihoodFullCovarianceRandomCoefficient.h>

Inheritance diagram for QUESO::GaussianLikelihoodFullCovarianceRandomCoefficient< V, M >:
QUESO::LikelihoodBase< V, M > QUESO::BaseScalarFunction< V, M >

Public Member Functions

virtual double lnValue (const V &domainVector) const
 Logarithm of the value of the scalar function. More...
 
Constructor/Destructor methods.
 GaussianLikelihoodFullCovarianceRandomCoefficient (const char *prefix, const VectorSet< V, M > &domainSet, const V &observations, const M &covariance)
 Default constructor. More...
 
virtual ~GaussianLikelihoodFullCovarianceRandomCoefficient ()
 Destructor. More...
 
- Public Member Functions inherited from QUESO::LikelihoodBase< V, M >
virtual void evaluateModel (const V &domainVector, const V *domainDirection, V &modelOutput, V *gradVector, M *hessianMatrix, V *hessianEffect) const
 Deprecated. Evaluates the user's model at the point domainVector. More...
 
virtual void evaluateModel (const V &domainVector, V &modelOutput) const
 Evaluates the user's model at the point domainVector. More...
 
virtual double actualValue (const V &domainVector, const V *, V *, M *, V *) const
 Actual value of the scalar function. More...
 
 LikelihoodBase (const char *prefix, const VectorSet< V, M > &domainSet, const V &observations)
 Default constructor. More...
 
virtual ~LikelihoodBase ()=0
 Destructor, pure to make this class abstract. 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 V *domainDirection, V *gradVector, M *hessianMatrix, V *hessianEffect) const
 Logarithm of the value of the scalar function. Deprecated. 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
 

Private Attributes

const M & m_covariance
 

Additional Inherited Members

- Protected Attributes inherited from QUESO::LikelihoodBase< V, M >
const V & m_observations
 
- Protected Attributes inherited from QUESO::BaseScalarFunction< V, M >
const BaseEnvironmentm_env
 
std::string m_prefix
 
const VectorSet< V, M > & m_domainSet
 Domain set of the scalar function. More...
 

Detailed Description

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

A class that represents a Gaussian likelihood with full covariance and random coefficient.

The random coefficient is a scalar that pre-multiplies the covariance matrix. This is treated as a hyperparameter to be inferred during the sampling procedure.

Definition at line 48 of file GaussianLikelihoodFullCovarianceRandomCoefficient.h.

Constructor & Destructor Documentation

template<class V , class M >
QUESO::GaussianLikelihoodFullCovarianceRandomCoefficient< V, M >::GaussianLikelihoodFullCovarianceRandomCoefficient ( const char *  prefix,
const VectorSet< V, M > &  domainSet,
const V &  observations,
const M &  covariance 
)

Default constructor.

Instantiates a Gaussian likelihood function, given a prefix, its domain, a set of observations and a full covariance matrix. The full covariance matrix is stored as a matrix in the covariance parameter.

The parameter covarianceCoefficient is a multiplying factor of covaraince and is treated as a random variable (i.e. it is solved for in a statistical inversion).

Definition at line 34 of file GaussianLikelihoodFullCovarianceRandomCoefficient.C.

37  : LikelihoodBase<V, M>(prefix, domainSet, observations),
38  m_covariance(covariance)
39 {
40  if (covariance.numRowsLocal() != observations.sizeLocal()) {
41  queso_error_msg("Covariance matrix not same size as observation vector");
42  }
43 }
const VectorSet< V, M > & domainSet() const
Access to the protected attribute m_domainSet: domain set of the scalar function. ...

Destructor.

Definition at line 46 of file GaussianLikelihoodFullCovarianceRandomCoefficient.C.

47 {
48 }

Member Function Documentation

template<class V , class M >
double QUESO::GaussianLikelihoodFullCovarianceRandomCoefficient< V, M >::lnValue ( const V &  domainVector) const
virtual

Logarithm of the value of the scalar function.

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

Definition at line 52 of file GaussianLikelihoodFullCovarianceRandomCoefficient.C.

53 {
54  V modelOutput(this->m_observations, 0, 0); // At least it's not a copy
55  V weightedMisfit(this->m_observations, 0, 0); // At least it's not a copy
56 
57  this->evaluateModel(domainVector, modelOutput);
58 
59  // Compute misfit G(x) - y
60  modelOutput -= this->m_observations;
61 
62  // Solve \Sigma u = G(x) - y for u
63  this->m_covariance.invertMultiply(modelOutput, weightedMisfit);
64 
65  // Compute (G(x) - y)^T \Sigma^{-1} (G(x) - y)
66  modelOutput *= weightedMisfit;
67 
68  // This is square of 2-norm
69  double norm2_squared = modelOutput.sumOfComponents();
70 
71  // Get the determinant of the covariance matrix |\Sigma|
72  double deter_cov = this->m_covariance.determinant();
73 
74  deter_cov = std::sqrt(deter_cov);
75 
76  // Set the right hyperparameter coefficient
77  // The last element of domainVector is the multiplicative coefficient of the
78  // covariance matrix
79  double cov_coeff = domainVector[domainVector.sizeLocal()-1];
80  cov_coeff = std::pow(std::sqrt(cov_coeff), this->m_observations.sizeLocal());
81 
82  return -0.5 * norm2_squared / cov_coeff - std::log(cov_coeff * deter_cov);
83 }
virtual void evaluateModel(const V &domainVector, const V *domainDirection, V &modelOutput, V *gradVector, M *hessianMatrix, V *hessianEffect) const
Deprecated. Evaluates the user&#39;s model at the point domainVector.

Member Data Documentation

template<class V = GslVector, class M = GslMatrix>
const M& QUESO::GaussianLikelihoodFullCovarianceRandomCoefficient< V, M >::m_covariance
private

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

Generated on Tue Jun 5 2018 19:49:26 for queso-0.57.1 by  doxygen 1.8.5