OpenStructure
distance_map.hh
Go to the documentation of this file.
1 //------------------------------------------------------------------------------
2 // This file is part of the OpenStructure project <www.openstructure.org>
3 //
4 // Copyright (C) 2008-2020 by the OpenStructure authors
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License as published by the Free
8 // Software Foundation; either version 3.0 of the License, or (at your option)
9 // any later version.
10 // This library is distributed in the hope that it will be useful, but WITHOUT
11 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
13 // details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with this library; if not, write to the Free Software Foundation, Inc.,
17 // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 //------------------------------------------------------------------------------
19 
20 /*
21  Author: Marco Biasini, Juergen Haas, Gerardo Tauriello
22  */
23 
24 #ifndef OST_SEQ_ALG_DISTANCE_MAP_HH
25 #define OST_SEQ_ALG_DISTANCE_MAP_HH
26 
27 #include <numeric>
28 #include <cmath>
29 #include <algorithm>
30 
31 #include <ost/tri_matrix.hh>
32 #include <ost/integrity_error.hh>
35 
36 namespace ost { namespace seq { namespace alg {
37 
39 // HELPERS
40 namespace {
41 
42 template <typename Pair>
43 struct sum_first {
44  typename Pair::first_type operator()(const typename Pair::first_type& x,
45  const Pair& y) const
46  {
47  return x + y.first;
48  }
49 };
50 
51 template <typename Pair>
52 struct CompareFirst {
53  bool operator()(const Pair& x, const Pair& y) const {
54  return (x.first < y.first);
55  }
56 };
57 
58 } // anon ns
60 
64 
65  typedef std::vector<std::pair<Real,int> > FloatArray;
66 
67  uint GetDataSize() const { return values_.size(); }
68 
69  void Add(Real dist, int index) {
70  values_.push_back(std::make_pair(dist,index));
71  }
72 
73  Real GetAverage() const {
74  if (values_.empty()) throw IntegrityError("List of distances empty!");
75  Real sum = std::accumulate(values_.begin(), values_.end(), 0.0,
76  sum_first<FloatArray::value_type>());
77  return sum/values_.size();
78  }
79 
80  std::pair<Real,int> GetMin() const {
81  if (values_.empty()) throw IntegrityError("List of distances empty!");
82  return *std::min_element(values_.begin(), values_.end(),
83  CompareFirst<FloatArray::value_type>());
84  }
85 
86  std::pair<Real,int> GetMax() const {
87  if (values_.empty()) throw IntegrityError("List of distances empty!");
88  return *std::max_element(values_.begin(), values_.end(),
89  CompareFirst<FloatArray::value_type>());
90  }
91 
92  std::pair<Real,int> GetDataElement(uint index) const {
93  if (values_.size() != 0) {
94  if (values_.size() > index) {
95  return values_[index];
96  } else {
97  std::stringstream error_message;
98  error_message << "size of distances data vector: " << values_.size()
99  << " is smaller than requested: " << index+1;
100  throw IntegrityError(error_message.str());
101  }
102  } else {
103  throw IntegrityError("List of distances empty!");
104  }
105  }
106 
107  Real GetStdDev() const {
108  if (values_.empty()) throw IntegrityError("List of distances empty!");
109  const Real avg = this->GetAverage();
110  return this->do_get_stddev(avg);
111  }
112 
113  Real GetWeightedStdDev(Real sigma) const {
114  if (values_.empty()) throw IntegrityError("List of distances empty!");
115  const Real avg = this->GetAverage();
116  const Real weight = std::exp(avg/(-2*sigma));
117  return this->do_get_stddev(avg) * weight;
118  }
119 
120  Real GetNormStdDev() const {
121  if (values_.empty()) throw IntegrityError("List of distances empty!");
122  const Real avg = this->GetAverage();
123  return this->do_get_stddev(avg) / avg;
124  }
125 
126  bool operator==(const Distances& rhs) const { return values_ == rhs.values_; }
127  bool operator!=(const Distances& rhs) const { return values_ != rhs.values_; }
128 
129 private:
130 
131  Real do_get_stddev(Real avg) const {
132  Real var = 0.0;
133  for (FloatArray::const_iterator i = values_.begin(),
134  e = values_.end(); i != e; ++i) {
135  const Real diff = i->first - avg;
136  var += diff*diff;
137  }
138  return std::sqrt(var/values_.size());
139  }
140 
141  FloatArray values_;
142 };
143 
147 class DistanceMap;
148 typedef boost::shared_ptr<DistanceMap> DistanceMapPtr;
149 
151 public:
152  DistanceMap(uint nresidues, uint num_structures)
153  : dist_(nresidues), num_structures_(num_structures) { }
154 
155  void AddDistance(uint i, uint j, Real distance, int index) {
156  dist_(i ,j).Add(distance,index);
157  }
158 
159  const Distances& GetDistances(uint i, uint j) const {
160  return dist_(i, j);
161  }
162 
163  uint GetSize() const { return dist_.GetSize(); }
164 
165  uint GetNumStructures() const { return num_structures_; }
166 
167 private:
168  TriMatrix<Distances> dist_;
169  uint num_structures_;
170 };
171 
172 
184 
185 }}}
186 
187 #endif
triangular matrix template
Definition: tri_matrix.hh:13
representation of a multiple sequence alignemnt consisting of two or more sequences
const Distances & GetDistances(uint i, uint j) const
void AddDistance(uint i, uint j, Real distance, int index)
DistanceMap(uint nresidues, uint num_structures)
unsigned int uint
Definition: base.hh:29
float Real
Definition: base.hh:44
def alignment(query_db, trg_db, resultDB, directory, resultDB_m8, sen=None, exe_path=None, start_sens=None, sens_steps=None, fmt=None)
Definition: mmseqs2.py:57
DistanceMapPtr DLLEXPORT_OST_SEQ_ALG CreateDistanceMap(const seq::AlignmentHandle &alignment)
create distance map from a multiple sequence alignment.
boost::shared_ptr< DistanceMap > DistanceMapPtr
Definition: base.dox:1
#define DLLEXPORT_OST_SEQ_ALG
Container for a pair wise distance for each structure. Each structure is identified by its index in t...
Definition: distance_map.hh:63
std::pair< Real, int > GetMax() const
Definition: distance_map.hh:86
void Add(Real dist, int index)
Definition: distance_map.hh:69
bool operator==(const Distances &rhs) const
Real GetWeightedStdDev(Real sigma) const
bool operator!=(const Distances &rhs) const
std::pair< Real, int > GetDataElement(uint index) const
Definition: distance_map.hh:92
std::pair< Real, int > GetMin() const
Definition: distance_map.hh:80
std::vector< std::pair< Real, int > > FloatArray
Definition: distance_map.hh:65