IMP logo
IMP Reference Guide  2.14.0
The Integrative Modeling Platform
GaussianProcessInterpolationRestraint.h
Go to the documentation of this file.
1 /**
2  * \file IMP/isd/GaussianProcessInterpolationRestraint.h
3  * \brief Restraint and ScoreState for GaussianProcessInterpolation
4  *
5  * Copyright 2007-2020 IMP Inventors. All rights reserved.
6  */
7 
8 #ifndef IMPISD_GAUSSIAN_PROCESS_INTERPOLATION_RESTRAINT_H
9 #define IMPISD_GAUSSIAN_PROCESS_INTERPOLATION_RESTRAINT_H
10 
11 #include <IMP/isd/isd_config.h>
12 #include <IMP/macros.h>
13 #include <boost/scoped_ptr.hpp>
14 #include <IMP/Restraint.h>
17 #include <IMP/Pointer.h>
18 #include <Eigen/Dense>
19 
20 #include <IMP/ScoreState.h>
21 
22 IMPISD_BEGIN_NAMESPACE
23 
24 class GaussianProcessInterpolationScoreState;
25 
26 //! gaussian process restraint
27 /* the restraint is a multivariate normal distribution on the vector of
28 * observations with mean and standard deviation given by the posterior of the
29 * gaussian process.
30 */
32  : public Restraint {
33  private:
34  // checks and makes necessary updates
35  void update_mean_and_covariance();
39  // number of observation points
40  unsigned M_;
41 
42  void create_score_state();
43 
44  public:
45  /** This is a restraint on other restraints. It first constructs the
46  necessary vectors from GaussianProcessInterpolation, then creates a
47  multivariate normal distribution around it. Upon evaluation, it
48  checks if parameters have changed, reconstructs the matrix if
49  necessary, changes the DA weight and passes it to the functions. */
52 
53  /** To call this, you need to update the scorestate before.
54  calling model.evaluate(False) is enough. */
55  double get_probability() const { return mvn_->density(); }
56 
57  //! Use conjugate gradients when possible (default false)
58  void set_use_cg(bool use, double tol) { mvn_->set_use_cg(use, tol); }
59 
60  //! Get minus log normalization and minus exponent separately
61  double get_minus_log_normalization() const;
62  double get_minus_exponent() const;
63 
64  //! Get hessian of the minus log likelihood
65  Eigen::MatrixXd get_hessian() const;
66 
67  //! Get log determinant of hessian
68  double get_logdet_hessian() const;
69 
70  //! call this one from Python
71  FloatsList get_hessian(bool unused) const;
72 
73  public:
74  double unprotected_evaluate(IMP::DerivativeAccumulator *accum) const
78 
79  // to allow the scorestate to get the restraint's objects
80  friend class GaussianProcessInterpolationScoreState;
81 };
82 
83 #if !defined(IMP_DOXYGEN) && !defined(SWIG)
84 class IMPISDEXPORT GaussianProcessInterpolationScoreState : public ScoreState {
85  private:
87 
88  private:
89  GaussianProcessInterpolationScoreState(
90  GaussianProcessInterpolationRestraint *gpir)
91  : ScoreState(gpir->get_model(),
92  "GaussianProcessInterpolationScoreState%1%"),
93  gpir_(gpir) {}
94 
95  public:
96  // only the GPIR can create this and add it to the model
97  friend class GaussianProcessInterpolationRestraint;
98  virtual void do_before_evaluate() IMP_OVERRIDE;
99  virtual void do_after_evaluate(DerivativeAccumulator *da) IMP_OVERRIDE;
102  IMP_OBJECT_METHODS(GaussianProcessInterpolationScoreState);
103 };
104 #endif
105 
106 IMPISD_END_NAMESPACE
107 
108 #endif /* IMPISD_GAUSSIAN_PROCESS_INTERPOLATION_RESTRAINT_H */
Normal distribution of Function.
Smart pointer to Object-derived classes that does not refcount.
Definition: WeakPointer.h:76
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Definition: object_macros.h:25
void set_use_cg(bool use, double tol)
Use conjugate gradients when possible (default false)
Various general useful macros for IMP.
Normal distribution of Function.
Class for storing model, its restraints, constraints, and particles.
Definition: Model.h:72
ScoreStates maintain invariants in the Model.
Definition: ScoreState.h:53
A smart pointer to a ref-counted Object that is a class member.
Definition: Pointer.h:146
Shared score state.
virtual ModelObjectsTemp do_get_outputs() const =0
A nullptr-initialized pointer to an IMP Object.
Abstract base class for all restraints.
virtual ModelObjectsTemp do_get_inputs() const =0
#define IMP_OVERRIDE
Cause a compile error if this method does not override a parent method.
Class for adding derivatives from restraints to the model.
A restraint is a term in an IMP ScoringFunction.
Definition: Restraint.h:54