IMP logo
IMP Reference Guide  2.8.0
The Integrative Modeling Platform
GaussianEMRestraint.h
Go to the documentation of this file.
1 /**
2  * \file IMP/isd/GaussianEMRestraint.h
3  * \brief Restrain two sets of Gaussians (model and GMM derived from EM map)
4  *
5  * Copyright 2007-2017 IMP Inventors. All rights reserved.
6  *
7  */
8 
9 #ifndef IMPISD_GAUSSIAN_EM_RESTRAINT_H
10 #define IMPISD_GAUSSIAN_EM_RESTRAINT_H
11 
12 #include "isd_config.h"
13 #include <IMP/PairContainer.h>
15 #include <IMP/container_macros.h>
16 #include <IMP/core/XYZ.h>
17 #include <IMP/core/Gaussian.h>
18 #include <IMP/algebra/Gaussian3D.h>
21 #include <IMP/atom/Mass.h>
22 #include <math.h>
23 #include <IMP/algebra/eigen3/Eigen/Dense>
24 #include <boost/unordered_map.hpp>
25 
26 IMPISD_BEGIN_NAMESPACE
27 
28 //! Restraint between two Gaussian Mixture Models, "model" and "density"
29 /** This restrains two sets of GMMs with a density overlap function.
30  The function is correlation of the two GMMs \f$f_M\f$ and \f$f_D\f$:
31  \f[
32  \frac{2\int{f_M(x)f_D(x)dx}}{\int{f_M^2(x)+f_D^2(x)dx}}
33  \f]
34  Where the integral is the "overlap function" given by:
35  \f[
36  ov(f_M,f_D) = \sum_{i=1}^{N_M} \sum_{j=1}^{N_D} \frac{1}{(2 \pi)^{3/2}|\Sigma_{Mi}+\Sigma_{Dj}|^{1/2}}\exp\left [-\frac{1}{2}(\boldsymbol\mu_{Mi} - \boldsymbol\mu_{Dj})^\top (\Sigma_{Mi}+\Sigma_{Dj})^{-1} (\boldsymbol\mu_{Mi} - \boldsymbol \mu_{Dj})\right ]
37  \f]
38  \note Source: Greenberg, Pellarin, Sali. In preparation.
39  */
40 class IMPISDEXPORT GaussianEMRestraint : public Restraint
41 {
42 
43  public:
44  //! Setup the GaussianEMRestraint
45  /**
46  \param[in] mdl the Model object to operate on
47  \param[in] model_ps particles for the model GMM
48  \param[in] density_ps particles for the density GMM
49  \param[in] global_sigma Particle to modulate the uncertainty
50  \param[in] model_cutoff_dist Cutoff for the model-model interactions
51  \param[in] density_cutoff_dist Cutoff for model-density interactions
52  \param[in] slope Gentle term to move all particles to the density
53  \param[in] update_model (DEPRECATED) update model each cycle
54  \param[in] backbone_slope Limit the slope only to backbone particles
55  \param[in] local Only consider density particles that are within the
56  specified cutoff of the model particles (experimental)
57  \param[in] name Name of this restraint
58  \note the model and density particles must be set up as Gaussian
59  */
61  ParticleIndexes model_ps, ParticleIndexes density_ps,
62  ParticleIndex global_sigma,
63  Float model_cutoff_dist, Float density_cutoff_dist,
64  Float slope,
65  bool update_model=true, bool backbone_slope=false,
66  bool local=false,
67  std::string name="GaussianEMRestraint%1%");
68 
69  //! Returns exp(score)
70  double get_probability() const {
71  return exp(-unprotected_evaluate(NULL));
72  }
73 
74  //! Get cross correlation between the model and the map
75  /** This CCC is that calculated from the last scoring function
76  evaluation; calling this function before the score is calculated
77  results in undefined behavior.
78  */
80  return cross_correlation_;
81  }
82 
83  //! Pre-calculate the density-density and model-model scores
84  /** This is automatically called by the constructor.
85  You only need to call it manually if you change Gaussian variances.
86  */
87  void compute_initial_scores();
88 
89  //! Set restraint slope
90  void set_slope(Float s){slope_=s;}
91 
92  //! Get restraint slope
93  Float get_slope(){return slope_;}
94  virtual double
95  unprotected_evaluate(IMP::DerivativeAccumulator *accum) const IMP_OVERRIDE;
97  void show(std::ostream &out) const { out << "GEM restraint"; }
99 
100  private:
101  ParticleIndexes model_ps_;
102  ParticleIndexes density_ps_;
103  ParticleIndex global_sigma_;
104  Float slope_;
105  bool update_model_,local_;
106  int msize_,dsize_;
107  Float normalization_;
108  Float dd_score_;
109  Float self_mm_score_;
112  ParticleIndexes slope_ps_; //experiment
113 
114  //variables needed to tabulate the exponential
115  Floats exp_grid_;
116  double invdx_;
117  double argmax_;
118 
119  mutable double cross_correlation_;
120 };
121 
122 IMPISD_END_NAMESPACE
123 
124 #endif /* IMPISD_GAUSSIAN_EM_RESTRAINT_H */
Declare an efficient stl-compatible map.
Decorator to hold Gaussian3D.
A decorator for particles with mass.
Store a list of ParticleIndexes.
Return all pairs from a SingletonContainer.
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Definition: object_macros.h:25
Float get_slope()
Get restraint slope.
void set_slope(Float s)
Set restraint slope.
Macros to define containers of objects.
Simple XYZ decorator.
Class for storing model, its restraints, constraints, and particles.
Definition: Model.h:72
double get_probability() const
Returns exp(score)
A container for Pairs.
double get_cross_correlation_coefficient() const
Get cross correlation between the model and the map.
Restraint between two Gaussian Mixture Models, "model" and "density".
Return all pairs from a SingletonContainer.
std::ostream & show(Hierarchy h, std::ostream &out=std::cout)
Print the hierarchy using a given decorator to display each node.
double Float
Basic floating-point value (could be float, double...)
Definition: types.h:20
Gaussian shape.
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:52