MUQ  0.4.3
muq::SamplingAlgorithms::DILIKernel Class Reference

An implementation of the Dimension Independent Likelihood Informed subspace (DILI) MCMC sampler. More...

#include <DILIKernel.h>

Inheritance diagram for muq::SamplingAlgorithms::DILIKernel:

Detailed Description

An implementation of the Dimension Independent Likelihood Informed subspace (DILI) MCMC sampler.

Cui, T., Law, K. J., & Marzouk, Y. M. (2016). Dimension-independent likelihood-informed MCMC. Journal of Computational Physics, 304, 109-137.

Configuration Parameters:

Parameter Key Type Default Value Description
"LIS Block" String - A string pointing to a block of kernel/proposal options for the Likelihood informed subspace.
"CS Block" String - A string pointing to a block of kernel/proposal options for the Complementary space. Typically this will be a Crank-Nicolson proposal.
"HessianType" String GaussNewton The type of posterior Hessian to use. Either "Exact" or "GaussNewton"
"Adapt Interval" int -1 How many MCMC steps to take before updating the average Hessian and the likelihood informed subspace. If negative, the LIS will be fixed to the initial value and will not be udpated.
"Adapt Start" int 1 The number of MCMC steps taken before updating the LIS.
"Adapt End" int -1 No LIS updates will occur after this number of MCMC steps. If negative, the LIS will continue to be updated until the end of the chain.
"Initial Weight" int 100 "Weight" or number of samples given to the to initial Hessian. The weight on the previous average Hessian estimate is given by $(N+W)/(N+W+1)$, where $N$ is the number of MCMC steps taken and $W$ is this parameter.
"Hess Tolerance" double 1e-4 The lower bound on eigenvalues used to construct the low-rank approximation of the global expected Hessian.
"LIS Tolerance" double 0.1 Lower bound on eigenvalues used to construct local and global likelihood informed subspaces. Must be larger than Hess Tolerance.
"Eigensolver Block" String - A string pointing to a block of eigensolver options for solving the generalized eigenvalue problems.

Definition at line 141 of file DILIKernel.h.

Public Member Functions

 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem)
 
 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem, Eigen::VectorXd const &genEigVals, Eigen::MatrixXd const &genEigVecs)
 
 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem, std::shared_ptr< muq::Modeling::GaussianBase > const &prior, std::shared_ptr< muq::Modeling::ModPiece > const &noiseModel, std::shared_ptr< muq::Modeling::ModPiece > const &forwardModelIn)
 
 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem, std::shared_ptr< muq::Modeling::GaussianBase > const &prior, std::shared_ptr< muq::Modeling::ModPiece > const &noiseModel, std::shared_ptr< muq::Modeling::ModPiece > const &forwardModelIn, Eigen::VectorXd const &genEigVals, Eigen::MatrixXd const &genEigVecs)
 
 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem, std::shared_ptr< muq::Modeling::GaussianBase > const &prior, std::shared_ptr< muq::Modeling::ModPiece > const &likelihood)
 
 DILIKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem, std::shared_ptr< muq::Modeling::GaussianBase > const &prior, std::shared_ptr< muq::Modeling::ModPiece > const &likelihood, Eigen::VectorXd const &genEigVals, Eigen::MatrixXd const &genEigVecs)
 
virtual ~DILIKernel ()=default
 
virtual std::shared_ptr< TransitionKernelLISKernel ()
 
virtual std::shared_ptr< TransitionKernelCSKernel ()
 
virtual void PostStep (unsigned int const t, std::vector< std::shared_ptr< SamplingState >> const &state) override
 Allow the kernel to adapt given a new state. More...
 
virtual std::vector< std::shared_ptr< SamplingState > > Step (unsigned int const t, std::shared_ptr< SamplingState > prevState) override
 
virtual void PrintStatus (std::string prefix) const override
 
Eigen::MatrixXd const & LISVecs () const
 
Eigen::VectorXd const & LISVals () const
 
Eigen::MatrixXd const & LISW () const
 
Eigen::VectorXd const & LISL () const
 
unsigned int LISDim () const
 
void CreateLIS (std::vector< Eigen::VectorXd > const &currState)
 
Eigen::VectorXd ToLIS (Eigen::VectorXd const &x) const
 
Eigen::VectorXd FromLIS (Eigen::VectorXd const &r) const
 
Eigen::VectorXd ToCS (Eigen::VectorXd const &x) const
 
- Public Member Functions inherited from muq::SamplingAlgorithms::TransitionKernel
 TransitionKernel (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem)
 
virtual ~TransitionKernel ()=default
 
virtual void SetCommunicator (std::shared_ptr< parcer::Communicator > newcomm)
 
virtual void PreStep (unsigned int const t, std::shared_ptr< SamplingState > state)
 Allow the kernel to preprocess the current step. More...
 
virtual void PrintStatus () const
 
virtual void SetBlockInd (int newBlockInd)
 
virtual int GetBlockInd () const
 
virtual std::shared_ptr< AbstractSamplingProblem > const & Problem () const
 

Static Public Member Functions

static std::shared_ptr< muq::Modeling::ModPieceExtractLikelihood (std::shared_ptr< AbstractSamplingProblem > const &problem, std::string const &nodeName)
 
static std::shared_ptr< muq::Modeling::GaussianBaseExtractPrior (std::shared_ptr< AbstractSamplingProblem > const &problem, std::string const &nodeName)
 
static std::shared_ptr< muq::Modeling::ModPieceExtractNoiseModel (std::shared_ptr< muq::Modeling::ModPiece > const &likelihood)
 
static std::shared_ptr< muq::Modeling::ModPieceExtractForwardModel (std::shared_ptr< muq::Modeling::ModPiece > const &likelihoodIn)
 
static std::shared_ptr< muq::Modeling::ModPieceCreateLikelihood (std::shared_ptr< muq::Modeling::ModPiece > const &forwardModel, std::shared_ptr< muq::Modeling::ModPiece > const &noiseDensity)
 
- Static Public Member Functions inherited from muq::SamplingAlgorithms::TransitionKernel
static std::shared_ptr< TransitionKernelConstruct (boost::property_tree::ptree const &pt, std::shared_ptr< AbstractSamplingProblem > problem)
 Static constructor for the transition kernel. More...
 
static std::shared_ptr< TransitionKernelMapGetTransitionKernelMap ()
 

Additional Inherited Members

- Public Types inherited from muq::SamplingAlgorithms::TransitionKernel
typedef boost::function< std::shared_ptr< TransitionKernel >boost::property_tree::ptree, std::shared_ptr< AbstractSamplingProblem >)> TransitionKernelConstructor
 
typedef std::map< std::string, TransitionKernelConstructorTransitionKernelMap
 

Constructor & Destructor Documentation

◆ DILIKernel() [1/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem 
)

Definition at line 71 of file DILIKernel.cpp.

◆ DILIKernel() [2/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem,
Eigen::VectorXd const &  genEigVals,
Eigen::MatrixXd const &  genEigVecs 
)

Definition at line 79 of file DILIKernel.cpp.

◆ DILIKernel() [3/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem,
std::shared_ptr< muq::Modeling::GaussianBase > const &  prior,
std::shared_ptr< muq::Modeling::ModPiece > const &  noiseModel,
std::shared_ptr< muq::Modeling::ModPiece > const &  forwardModelIn 
)

Definition at line 92 of file DILIKernel.cpp.

References eigOpts.

◆ DILIKernel() [4/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem,
std::shared_ptr< muq::Modeling::GaussianBase > const &  prior,
std::shared_ptr< muq::Modeling::ModPiece > const &  noiseModel,
std::shared_ptr< muq::Modeling::ModPiece > const &  forwardModelIn,
Eigen::VectorXd const &  genEigVals,
Eigen::MatrixXd const &  genEigVecs 
)

Definition at line 119 of file DILIKernel.cpp.

References eigOpts, and SetLIS().

◆ DILIKernel() [5/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem,
std::shared_ptr< muq::Modeling::GaussianBase > const &  prior,
std::shared_ptr< muq::Modeling::ModPiece > const &  likelihood 
)

Definition at line 150 of file DILIKernel.cpp.

References eigOpts.

◆ DILIKernel() [6/6]

DILIKernel::DILIKernel ( boost::property_tree::ptree const &  pt,
std::shared_ptr< AbstractSamplingProblem problem,
std::shared_ptr< muq::Modeling::GaussianBase > const &  prior,
std::shared_ptr< muq::Modeling::ModPiece > const &  likelihood,
Eigen::VectorXd const &  genEigVals,
Eigen::MatrixXd const &  genEigVecs 
)

Definition at line 176 of file DILIKernel.cpp.

References eigOpts, and SetLIS().

◆ ~DILIKernel()

virtual muq::SamplingAlgorithms::DILIKernel::~DILIKernel ( )
virtualdefault

Member Function Documentation

◆ ComputeLocalLIS()

std::pair< Eigen::MatrixXd, Eigen::VectorXd > DILIKernel::ComputeLocalLIS ( std::vector< Eigen::VectorXd > const &  currState)
protected

◆ CreateLikelihood()

std::shared_ptr< muq::Modeling::ModPiece > DILIKernel::CreateLikelihood ( std::shared_ptr< muq::Modeling::ModPiece > const &  forwardModel,
std::shared_ptr< muq::Modeling::ModPiece > const &  noiseDensity 
)
static

◆ CreateLIS()

void DILIKernel::CreateLIS ( std::vector< Eigen::VectorXd > const &  currState)

Create the likelihood informed subspace for the first time.

Definition at line 372 of file DILIKernel.cpp.

References ComputeLocalLIS(), numLisUpdates, and SetLIS().

Referenced by Step().

◆ CSKernel()

virtual std::shared_ptr<TransitionKernel> muq::SamplingAlgorithms::DILIKernel::CSKernel ( )
inlinevirtual

Definition at line 182 of file DILIKernel.h.

References csKernel.

◆ ExtractForwardModel()

std::shared_ptr< muq::Modeling::ModPiece > DILIKernel::ExtractForwardModel ( std::shared_ptr< muq::Modeling::ModPiece > const &  likelihoodIn)
static

Definition at line 481 of file DILIKernel.cpp.

◆ ExtractLikelihood()

std::shared_ptr< muq::Modeling::ModPiece > DILIKernel::ExtractLikelihood ( std::shared_ptr< AbstractSamplingProblem > const &  problem,
std::string const &  nodeName 
)
static

From a ModGraphPiece defining the posterior log density, this function extracts a ModGraphPiece defining the likelihood function.

Definition at line 498 of file DILIKernel.cpp.

References muq::SamplingAlgorithms::TransitionKernel::problem.

◆ ExtractNoiseModel()

std::shared_ptr< muq::Modeling::ModPiece > DILIKernel::ExtractNoiseModel ( std::shared_ptr< muq::Modeling::ModPiece > const &  likelihood)
static

Definition at line 471 of file DILIKernel.cpp.

◆ ExtractPrior()

std::shared_ptr< muq::Modeling::GaussianBase > DILIKernel::ExtractPrior ( std::shared_ptr< AbstractSamplingProblem > const &  problem,
std::string const &  nodeName 
)
static

From a ModGraphPiece defining the posterior log density, this function extracts a GaussianBase instance defining the prior.

Definition at line 441 of file DILIKernel.cpp.

References muq::SamplingAlgorithms::TransitionKernel::problem.

◆ FromLIS()

Eigen::VectorXd DILIKernel::FromLIS ( Eigen::VectorXd const &  r) const

Returns a vector in the full space given a vector \(r\) in the low dimensional LIS.

Definition at line 223 of file DILIKernel.cpp.

References hessU, lisDim, and lisL.

Referenced by Step().

◆ LISDim()

unsigned int muq::SamplingAlgorithms::DILIKernel::LISDim ( ) const
inline

Returns the dimension of the LIS.

Definition at line 219 of file DILIKernel.h.

References lisDim.

◆ LISKernel()

virtual std::shared_ptr<TransitionKernel> muq::SamplingAlgorithms::DILIKernel::LISKernel ( )
inlinevirtual

Definition at line 181 of file DILIKernel.h.

References lisKernel.

◆ LISL()

Eigen::VectorXd const& muq::SamplingAlgorithms::DILIKernel::LISL ( ) const
inline

Definition at line 216 of file DILIKernel.h.

References lisL.

◆ LISVals()

Eigen::VectorXd const& muq::SamplingAlgorithms::DILIKernel::LISVals ( ) const
inline

Definition at line 214 of file DILIKernel.h.

References hessEigVals.

◆ LISVecs()

Eigen::MatrixXd const& muq::SamplingAlgorithms::DILIKernel::LISVecs ( ) const
inline

Definition at line 213 of file DILIKernel.h.

References hessU.

◆ LISW()

Eigen::MatrixXd const& muq::SamplingAlgorithms::DILIKernel::LISW ( ) const
inline

Definition at line 215 of file DILIKernel.h.

References hessW.

◆ PostStep()

void DILIKernel::PostStep ( unsigned int const  t,
std::vector< std::shared_ptr< SamplingState >> const &  state 
)
overridevirtual

Allow the kernel to adapt given a new state.

By default this function does nothing but children can override it to adapt the kernel

Parameters
[in]tThe current step
[in]stateThe current state

Reimplemented from muq::SamplingAlgorithms::TransitionKernel.

Definition at line 206 of file DILIKernel.cpp.

References adaptEnd, adaptStart, numLisUpdates, updateInterval, and UpdateLIS().

◆ PrintStatus()

void DILIKernel::PrintStatus ( std::string  prefix) const
overridevirtual

Reimplemented from muq::SamplingAlgorithms::TransitionKernel.

Definition at line 295 of file DILIKernel.cpp.

References csKernel, lisDim, and lisKernel.

◆ SetLIS()

void DILIKernel::SetLIS ( Eigen::VectorXd const &  eigVals,
Eigen::MatrixXd const &  eigVecs 
)
protected

Sets up the LIS based on the eigenvalues and eigenvectors solving the problem $Hu = \lambda \Gamma^{-1}u$

Definition at line 305 of file DILIKernel.cpp.

References hessEigVals, hessU, hessUQR, hessW, lisDim, lisL, lisValTol, prior, and UpdateKernels().

Referenced by CreateLIS(), DILIKernel(), and UpdateLIS().

◆ Step()

std::vector< std::shared_ptr< SamplingState > > DILIKernel::Step ( unsigned int const  t,
std::shared_ptr< SamplingState prevState 
)
overridevirtual
Parameters
[in]tThe current step
[in]stateThe current state

Implements muq::SamplingAlgorithms::TransitionKernel.

Definition at line 235 of file DILIKernel.cpp.

References muq::SamplingAlgorithms::TransitionKernel::blockInd, CreateLIS(), csKernel, FromLIS(), hessU, lisDim, lisKernel, ToCS(), and ToLIS().

◆ ToCS()

Eigen::VectorXd DILIKernel::ToCS ( Eigen::VectorXd const &  x) const

Projects a point \(x\in\mathbb{R}^N\) in the full parameter space to a point \(y\in\mathbb{R}^N\) that has the same dimension as \(x\) but lies on the complementary space.

Definition at line 229 of file DILIKernel.cpp.

References hessU, hessW, and lisDim.

Referenced by Step().

◆ ToLIS()

Eigen::VectorXd DILIKernel::ToLIS ( Eigen::VectorXd const &  x) const

Returns the low dimensional LIS representation of a vector \(x\) defined in the full parameter space.

Definition at line 215 of file DILIKernel.cpp.

References hessW, lisDim, and lisL.

Referenced by Step().

◆ UpdateKernels()

void DILIKernel::UpdateKernels ( )
protected

◆ UpdateLIS()

void DILIKernel::UpdateLIS ( unsigned int  numSamps,
std::vector< Eigen::VectorXd > const &  currState 
)
protected

Update the likelihood informed subspace using Hessian information at a new point.

Parameters
[in]numSampsThe number of samples that have been used to estimate the current LIS.
[in]currStateThe current state of the chain.

Definition at line 413 of file DILIKernel.cpp.

References muq::SamplingAlgorithms::TransitionKernel::blockInd, muq::Modeling::StochasticEigenSolver::compute(), muq::Modeling::GeneralizedEigenSolver::eigenvalues(), muq::Modeling::GeneralizedEigenSolver::eigenvectors(), eigOpts, forwardModel, hessEigVals, hessType, hessUQR, hessValTol, hessW, initialHessSamps, logLikelihood, noiseDensity, prior, and SetLIS().

Referenced by PostStep().

Member Data Documentation

◆ adaptEnd

const int muq::SamplingAlgorithms::DILIKernel::adaptEnd
protected

Definition at line 307 of file DILIKernel.h.

Referenced by PostStep().

◆ adaptStart

const int muq::SamplingAlgorithms::DILIKernel::adaptStart
protected

Definition at line 306 of file DILIKernel.h.

Referenced by PostStep().

◆ csKernel

std::shared_ptr<TransitionKernel> muq::SamplingAlgorithms::DILIKernel::csKernel
protected

Definition at line 297 of file DILIKernel.h.

Referenced by CSKernel(), PrintStatus(), Step(), and UpdateKernels().

◆ csKernelOpts

boost::property_tree::ptree muq::SamplingAlgorithms::DILIKernel::csKernelOpts
protected

Definition at line 264 of file DILIKernel.h.

Referenced by UpdateKernels().

◆ eigOpts

boost::property_tree::ptree muq::SamplingAlgorithms::DILIKernel::eigOpts
protected

Definition at line 303 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), DILIKernel(), and UpdateLIS().

◆ forwardModel

std::shared_ptr<muq::Modeling::ModPiece> muq::SamplingAlgorithms::DILIKernel::forwardModel
protected

Definition at line 269 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), CreateLikelihood(), and UpdateLIS().

◆ fullToCS

std::shared_ptr<muq::Modeling::LinearOperator> muq::SamplingAlgorithms::DILIKernel::fullToCS
protected

Definition at line 291 of file DILIKernel.h.

Referenced by UpdateKernels().

◆ hessEigVals

std::shared_ptr<Eigen::VectorXd> muq::SamplingAlgorithms::DILIKernel::hessEigVals
protected

Definition at line 279 of file DILIKernel.h.

Referenced by LISVals(), SetLIS(), and UpdateLIS().

◆ hessType

std::string muq::SamplingAlgorithms::DILIKernel::hessType
protected

Definition at line 300 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), and UpdateLIS().

◆ hessU

std::shared_ptr<Eigen::MatrixXd> muq::SamplingAlgorithms::DILIKernel::hessU
protected

Definition at line 273 of file DILIKernel.h.

Referenced by FromLIS(), LISVecs(), SetLIS(), Step(), ToCS(), and UpdateKernels().

◆ hessUQR

std::shared_ptr<Eigen::ColPivHouseholderQR<Eigen::MatrixXd> > muq::SamplingAlgorithms::DILIKernel::hessUQR
protected

Definition at line 276 of file DILIKernel.h.

Referenced by SetLIS(), and UpdateLIS().

◆ hessValTol

const double muq::SamplingAlgorithms::DILIKernel::hessValTol
protected

Definition at line 321 of file DILIKernel.h.

Referenced by UpdateLIS().

◆ hessW

std::shared_ptr<Eigen::MatrixXd> muq::SamplingAlgorithms::DILIKernel::hessW
protected

Definition at line 282 of file DILIKernel.h.

Referenced by LISW(), SetLIS(), ToCS(), ToLIS(), UpdateKernels(), and UpdateLIS().

◆ initialHessSamps

const int muq::SamplingAlgorithms::DILIKernel::initialHessSamps
protected

Definition at line 308 of file DILIKernel.h.

Referenced by UpdateLIS().

◆ lisDim

unsigned int muq::SamplingAlgorithms::DILIKernel::lisDim =0
protected

◆ lisKernel

std::shared_ptr<TransitionKernel> muq::SamplingAlgorithms::DILIKernel::lisKernel
protected

Definition at line 294 of file DILIKernel.h.

Referenced by LISKernel(), PrintStatus(), Step(), and UpdateKernels().

◆ lisKernelOpts

boost::property_tree::ptree muq::SamplingAlgorithms::DILIKernel::lisKernelOpts
protected

Definition at line 263 of file DILIKernel.h.

Referenced by UpdateKernels().

◆ lisL

std::shared_ptr<Eigen::VectorXd> muq::SamplingAlgorithms::DILIKernel::lisL
protected

Definition at line 285 of file DILIKernel.h.

Referenced by FromLIS(), LISL(), SetLIS(), ToLIS(), and UpdateKernels().

◆ lisToFull

std::shared_ptr<muq::Modeling::LinearOperator> muq::SamplingAlgorithms::DILIKernel::lisToFull
protected

Definition at line 288 of file DILIKernel.h.

Referenced by UpdateKernels().

◆ lisValTol

const double muq::SamplingAlgorithms::DILIKernel::lisValTol
protected

Definition at line 327 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), and SetLIS().

◆ logLikelihood

std::shared_ptr<muq::Modeling::ModPiece> muq::SamplingAlgorithms::DILIKernel::logLikelihood
protected

Definition at line 266 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), UpdateKernels(), and UpdateLIS().

◆ noiseDensity

std::shared_ptr<muq::Modeling::ModPiece> muq::SamplingAlgorithms::DILIKernel::noiseDensity
protected

Definition at line 270 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), CreateLikelihood(), and UpdateLIS().

◆ numLisUpdates

unsigned int muq::SamplingAlgorithms::DILIKernel::numLisUpdates
protected

Definition at line 329 of file DILIKernel.h.

Referenced by CreateLIS(), and PostStep().

◆ prior

std::shared_ptr<muq::Modeling::GaussianBase> muq::SamplingAlgorithms::DILIKernel::prior
protected

Definition at line 267 of file DILIKernel.h.

Referenced by ComputeLocalLIS(), SetLIS(), UpdateKernels(), and UpdateLIS().

◆ updateInterval

const int muq::SamplingAlgorithms::DILIKernel::updateInterval
protected

Definition at line 305 of file DILIKernel.h.

Referenced by PostStep().


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