Added score container to keep track of different loop candidate scores.

modelling/pymod/CMakeLists.txt

@@ -7,6 +7,7 @@ set(MODELLING_CPP
   export_model.cc
   export_monte_carlo.cc
   export_rigid_blocks.cc
+  export_score_container.cc
   export_scoring_weights.cc
   export_el_enumerador.cc
   wrap_modelling.cc

modelling/pymod/_closegaps.py

@@ -37,7 +37,7 @@ def _CombineScores(all_scores, weights):
     return combined_scores
 
 def _GetBestLC(mhandle, loop_candidates, start_resnum, chain_idx,
-               db_scores={}, max_num_all_atom=0):
+               db_scores, max_num_all_atom=0):
     # returns (min_idx, min_score)
 
     # dummy check
@@ -45,45 +45,44 @@ def _GetBestLC(mhandle, loop_candidates, start_resnum, chain_idx,
         raise RuntimeError("Need at least 1 candidate to choose from!")
     
     # setup weights
-    all_atom = (max_num_all_atom > 0)
-    weights = ScoringWeights.GetWeights(bool(db_scores), all_atom)
+    with_db = (not db_scores.IsEmpty())
+    with_aa = (max_num_all_atom > 0)
+    weights = ScoringWeights.GetWeights(with_db, with_aa)
 
     # setup all scores object (we collect needed scores there)
-    all_scores = dict(db_scores) # copy to avoid changing db_scores
+    all_scores = db_scores.Copy()  # copy to avoid changing db_scores
 
     # add backbone scores
     for key in ScoringWeights.GetBackboneScoringKeys():
         scores = loop_candidates.CalculateScores(mhandle.backbone_scorer,
                                                  key, start_resnum, chain_idx)
-        all_scores[key] = scores
+        all_scores.Set(key, scores)
 
     # add all atom scores?
-    if all_atom:
+    if with_aa:
         # get stuff from mhandle
         aa_scorer = mhandle.all_atom_scorer
         aa_score_env = mhandle.all_atom_scorer_env
         aa_sc_env = mhandle.all_atom_sidechain_env
         sc_rec = mhandle.sidechain_reconstructor
         # get best ones based on previous scores
-        non_aa_weights = ScoringWeights.GetWeights(bool(db_scores), False)
-        non_aa_scores = _CombineScores(all_scores, non_aa_weights)
-        arg_sorted_scores = []
-        for i,v in enumerate(non_aa_scores):
-            arg_sorted_scores.append((v,i))
-        arg_sorted_scores.sort()
+        non_aa_weights = ScoringWeights.GetWeights(with_db, False)
+        non_aa_scores = all_scores.LinearCombine(non_aa_weights)
+        arg_sorted_scores = sorted([(v,i) for i,v in enumerate(non_aa_scores)])
         min_indices = [i for v,i in arg_sorted_scores[:max_num_all_atom]]
         # extract relevant loop candidates and scores
         aa_loop_candidates = loop_candidates.Extract(min_indices)
-        for key in all_scores:
-            all_scores[key] = [all_scores[key][i] for i in min_indices]
+        all_scores = all_scores.Extract(min_indices)
+        
+        # TODO: use LoopCandidates feature to do all of the below...
+
         # keep backup of loop of interest
         loop_length = len(loop_candidates[0])
         loop_backup = aa_sc_env.GetEnvironment(start_resnum, loop_length,
                                                chain_idx)
         # setup empty lists for AA results
         aa_scoring_keys = ScoringWeights.GetAllAtomScoringKeys()
-        for key in aa_scoring_keys:
-            all_scores[key] = list()
+        aa_scores = {key: list() for key in aa_scoring_keys}
         # score them w/o relaxation
         for lc_bb_list in aa_loop_candidates:
             # setup sc env. and reconstruct
@@ -94,19 +93,21 @@ def _GetBestLC(mhandle, loop_candidates, start_resnum, chain_idx,
             for key in aa_scoring_keys:
                 score = aa_scorer[key].CalculateScore(start_resnum,
                                                       loop_length, chain_idx)
-                all_scores[key].append(score)
+                aa_scores[key].append(score)
         # reset aa_sc_env
         aa_sc_env.SetEnvironment(loop_backup)
         # NOTE: we leave the aa_score_env in a modified state
         # -> this is ok for aa_score_env as we always set proper surroundings
         # -> this is not ok for SC (due to frame approach), so we back it up
+        for key in aa_scoring_keys:
+            all_scores.Set(key, aa_scores[key])
 
     # get / return best
-    scores = _CombineScores(all_scores, weights)
+    scores = all_scores.LinearCombine(weights)
     min_score = min(scores)
     min_idx = scores.index(min_score)
     # translate to original indexing for all atom
-    if all_atom:
+    if with_aa:
         min_idx = min_indices[min_idx]
     return (min_idx, min_score)
 
@@ -198,14 +199,9 @@ def _CloseLoopFrame(mhandle, gap_orig, actual_candidates, actual_extended_gaps,
             back_mapper.append((i, j))
 
     # db scores requested
-    db_scores = {}
-    if actual_db_scores:
-        # build up lists matching loop candidates
-        for key in _GetDbScoringKeys():
-            scores = []
-            for actual_scores in actual_db_scores:
-                scores.extend(actual_scores[key])
-            db_scores[key] = scores
+    db_scores = ScoreContainer()
+    for actual_scores in actual_db_scores:
+        db_scores.Extend(actual_scores)
     
     if len(final_loop_candidates) > 0:
         # get best
@@ -249,7 +245,7 @@ def _CloseLoopBare(mhandle, gap_orig, actual_candidates, actual_extended_gaps,
         n_candidates += len(loop_candidates)
         # get current best
         before_resnum = actual_extended_gaps[i].before.GetNumber().GetNum()
-        db_scores = {}
+        db_scores = ScoreContainer()
         if actual_db_scores:
             db_scores = actual_db_scores[i]
         best_idx, score = _GetBestLC(mhandle, loop_candidates, before_resnum,
@@ -771,7 +767,9 @@ def FillLoopsByDatabase(mhandle, fragment_db, structure_db,
                 continue
             # check for stem rmsd before ApplyCCD
             if add_db_features:
-                stem_rmsds = candidates.CalculateStemRMSDs(n_stem, c_stem)
+                db_scores = ScoreContainer()
+                db_scores.Set(ScoringWeights.GetStemRMSDsKey(),
+                              candidates.CalculateStemRMSDs(n_stem, c_stem))
             # try to close loops
             try:
                 #pylint: disable=broad-except
@@ -800,22 +798,19 @@ def FillLoopsByDatabase(mhandle, fragment_db, structure_db,
             if add_db_features:
                 # setup
                 assert(len(orig_indices) == len(candidates))
-                stem_rmsd_key = ScoringWeights.GetStemRMSDsKey()
                 seq_prof_key = ScoringWeights.GetSequenceProfileScoresKey()
                 str_prof_key = ScoringWeights.GetStructureProfileScoresKey()
                 # get subset of stem rmsd scores
-                db_scores = {stem_rmsd_key: []}
-                for i in orig_indices:
-                    db_scores[stem_rmsd_key].append(stem_rmsds[i])
+                db_scores = db_scores.Extract(orig_indices)
                 # add profile scores
                 prof = mhandle.profiles[actual_gap.GetChainIndex()]
                 start_pos = n_stem.GetNumber().GetNum() - 1
                 scores = candidates.CalculateSequenceProfileScores(\
                                structure_db, prof, start_pos)
-                db_scores[seq_prof_key] = scores
+                db_scores.Set(seq_prof_key, scores)
                 scores = candidates.CalculateStructureProfileScores(\
                                structure_db, prof, start_pos)
-                db_scores[str_prof_key] = scores
+                db_scores.Set(str_prof_key, scores)
                 # update list
                 actual_db_scores.append(db_scores)
 

modelling/pymod/export_score_container.cc

+#include <boost/python.hpp>
+#include <promod3/core/message.hh>
+#include <promod3/core/export_helper.hh>
+
+#include <promod3/modelling/score_container.hh>
+
+using namespace promod3;
+using namespace promod3::modelling;
+using namespace boost::python;
+
+// WRAPPERS
+namespace {
+
+boost::python::list WrapGet(const ScoreContainer& sc, const String& key) {
+  boost::python::list result;
+  core::AppendVectorToList(sc.Get(key), result);
+  return result; 
+}
+void WrapSet(ScoreContainer& sc, const String& key,
+             const boost::python::list& scores) {
+  ScoreContainer::ScoreVector v;
+  core::ConvertListToVector(scores, v);
+  sc.Set(key, v);
+}
+
+boost::python::list WrapLinearCombine(const ScoreContainer& sc,
+                                      const boost::python::dict& weights) {
+  // get scores
+  std::map<String, Real> w;
+  core::ConvertDictToMap(weights, w);
+  ScoreContainer::ScoreVector scores = sc.LinearCombine(w);
+  // convert and return
+  boost::python::list return_list;
+  core::AppendVectorToList(scores, return_list);
+  return return_list;
+}
+
+ScoreContainerPtr WrapExtract(const ScoreContainer& sc,
+                              const boost::python::list indices) {
+  std::vector<uint> v;
+  core::ConvertListToVector(indices, v);
+  return sc.Extract(v);
+}
+
+} // anon wrapper ns
+
+void export_score_container() {
+
+  class_<ScoreContainer, ScoreContainerPtr>("ScoreContainer", init<>())
+    .def("IsEmpty", &ScoreContainer::IsEmpty)
+    .def("Contains", &ScoreContainer::Contains, (arg("key")))
+    .def("Get", WrapGet, (arg("key")))
+    .def("Set", WrapSet, (arg("key"), arg("scores")))
+    .def("GetNumCandidates", &ScoreContainer::GetNumCandidates)
+    .def("LinearCombine", WrapLinearCombine, (arg("linear_weights")))
+    .def("Copy", &ScoreContainer::Copy)
+    .def("Extract", WrapExtract, (arg("indices")))
+    .def("Extend", &ScoreContainer::Extend, (arg("other")))
+  ;
+
+}

modelling/pymod/wrap_modelling.cc

@@ -8,6 +8,7 @@ void export_loop_candidate();
 void export_loop_closure();
 void export_monte_carlo();
 void export_rigid_blocks();
+void export_score_container();
 void export_scoring_weights();
 void export_el_enumerador();
 
@@ -21,6 +22,7 @@ BOOST_PYTHON_MODULE(_modelling)
   export_loop_closure();
   export_monte_carlo();
   export_rigid_blocks();
+  export_score_container();
   export_scoring_weights();
   export_el_enumerador();
 }

modelling/src/CMakeLists.txt

@@ -14,6 +14,7 @@ set(MODELLING_SOURCES
   monte_carlo_sampler.cc
   monte_carlo_scorer.cc
   rigid_blocks.cc
+  score_container.cc
   scoring_weights.cc
   el_enumerador.cc
 )
@@ -34,6 +35,7 @@ set(MODELLING_HEADERS
   monte_carlo_sampler.hh
   monte_carlo_scorer.hh
   rigid_blocks.hh
+  score_container.hh
   scoring_weights.hh
   el_enumerador.hh
 )

modelling/src/score_container.cc

+#include <promod3/modelling/score_container.hh>
+
+namespace promod3 { namespace modelling {
+
+ScoreContainer::ScoreVector ScoreContainer::LinearCombine(
+                    const std::map<String, Real>& linear_weights) const {
+  ScoreVector result(num_candidates_, 0);
+  for (std::map<String, Real>::const_iterator i = linear_weights.begin();
+       i != linear_weights.end(); ++i) {
+    const ScoreVector& scores = Get(i->first);
+    const Real weight = i->second;
+    for (uint j = 0; j < num_candidates_; ++j) result[j] += weight * scores[j];
+  }
+  return result;
+}
+
+ScoreContainerPtr
+ScoreContainer::Extract(const std::vector<uint>& indices) const {
+  // check indices
+  for (uint j = 0; j < indices.size(); ++j) {
+    if (indices[j] >= num_candidates_) {
+      throw promod3::Error("Index out-of-bounds in ScoreContainer::Extract!");
+    }
+  }
+  // get results
+  ScoreContainerPtr result(new ScoreContainer);
+  for (std::map<String, ScoreVector>::const_iterator i = map_.begin();
+       i != map_.end(); ++i) {
+    // setup scores directly in map of result
+    ScoreVector& target_scores = result->map_[i->first];
+    target_scores.resize(indices.size());
+    for (uint j = 0; j < indices.size(); ++j) {
+      target_scores[j] = i->second[indices[j]];
+    }
+  }
+  result->num_candidates_ = indices.size();
+  return result;
+}
+
+void ScoreContainer::Extend(const ScoreContainer& other) {
+  // check consistency
+  if (!IsEmpty()) {
+    if (map_.size() != other.map_.size()) {
+      throw promod3::Error("Inconsistent keys in ScoreContainer::Extend!");
+    }
+    for (std::map<String, ScoreVector>::const_iterator i = map_.begin();
+         i != map_.end(); ++i) {
+      if (!other.Contains(i->first)) {
+        throw promod3::Error("Inconsistent keys in ScoreContainer::Extend!");
+      }
+    }
+  }
+  // extend this
+  for (std::map<String, ScoreVector>::const_iterator i = other.map_.begin();
+       i != other.map_.end(); ++i) {
+    ScoreVector& scores = map_[i->first];
+    const ScoreVector& scores_to_add = i->second;
+    scores.insert(scores.end(), scores_to_add.begin(), scores_to_add.end());
+  }
+  num_candidates_ += other.num_candidates_;
+}
+
+}} // ns

modelling/src/score_container.hh

+#ifndef PROMOD3_MODELLING_SCORE_CONTAINER_HH
+#define PROMOD3_MODELLING_SCORE_CONTAINER_HH
+
+#include <boost/shared_ptr.hpp>
+#include <map>
+#include <vector>
+#include <ost/base.hh>
+#include <promod3/core/message.hh>
+
+namespace promod3 { namespace modelling {
+
+class ScoreContainer;
+typedef boost::shared_ptr<ScoreContainer> ScoreContainerPtr;
+
+// convenience class to keep multiple scores for loop candidates
+class ScoreContainer {
+public:
+  // score vector type
+  typedef std::vector<Real> ScoreVector;
+
+  // construct empty container
+  ScoreContainer(): num_candidates_(0) { }
+
+  // check if there's anything
+  bool IsEmpty() { return map_.empty(); }
+
+  // check existance of scores
+  bool Contains(const String& key) const { return map_.find(key) != map_.end(); }
+
+  // r/o get access (throws exc. if nothing found)
+  const ScoreVector& Get(const String& key) const {
+    std::map<String, ScoreVector>::const_iterator found = map_.find(key);
+    if (found == map_.end()) {
+      throw promod3::Error("ScoreContainer does not contain scores for key '"
+                           + key + "'!");
+    } else {
+      return found->second;
+    }
+  }
+
+  // set new scores (consistency enforced)
+  void Set(const String& key, const ScoreVector& scores) {
+    if (!IsEmpty() && num_candidates_ != scores.size()) {
+      throw promod3::Error("Inconsistent scores size in ScoreContainer::Set!");
+    } else {
+      num_candidates_ = scores.size();
+      map_[key] = scores;
+    }
+  }
+
+  // get number of candidates that are scored
+  uint GetNumCandidates() const { return num_candidates_; }
+
+  // calculate linear combination of scores
+  ScoreVector LinearCombine(const std::map<String, Real>& linear_weights) const;
+
+  // extract copy
+  ScoreContainerPtr Copy() const {
+    return ScoreContainerPtr(new ScoreContainer(*this));
+  }
+  // extract copy with subset of data
+  ScoreContainerPtr Extract(const std::vector<uint>& indices) const;
+
+  // extend scores
+  void Extend(const ScoreContainer& other);
+
+private:
+  // DATA
+  std::map<String, ScoreVector> map_;
+  uint num_candidates_;
+};
+
+}} // ns
+
+#endif