IMP  2.2.1
The Integrative Modeling Platform
GaussianProcessInterpolationRestraint.h
Go to the documentation of this file.
1 /**
2  * \file IMP/isd/GaussianProcessInterpolationRestraint.h
3  * \brief kernel::Restraint and ScoreState for GaussianProcessInterpolation
4  *
5  * Copyright 2007-2014 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/kernel/Restraint.h>
17 #include <IMP/base/Pointer.h>
18 #include <IMP/algebra/eigen3/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 kernel::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  void stats() const { return mvn_->stats(); }
58 
59  // use conjugate gradients when possible (default false)
60  void set_use_cg(bool use, double tol) { mvn_->set_use_cg(use, tol); }
61 
62  // get minus log normalization and minus exponent separately
63  double get_minus_log_normalization() const;
64  double get_minus_exponent() const;
65 
66  // get hessian of the minus log likelihood
67  IMP_Eigen::MatrixXd get_hessian() const;
68 
69  // get log determinant of hessian
70  double get_logdet_hessian() const;
71 
72  // call this one from python
73  FloatsList get_hessian(bool unused) const;
74 
75  public:
76  double unprotected_evaluate(IMP::DerivativeAccumulator *accum) const
77  IMP_OVERRIDE;
78  IMP::kernel::ModelObjectsTemp do_get_inputs() const IMP_OVERRIDE;
80  ;
81 
82  // to allow the scorestate to get the restraint's objects
83  friend class GaussianProcessInterpolationScoreState;
84 };
85 
86 #if !defined(IMP_DOXYGEN) && !defined(SWIG)
87 class IMPISDEXPORT GaussianProcessInterpolationScoreState : public ScoreState {
88  private:
90 
91  private:
92  GaussianProcessInterpolationScoreState(
94  : ScoreState(gpir->get_model(),
95  "GaussianProcessInterpolationScoreState%1%"),
96  gpir_(gpir) {}
97 
98  public:
99  // only the GPIR can create this and add it to the model
100  friend class GaussianProcessInterpolationRestraint;
101  virtual void do_before_evaluate() IMP_OVERRIDE;
102  virtual void do_after_evaluate(DerivativeAccumulator *da) IMP_OVERRIDE;
103  virtual kernel::ModelObjectsTemp do_get_inputs() const IMP_OVERRIDE;
104  virtual kernel::ModelObjectsTemp do_get_outputs() const IMP_OVERRIDE;
105  IMP_OBJECT_METHODS(GaussianProcessInterpolationScoreState);
106 };
107 #endif
108 
109 IMPISD_END_NAMESPACE
110 
111 #endif /* IMPISD_GAUSSIAN_PROCESS_INTERPOLATION_RESTRAINT_H */
Normal distribution of Function.
Class for adding derivatives from restraints to the model.
IMP::base::Vector< IMP::base::WeakPointer< kernel::ModelObject > > ModelObjectsTemp
A nullptr-initialized pointer to an IMP Object.
A smart pointer to a ref-counted Object that is a class memeber.
Definition: base/Pointer.h:147
A smart pointer to a reference counted object.
Definition: base/Pointer.h:87
Import IMP/kernel/macros.h in the namespace.
ScoreStates maintian invariants in the Model.
Normal distribution of Function.
Abstract base class for all restraints.
A restraint is a term in an IMP ScoringFunction.
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Import IMP/kernel/ScoreState.h in the namespace.
virtual ModelObjectsTemp do_get_inputs() const =0
Class for storing model, its restraints, constraints, and particles.
Definition: kernel/Model.h:72