IsoSpec  1.95
marginalTrek++.h
1 /*
2  * Copyright (C) 2015-2019 Mateusz Łącki and Michał Startek.
3  *
4  * This file is part of IsoSpec.
5  *
6  * IsoSpec is free software: you can redistribute it and/or modify
7  * it under the terms of the Simplified ("2-clause") BSD licence.
8  *
9  * IsoSpec is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
12  *
13  * You should have received a copy of the Simplified BSD Licence
14  * along with IsoSpec. If not, see <https://opensource.org/licenses/BSD-2-Clause>.
15  */
16 
17 #pragma once
18 
19 #include <tuple>
20 #include <unordered_map>
21 #include <queue>
22 #include "conf.h"
23 #include "allocator.h"
24 #include "operators.h"
25 #include "summator.h"
26 
27 
28 namespace IsoSpec
29 {
30 
31 Conf initialConfigure(int atomCnt, int isotopeNo, const double* probs);
32 
33 
34 void printMarginal(const std::tuple<double*,double*,int*,int>& results, int dim);
35 
37 
44 class Marginal
45 {
46 private:
47  bool disowned;
48 protected:
49  const unsigned int isotopeNo;
50  const unsigned int atomCnt;
51  const double* const atom_lProbs;
52  const double* const atom_masses;
53  const double loggamma_nominator;
54  const Conf mode_conf;
55  const double mode_lprob;
56  const double mode_mass;
57  const double mode_prob;
58  const double smallest_lprob;
61 public:
63 
70  Marginal(
71  const double* _masses, // masses size = logProbs size = isotopeNo
72  const double* _probs,
73  int _isotopeNo,
74  int _atomCnt
75  );
76 
77  // Get rid of the C++ generated copy and assignment constructors.
78  Marginal(Marginal& other) = delete;
79  Marginal& operator= (const Marginal& other) = delete;
80 
82  Marginal(Marginal&& other);
83 
85  virtual ~Marginal();
86 
88 
91  inline int get_isotopeNo() const { return isotopeNo; };
92 
93  inline const double* get_lProbs() const { return atom_lProbs; };
94 
96 
99  double getLightestConfMass() const;
100 
102 
105  double getHeaviestConfMass() const;
106 
108 
112  double getMonoisotopicConfMass() const;
113 
115 
118  inline double getModeLProb() const { return mode_lprob; };
119 
121 
124  inline double getModeMass() const { return mode_mass; };
125 
127 
130  inline double getModeProb() const { return mode_prob; };
131 
133 
136  inline double getSmallestLProb() const { return smallest_lprob; };
137 
139 
142  double getTheoreticalAverageMass() const;
143 
145 
149  inline double logProb(Conf conf) const { return loggamma_nominator + unnormalized_logProb(conf, atom_lProbs, isotopeNo); };
150 };
151 
152 
154 class MarginalTrek : public Marginal
155 {
156 private:
157  int current_count;
158  const KeyHasher keyHasher;
159  const ConfEqual equalizer;
160  const ConfOrderMarginal orderMarginal;
161  std::unordered_map<Conf,int,KeyHasher,ConfEqual> visited;
162  std::priority_queue<Conf,std::vector<Conf>,ConfOrderMarginal> pq;
163  Summator totalProb;
164  Conf candidate;
165  Allocator<int> allocator;
166  std::vector<double> _conf_lprobs;
167  std::vector<double> _conf_masses;
168  std::vector<int*> _confs;
169 
171  bool add_next_conf();
172 
173 public:
175 
179  MarginalTrek(
180  Marginal&& m,
181  int tabSize = 1000,
182  int hashSize = 1000
183  );
184 
186 
192  inline bool probeConfigurationIdx(int idx)
193  {
194  while(current_count <= idx)
195  if(!add_next_conf())
196  return false;
197  return true;
198  }
199 
200 
202 
206  int processUntilCutoff(double cutoff);
207 
208  inline const std::vector<double>& conf_lprobs() const { return _conf_lprobs; };
209  inline const std::vector<double>& conf_masses() const { return _conf_masses; };
210  inline const std::vector<int*>& confs() const { return _confs; };
211 
212 
213  virtual ~MarginalTrek();
214 };
215 
216 
218 
227 {
228 protected:
229  std::vector<Conf> configurations;
230  Conf* confs;
231  unsigned int no_confs;
232  double* masses;
233  double* lProbs;
234  double* probs;
235  Allocator<int> allocator;
236 public:
238 
246  Marginal&& m,
247  double lCutOff,
248  bool sort = true,
249  int tabSize = 1000,
250  int hashSize = 1000
251  );
252 
254  virtual ~PrecalculatedMarginal();
255 
257 
260  inline bool inRange(unsigned int idx) const { return idx < no_confs; };
261 
263 
267  inline const double& get_lProb(int idx) const { return lProbs[idx]; };
268 
270 
274  inline const double& get_prob(int idx) const { return probs[idx]; };
275 
277 
281  inline const double& get_mass(int idx) const { return masses[idx]; };
282 
284 
287  inline const double* get_lProbs_ptr() const { return lProbs; };
288 
290 
293  inline const double* get_masses_ptr() const { return masses; };
294 
295 
297 
301  inline const Conf& get_conf(int idx) const { return confs[idx]; };
302 
304 
307  inline unsigned int get_no_confs() const { return no_confs; };
308 };
309 
310 
312 
315 class LayeredMarginal : public Marginal
316 {
317 private:
318  double current_threshold;
319  std::vector<Conf> configurations;
320  std::vector<Conf> fringe;
321  Allocator<int> allocator;
322  unsigned int sorted_up_to_idx;
323  const ConfEqual equalizer;
324  const KeyHasher keyHasher;
325  const ConfOrderMarginalDescending orderMarginal;
326  std::vector<double> lProbs;
327  std::vector<double> probs;
328  std::vector<double> masses;
329  double* guarded_lProbs;
330  const int hashSize;
331 
332 public:
334 
338  LayeredMarginal(Marginal&& m, int tabSize = 1000, int hashSize = 1000);
339 
341 
345  bool extend(double new_threshold);
346 
348  inline double get_lProb(int idx) const { return guarded_lProbs[idx]; }; // access to idx == -1 is valid and gives a guardian of +inf
349 
351  inline double get_prob(int idx) const { return probs[idx]; };
352 
354  inline double get_mass(int idx) const { return masses[idx]; };
355 
357  inline const double* get_lProbs_ptr() const { return lProbs.data()+1; };
358 
360  inline const Conf& get_conf(int idx) const { return configurations[idx]; };
361 
363  inline unsigned int get_no_confs() const { return configurations.size(); };
364 };
365 
366 
367 } // namespace IsoSpec
368 
IsoSpec::Marginal::getMonoisotopicConfMass
double getMonoisotopicConfMass() const
Get the mass of the monoisotopic subisotopologue.
Definition: marginalTrek++.cpp:262
IsoSpec::PrecalculatedMarginal::PrecalculatedMarginal
PrecalculatedMarginal(Marginal &&m, double lCutOff, bool sort=true, int tabSize=1000, int hashSize=1000)
The move constructor (disowns the Marginal).
Definition: marginalTrek++.cpp:390
IsoSpec::Marginal::Marginal
Marginal(const double *_masses, const double *_probs, int _isotopeNo, int _atomCnt)
Class constructor.
Definition: marginalTrek++.cpp:197
IsoSpec::PrecalculatedMarginal::get_masses_ptr
const double * get_masses_ptr() const
Get the table of the masses of subisotopologues.
Definition: marginalTrek++.h:293
IsoSpec::PrecalculatedMarginal
Precalculated Marginal class.
Definition: marginalTrek++.h:226
IsoSpec::Marginal
The marginal distribution class (a subisotopologue).
Definition: marginalTrek++.h:44
IsoSpec::Marginal::mode_lprob
const double mode_lprob
Definition: marginalTrek++.h:55
IsoSpec::LayeredMarginal::get_no_confs
unsigned int get_no_confs() const
Get the number of precomputed subisotopologues, see details in PrecalculatedMarginal::get_no_confs.
Definition: marginalTrek++.h:363
IsoSpec::LayeredMarginal::extend
bool extend(double new_threshold)
Extend the set of computed subisotopologues to those above the new threshold.
Definition: marginalTrek++.cpp:492
IsoSpec::PrecalculatedMarginal::get_lProbs_ptr
const double * get_lProbs_ptr() const
Get the table of the log-probabilities of subisotopologues.
Definition: marginalTrek++.h:287
IsoSpec::PrecalculatedMarginal::get_mass
const double & get_mass(int idx) const
Get the mass of the idx-th subisotopologue.
Definition: marginalTrek++.h:281
IsoSpec::Marginal::logProb
double logProb(Conf conf) const
Calculate the log-probability of a given subisotopologue.
Definition: marginalTrek++.h:149
IsoSpec::LayeredMarginal::LayeredMarginal
LayeredMarginal(Marginal &&m, int tabSize=1000, int hashSize=1000)
Move constructor: specializes the Marginal class.
Definition: marginalTrek++.cpp:482
IsoSpec::Allocator< int >
IsoSpec::MarginalTrek::processUntilCutoff
int processUntilCutoff(double cutoff)
Calculate subisotopologues with probability above or equal to the cut-off.
Definition: marginalTrek++.cpp:361
IsoSpec::Marginal::getHeaviestConfMass
double getHeaviestConfMass() const
Get the mass of the heaviest subisotopologue.
Definition: marginalTrek++.cpp:253
IsoSpec::Marginal::smallest_lprob
const double smallest_lprob
Definition: marginalTrek++.h:58
IsoSpec::Marginal::getSmallestLProb
double getSmallestLProb() const
The the log-probability of the lightest subisotopologue.
Definition: marginalTrek++.h:136
IsoSpec::initialConfigure
Conf initialConfigure(const int atomCnt, const int isotopeNo, const double *probs, const double *lprobs)
Find one of the most probable subisotopologues.
Definition: marginalTrek++.cpp:55
IsoSpec
Definition: allocator.cpp:21
IsoSpec::Marginal::getModeLProb
double getModeLProb() const
Get the log-probability of the mode subisotopologue.
Definition: marginalTrek++.h:118
IsoSpec::PrecalculatedMarginal::get_conf
const Conf & get_conf(int idx) const
Get the counts of isotopes that define the subisotopologue.
Definition: marginalTrek++.h:301
IsoSpec::PrecalculatedMarginal::get_lProb
const double & get_lProb(int idx) const
Get the log-probability of the idx-th subisotopologue.
Definition: marginalTrek++.h:267
IsoSpec::PrecalculatedMarginal::~PrecalculatedMarginal
virtual ~PrecalculatedMarginal()
Destructor.
Definition: marginalTrek++.cpp:466
IsoSpec::LayeredMarginal::get_mass
double get_mass(int idx) const
get the mass of the idx-th subisotopologue, see details in PrecalculatedMarginal::get_mass.
Definition: marginalTrek++.h:354
IsoSpec::Marginal::mode_prob
const double mode_prob
Definition: marginalTrek++.h:57
IsoSpec::Marginal::getLightestConfMass
double getLightestConfMass() const
Get the mass of the lightest subisotopologue.
Definition: marginalTrek++.cpp:244
IsoSpec::ConfOrderMarginal
Definition: operators.h:78
IsoSpec::KeyHasher
Definition: operators.h:27
IsoSpec::Marginal::atom_masses
const double *const atom_masses
Definition: marginalTrek++.h:52
IsoSpec::Marginal::mode_mass
const double mode_mass
Definition: marginalTrek++.h:56
IsoSpec::Summator
Definition: summator.h:75
IsoSpec::Marginal::atomCnt
const unsigned int atomCnt
Definition: marginalTrek++.h:50
IsoSpec::ConfEqual
Definition: operators.h:45
IsoSpec::LayeredMarginal::get_conf
const Conf & get_conf(int idx) const
get the counts of isotopes that define the subisotopologue, see details in PrecalculatedMarginal::get...
Definition: marginalTrek++.h:360
IsoSpec::PrecalculatedMarginal::get_no_confs
unsigned int get_no_confs() const
Get the number of precomputed subisotopologues.
Definition: marginalTrek++.h:307
IsoSpec::Marginal::getModeMass
double getModeMass() const
The the mass of the mode subisotopologue.
Definition: marginalTrek++.h:124
IsoSpec::Marginal::getTheoreticalAverageMass
double getTheoreticalAverageMass() const
The theoretical average mass of the molecule.
Definition: marginalTrek++.cpp:275
IsoSpec::LayeredMarginal
LayeredMarginal class.
Definition: marginalTrek++.h:315
IsoSpec::LayeredMarginal::get_lProbs_ptr
const double * get_lProbs_ptr() const
get the pointer to lProbs array. Accessing index -1 is legal and returns a guardian of -inf....
Definition: marginalTrek++.h:357
IsoSpec::MarginalTrek::probeConfigurationIdx
bool probeConfigurationIdx(int idx)
Check if the table of computed subisotopologues does not have to be extended.
Definition: marginalTrek++.h:192
IsoSpec::Marginal::get_isotopeNo
int get_isotopeNo() const
Get the number of isotopes of the investigated element.
Definition: marginalTrek++.h:91
IsoSpec::Marginal::atom_lProbs
const double *const atom_lProbs
Definition: marginalTrek++.h:51
IsoSpec::Marginal::~Marginal
virtual ~Marginal()
Destructor.
Definition: marginalTrek++.cpp:233
IsoSpec::LayeredMarginal::get_prob
double get_prob(int idx) const
get the probability of the idx-th subisotopologue, see details in PrecalculatedMarginal::get_eProb.
Definition: marginalTrek++.h:351
IsoSpec::Marginal::getModeProb
double getModeProb() const
The the probability of the mode subisotopologue.
Definition: marginalTrek++.h:130
IsoSpec::MarginalTrek::MarginalTrek
MarginalTrek(Marginal &&m, int tabSize=1000, int hashSize=1000)
Move constructor: specializes the Marginal class.
Definition: marginalTrek++.cpp:284
IsoSpec::PrecalculatedMarginal::inRange
bool inRange(unsigned int idx) const
Is there a subisotopologue with a given number?
Definition: marginalTrek++.h:260
IsoSpec::MarginalTrek
The marginal distribution class (a subisotopologue).
Definition: marginalTrek++.h:154
IsoSpec::LayeredMarginal::get_lProb
double get_lProb(int idx) const
get the log-probability of the idx-th subisotopologue, see details in PrecalculatedMarginal::get_lPro...
Definition: marginalTrek++.h:348
IsoSpec::PrecalculatedMarginal::get_prob
const double & get_prob(int idx) const
Get the probability of the idx-th subisotopologue.
Definition: marginalTrek++.h:274
IsoSpec::Marginal::isotopeNo
const unsigned int isotopeNo
Definition: marginalTrek++.h:49
IsoSpec::Marginal::loggamma_nominator
const double loggamma_nominator
Definition: marginalTrek++.h:53
IsoSpec::Marginal::mode_conf
const Conf mode_conf
Definition: marginalTrek++.h:54
IsoSpec::ConfOrderMarginalDescending
Definition: operators.h:92