Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages
FeatureFinderAlgorithmSimple.h
Go to the documentation of this file.
1 // --------------------------------------------------------------------------
2 // OpenMS -- Open-Source Mass Spectrometry
3 // --------------------------------------------------------------------------
4 // Copyright The OpenMS Team -- Eberhard Karls University Tuebingen,
5 // ETH Zurich, and Freie Universitaet Berlin 2002-2013.
6 //
7 // This software is released under a three-clause BSD license:
8 // * Redistributions of source code must retain the above copyright
9 // notice, this list of conditions and the following disclaimer.
10 // * Redistributions in binary form must reproduce the above copyright
11 // notice, this list of conditions and the following disclaimer in the
12 // documentation and/or other materials provided with the distribution.
13 // * Neither the name of any author or any participating institution
14 // may be used to endorse or promote products derived from this software
15 // without specific prior written permission.
16 // For a full list of authors, refer to the file AUTHORS.
17 // --------------------------------------------------------------------------
18 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
19 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
21 // ARE DISCLAIMED. IN NO EVENT SHALL ANY OF THE AUTHORS OR THE CONTRIBUTING
22 // INSTITUTIONS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
23 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
24 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
25 // OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
26 // WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
27 // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
28 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 //
30 // --------------------------------------------------------------------------
31 // $Maintainer: Clemens Groepl $
32 // $Authors: $
33 // --------------------------------------------------------------------------
34 
35 #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMSIMPLE_H
36 #define OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMSIMPLE_H
37 
39 
43 
44 namespace OpenMS
45 {
61  template <class PeakType, class FeatureType>
63  public FeatureFinderAlgorithm<PeakType, FeatureType>,
64  public FeatureFinderDefs
65  {
66 
67 public:
70  FeatureFinderAlgorithm<PeakType, FeatureType>()
71  {
73  this->check_defaults_ = false;
74  }
75 
76  virtual Param getDefaultParameters() const
77  {
78  Param tmp;
79 
80  SimpleSeeder<PeakType, FeatureType> seeder(this->map_, this->features_, this->ff_);
81  tmp.insert("seeder:", seeder.getParameters());
82  tmp.setSectionDescription("seeder", "Settings for the seeder (Determines potential feature regions)");
83 
84  SimpleExtender<PeakType, FeatureType> extender(this->map_, this->features_, this->ff_);
85  tmp.insert("extender:", extender.getParameters());
86  tmp.setSectionDescription("extender", "Settings for the extender (Collects all peaks belonging to a feature)");
87 
88  ModelFitter<PeakType, FeatureType> fitter(this->map_, this->features_, this->ff_);
89  tmp.insert("fitter:", fitter.getParameters());
90  tmp.setSectionDescription("fitter", "Settings for the modefitter (Fits a model to the data determinging the probapility that they represent a feature.)");
91 
92  return tmp;
93  }
94 
95  virtual void run()
96  {
97 #ifdef DEBUG_FEATUREFINDER
98  UInt seed_nr = 0;
99 #endif
100  SimpleSeeder<PeakType, FeatureType> seeder(this->map_, this->features_, this->ff_);
101  seeder.setParameters(this->getParameters().copy("seeder:", true));
102 
103  SimpleExtender<PeakType, FeatureType> extender(this->map_, this->features_, this->ff_);
104  extender.setParameters(this->getParameters().copy("extender:", true));
105 
106  ModelFitter<PeakType, FeatureType> fitter(this->map_, this->features_, this->ff_);
107  Param params;
108  params.setDefaults(this->getParameters().copy("fitter:", true));
109  params.setValue("fit_algorithm", "simple");
110  fitter.setParameters(params);
111 
113  Summary summary;
114 
115  try
116  {
117  for (;; )
118  {
119 #ifdef DEBUG_FEATUREFINDER
120  std::cout << "===============================" << std::endl;
121  std::cout << "### Seeder (seed # " << ++seed_nr << ")..." << std::endl;
122 #endif
123  IndexPair seed = seeder.nextSeed();
124 
125 #ifdef DEBUG_FEATUREFINDER
126  std::cout << "seed ... " << seed.first << " - " << seed.second << std::endl;
127  std::cout << "### Extender..." << std::endl;
128 #endif
129  ChargedIndexSet index_set;
130  index_set.insert(seed);
131  ChargedIndexSet region;
132  extender.extend(index_set, region);
133 
134 #ifdef DEBUG_FEATUREFINDER
135  std::cout << "### ModelFitter..." << std::endl;
136 #endif
137  try
138  {
139  this->features_->push_back(fitter.fit(region));
140 
141  // gather information for fitting summary
142  {
143  const Feature & f = this->features_->back();
144 
145  // quality, correlation
146  DoubleReal corr = f.getOverallQuality();
147  summary.corr_mean += corr;
148  if (corr < summary.corr_min) summary.corr_min = corr;
149  if (corr > summary.corr_max) summary.corr_max = corr;
150 
151  // charge
152  UInt ch = f.getCharge();
153  if (ch >= summary.charge.size())
154  {
155  summary.charge.resize(ch + 1);
156  }
157  summary.charge[ch]++;
158 
159  // MZ model type
160  const Param & p = f.getModelDescription().getParam();
161  ++summary.mz_model[p.getValue("MZ")];
162 
163  // standard deviation of isotopic peaks
164  if (p.exists("MZ:isotope:stdev") && p.getValue("MZ:isotope:stdev") != DataValue::EMPTY)
165  {
166  ++summary.mz_stdev[p.getValue("MZ:isotope:stdev")];
167  }
168  }
169  }
170  catch (Exception::UnableToFit ex)
171  {
172  //std::cout << "UnableToFit: " << ex.what() << std::endl;
173 
174  // set unused flag for all data points
175  for (IndexSet::const_iterator it = region.begin(); it != region.end(); ++it)
176  {
177  this->ff_->getPeakFlag(*it) = UNUSED;
178  }
179 
180  // gather information for fitting summary
181  {
182  ++summary.no_exceptions;
183  ++summary.exception[ex.getName()];
184  }
185  }
186  } // for
187  } // try
188  catch (NoSuccessor ex)
189  {
190  }
191 
192  this->ff_->endProgress();
193 
194  // print fitting summary
195  {
196  Size size = this->features_->size();
197  std::cout << size << " features were found. " << std::endl;
198 
199  // compute corr_mean
200  summary.corr_mean /= size;
201 
202  std::cout << "FeatureFinder summary:\n"
203  << "Correlation:\n\tminimum: " << summary.corr_min << "\n\tmean: " << summary.corr_mean
204  << "\n\tmaximum: " << summary.corr_max << std::endl;
205 
206  std::cout << "Exceptions:\n";
207  for (std::map<String, UInt>::const_iterator it = summary.exception.begin(); it != summary.exception.end(); ++it)
208  {
209  std::cout << "\t" << it->first << ": " << it->second * 100 / summary.no_exceptions << "% (" << it->second << ")\n";
210  }
211 
212  std::cout << "Chosen mz models:\n";
213  for (std::map<String, UInt>::const_iterator it = summary.mz_model.begin(); it != summary.mz_model.end(); ++it)
214  {
215  std::cout << "\t" << it->first << ": " << it->second * 100 / size << "% (" << it->second << ")\n";
216  }
217 
218  std::cout << "Chosen mz stdevs:\n";
219  for (std::map<float, UInt>::const_iterator it = summary.mz_stdev.begin(); it != summary.mz_stdev.end(); ++it)
220  {
221  std::cout << "\t" << it->first << ": " << it->second * 100 / (size - summary.charge[0]) << "% (" << it->second << ")\n";
222  }
223 
224  std::cout << "Charges:\n";
225  for (Size i = 1; i < summary.charge.size(); ++i)
226  {
227  if (summary.charge[i] != 0)
228  {
229  std::cout << "\t+" << i << ": " << summary.charge[i] * 100 / (size - summary.charge[0]) << "% (" << summary.charge[i] << ")\n";
230  }
231  }
232  }
233  } // run
234 
236  {
237  return new FeatureFinderAlgorithmSimple();
238  }
239 
240  static const String getProductName()
241  {
242  return "simple";
243  }
244 
245 private:
250 
251  }; // FeatureFinderAlgorithmSimple
252 
253 } // namespace OpenMS
254 
255 #endif // OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMSIMPLE_H
UInt no_exceptions
Definition: FeatureFinderAlgorithm.h:53
const ChargeType & getCharge() const
Non-mutable access to charge state.
void extend(const ChargedIndexSet &seed_region, ChargedIndexSet &result_region)
return next seed
Definition: SimpleExtender.h:127
virtual Param getDefaultParameters() const
Returns the default parameters. Reimplment.
Definition: FeatureFinderAlgorithmSimple.h:76
IsotopeCluster::IndexPair IndexPair
Index to peak consisting of two UInts (scan index / peak index)
Definition: FeatureFinderDefs.h:54
FeatureMapType * features_
Output data pointer.
Definition: FeatureFinderAlgorithm.h:144
Param defaults_
Container for default parameters. This member should be filled in the constructor of derived classes!...
Definition: DefaultParamHandler.h:155
void setValue(const String &key, const DataValue &value, const String &description="", const StringList &tags=StringList())
Sets a value.
A more convenient string class.
Definition: String.h:56
const ModelDescription< 2 > & getModelDescription() const
Non-mutable access to the model description.
void insert(const String &prefix, const Param &param)
FeatureFinder * ff_
Pointer to the calling FeatureFinder that is used to access the feature flags.
Definition: FeatureFinderAlgorithm.h:147
A 2-dimensional raw data point or peak.
Definition: Peak2D.h:55
QualityType getOverallQuality() const
Non-mutable access to the overall quality.
Definition: FeatureFinderDefs.h:63
DoubleReal corr_mean
Definition: FeatureFinderAlgorithm.h:57
const MapType * map_
Input data pointer.
Definition: FeatureFinderAlgorithm.h:141
std::map< String, UInt > mz_model
Definition: FeatureFinderAlgorithm.h:54
static const DataValue EMPTY
Empty data value for comparisons.
Definition: DataValue.h:63
Abstract base class for FeatureFinder algorithms.
Definition: FeatureFinderAlgorithm.h:74
void setSectionDescription(const String &key, const String &description)
Sets a description for an existing section.
std::map< float, UInt > mz_stdev
Definition: FeatureFinderAlgorithm.h:55
IndexPair nextSeed()
return the next seed
Definition: SimpleSeeder.h:94
void setParameters(const Param &param)
Sets the parameters.
The purpose of this struct is to provide definitions of classes and typedefs which are used throughou...
Definition: FeatureFinderDefs.h:51
static FeatureFinderAlgorithm< PeakType, FeatureType > * create()
Definition: FeatureFinderAlgorithmSimple.h:235
DoubleReal corr_max
Definition: FeatureFinderAlgorithm.h:57
static const String getProductName()
Definition: FeatureFinderAlgorithmSimple.h:240
void endProgress() const
Ends the progress display.
const DataValue & getValue(const String &key) const
Returns a value of a parameter.
FeatureFinderAlgorithm implementation using the Simple* modules.
Definition: FeatureFinderAlgorithmSimple.h:62
std::map< String, UInt > exception
Definition: FeatureFinderAlgorithm.h:52
const char * getName() const
Returns the name of the exception.
void setDefaults(const Param &defaults, const String &prefix="", bool showMessage=false)
Insert all values of defaults and adds the prefix prefix, if the values are not already set...
std::vector< UInt > charge
Definition: FeatureFinderAlgorithm.h:56
Summary of fitting results.
Definition: FeatureFinderAlgorithm.h:50
Simple feature extension algorithm.
Definition: SimpleExtender.h:79
Tests a group of data points in an LC-MS map for goodness-of-fit with a 2D averagine model...
Definition: ModelFitter.h:77
bool exists(const String &key) const
Tests if a parameter is set (expecting its fully qualified name, e.g., TextExporter:1:proteins_only) ...
Exception that is thrown if a method an invalid IndexPair is given.
Definition: FeatureFinderDefs.h:66
An LC-MS feature.
Definition: Feature.h:66
FeatureFinderAlgorithmSimple()
default constructor
Definition: FeatureFinderAlgorithmSimple.h:69
Management and storage of parameters / INI files.
Definition: Param.h:69
DoubleReal corr_min
Definition: FeatureFinderAlgorithm.h:57
Exception used if an error occurred while fitting a model to a given dataset.
Definition: Exception.h:662
index set with associated charge estimate
Definition: IsotopeCluster.h:53
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:144
const Param & getParam() const
Non-mutable access to model parameters.
Definition: ModelDescription.h:127
bool check_defaults_
If this member is set to false no checking if parameters in done;.
Definition: DefaultParamHandler.h:174
const Flag & getPeakFlag(const IndexPair &index) const
Returns a non-mutable reference to a peak flag.
Definition: FeatureFinder.h:91
Simple seeding class that uses the strongest peak as next seed.
Definition: SimpleSeeder.h:61
Feature fit(const ChargedIndexSet &index_set)
Return next feature.
Definition: ModelFitter.h:226
const Param & getParameters() const
Non-mutable access to the parameters.
FeatureFinderAlgorithmSimple & operator=(const FeatureFinderAlgorithmSimple &)
Not implemented.
virtual void run()
Main method that implements the actual algorithm.
Definition: FeatureFinderAlgorithmSimple.h:95

OpenMS / TOPP release 1.11.1 Documentation generated on Thu Nov 14 2013 11:19:13 using doxygen 1.8.5