From a4faf9cb1288bd1d95d112788f1a3bf3ac41edf6 Mon Sep 17 00:00:00 2001
From: Gabriel Studer <gabriel.studer@unibas.ch>
Date: Mon, 10 Aug 2020 22:53:35 +0200
Subject: [PATCH] allow to set linear VINA weights individually

---
 sidechain/doc/rotamer.rst                     | 39 ++++++++++++++++++-
 sidechain/pymod/export_rotamer_constructor.cc | 12 ++++++
 sidechain/src/vina_particle_scoring.cc        | 36 +++++++++++------
 sidechain/src/vina_particle_scoring.hh        | 17 ++++++++
 4 files changed, 91 insertions(+), 13 deletions(-)

diff --git a/sidechain/doc/rotamer.rst b/sidechain/doc/rotamer.rst
index 7e87093f..c729b48e 100644
--- a/sidechain/doc/rotamer.rst
+++ b/sidechain/doc/rotamer.rst
@@ -254,7 +254,7 @@ Details can be found in the relevant publication [canutescu2003]_.
 The VINA scoring function
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
-The VINA scoring function is a combination of scores that are named 
+The VINA scoring function is a combination of functional forms that are named 
 gaussian1, gaussian2, repulsion, hydrophobic and hbond in the Autodock Vina
 software [trott2010]_. VINA only evaluates heavy atoms. Gaussian1, gaussian2
 and repulsion are evaluated for all pairs of particles. Hydrophobic is only
@@ -263,6 +263,43 @@ evaluated between hydrogen bond donor/acceptor pairs. While SCWRL3 and SCWRL4
 are intended to evaluate sidechain-sidechain interactions in proteins, 
 VINA is mainly targeted at interactions between sidechains and ligands.
 
+The functional forms are linearly combined with the default weights from
+the Autodock Vina Software. They're set as global variables and can be
+extracted with:
+
+.. method:: GetVINAWeightGaussian1
+
+
+.. method:: GetVINAWeightGaussian2
+
+
+.. method:: GetVINAWeightRepulsion
+
+
+.. method:: GetVINAWeightHydrophobic
+
+
+.. method:: GetVINAWeightHBond
+
+You can set custom weights. A call to the following functions overwrites 
+according weights globally which affects any subsequent score evaluation:
+
+.. method:: SetVINAWeightGaussian1(w)
+
+
+.. method:: SetVINAWeightGaussian2(w)
+
+
+.. method:: SetVINAWeightRepulsion(w)
+
+
+.. method:: SetVINAWeightHydrophobic(w)
+
+
+.. method:: SetVINAWeightHBond(w)
+
+
+
 The VINA scoring function differentiates between the following particle types:
 
 .. class:: VINAParticleType
diff --git a/sidechain/pymod/export_rotamer_constructor.cc b/sidechain/pymod/export_rotamer_constructor.cc
index 0950dddb..ae7853a0 100644
--- a/sidechain/pymod/export_rotamer_constructor.cc
+++ b/sidechain/pymod/export_rotamer_constructor.cc
@@ -23,6 +23,7 @@
 #include <promod3/sidechain/scwrl4_rotamer_constructor.hh>
 #include <promod3/sidechain/scwrl3_rotamer_constructor.hh>
 #include <promod3/sidechain/vina_rotamer_constructor.hh>
+#include <promod3/sidechain/vina_particle_scoring.hh>
 
 using namespace boost::python;
 using namespace promod3::sidechain;
@@ -475,4 +476,15 @@ void export_RotamerConstructor(){
          &VINARotamerConstructor::ConstructFRMRotamerHeuristic, (arg("res")))
   ;
 
+  def("SetVINAWeightGaussian1", &SetVINAWeightGaussian1, (arg("weight")));
+  def("SetVINAWeightGaussian2", &SetVINAWeightGaussian2, (arg("weight")));
+  def("SetVINAWeightRepulsion", &SetVINAWeightRepulsion, (arg("weight")));
+  def("SetVINAWeightHydrophobic", &SetVINAWeightHydrophobic, (arg("weight")));
+  def("SetVINAWeightHBond", &SetVINAWeightHBond, (arg("weight")));
+
+  def("GetVINAWeightGaussian1", &GetVINAWeightGaussian1);
+  def("GetVINAWeightGaussian2", &GetVINAWeightGaussian2);
+  def("GetVINAWeightRepulsion", &GetVINAWeightRepulsion);
+  def("GetVINAWeightHydrophobic", &GetVINAWeightHydrophobic);
+  def("GetVINAWeightHBond", &GetVINAWeightHBond);
 }
diff --git a/sidechain/src/vina_particle_scoring.cc b/sidechain/src/vina_particle_scoring.cc
index c048e214..c5a52d6f 100644
--- a/sidechain/src/vina_particle_scoring.cc
+++ b/sidechain/src/vina_particle_scoring.cc
@@ -20,6 +20,25 @@
 
 namespace promod3 { namespace sidechain {
 
+Real VINAParam::gaussian1_weight_ = -0.035579;
+Real VINAParam::gaussian2_weight_ = -0.005156;
+Real VINAParam::repulsion_weight_ = 0.840245;
+Real VINAParam::hydrophobic_weight_ = -0.035069;
+Real VINAParam::hbond_weight_ = -0.587439;
+
+void SetVINAWeightGaussian1(float w) { VINAParam::gaussian1_weight_ = w; }
+void SetVINAWeightGaussian2(float w) { VINAParam::gaussian2_weight_ = w; }
+void SetVINAWeightRepulsion(float w) { VINAParam::repulsion_weight_ = w; }
+void SetVINAWeightHydrophobic(float w) { VINAParam::hydrophobic_weight_ = w; }
+void SetVINAWeightHBond(float w) { VINAParam::hbond_weight_ = w; }
+
+Real GetVINAWeightGaussian1() { return VINAParam::gaussian1_weight_; }
+Real GetVINAWeightGaussian2() { return VINAParam::gaussian2_weight_; }
+Real GetVINAWeightRepulsion() { return VINAParam::repulsion_weight_; }
+Real GetVINAWeightHydrophobic() { return VINAParam::hydrophobic_weight_; }
+Real GetVINAWeightHBond() { return VINAParam::hbond_weight_; }
+
+
 VINAParam::VINAParam(VINAParticleType t, 
                      const geom::Vec3& pos): pos_(pos),
                                              radius_(VINAVDWRadii[t]),
@@ -82,13 +101,6 @@ bool VINAParam::EqualTo(PScoringParam* other) const {
 
 Real VINAPairwiseScore(VINAParam* p_one, VINAParam* p_two) {
 
-  // linear weights of the different terms
-  const Real w_gaussian1 = -0.035579;
-  const Real w_gaussian2 = -0.005156;
-  const Real w_repulsion = 0.840245;
-  const Real w_hydrophobic = -0.035069;
-  const Real w_hbond = -0.587439; 
-
   Real r = geom::Distance(p_one->pos_, p_two->pos_);
   Real d = r - (p_one->radius_ + p_two->radius_);
 
@@ -120,11 +132,11 @@ Real VINAPairwiseScore(VINAParam* p_one, VINAParam* p_two) {
     e_hbond = (d < Real(-0.7) ? 1.0 : (d/Real(-0.7)));
   }
 
-  Real e = w_gaussian1 * e_gaussian1 +
-           w_gaussian2 * e_gaussian2 +
-           w_repulsion * e_repulsion +
-           w_hydrophobic * e_hydrophobic +
-           w_hbond * e_hbond;
+  Real e = VINAParam::gaussian1_weight_ * e_gaussian1 +
+           VINAParam::gaussian2_weight_ * e_gaussian2 +
+           VINAParam::repulsion_weight_ * e_repulsion +
+           VINAParam::hydrophobic_weight_ * e_hydrophobic +
+           VINAParam::hbond_weight_ * e_hbond;
 
   return e;
 }
diff --git a/sidechain/src/vina_particle_scoring.hh b/sidechain/src/vina_particle_scoring.hh
index 204cfb61..b2df5010 100644
--- a/sidechain/src/vina_particle_scoring.hh
+++ b/sidechain/src/vina_particle_scoring.hh
@@ -64,6 +64,18 @@ const Real VINAVDWRadii[18] = {
   0.0 // INVALID
 };
 
+void SetVINAWeightGaussian1(float w);
+void SetVINAWeightGaussian2(float w);
+void SetVINAWeightRepulsion(float w);
+void SetVINAWeightHydrophobic(float w);
+void SetVINAWeightHBond(float w);
+
+Real GetVINAWeightGaussian1();
+Real GetVINAWeightGaussian2();
+Real GetVINAWeightRepulsion();
+Real GetVINAWeightHydrophobic();
+Real GetVINAWeightHBond();
+
 struct VINAParam : public PScoringParam {
 
   VINAParam(VINAParticleType t, const geom::Vec3& pos);
@@ -95,6 +107,11 @@ struct VINAParam : public PScoringParam {
   bool hydrophobic_;
   bool hbond_acceptor_;
   bool hbond_donor_;
+  static Real gaussian1_weight_;
+  static Real gaussian2_weight_;
+  static Real repulsion_weight_;
+  static Real hydrophobic_weight_;
+  static Real hbond_weight_;
 };
 
 Real VINAPairwiseScore(VINAParam* p_one, VINAParam* p_two);
-- 
GitLab