OpenMS  2.4.0
MapAlignmentAlgorithmIdentification.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-2018.
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: Hendrik Weisser $
32 // $Authors: Eva Lange, Clemens Groepl, Hendrik Weisser $
33 // --------------------------------------------------------------------------
34 
35 #pragma once
36 
46 
47 #include <cmath> // for "abs"
48 #include <limits> // for "max"
49 #include <map>
50 
51 namespace OpenMS
52 {
71  public DefaultParamHandler,
72  public ProgressLogger
73  {
74 public:
77 
80 
81  // Set a reference for the alignment
82  template <typename DataType> void setReference(DataType& data)
83  {
84  reference_.clear();
85  if (data.empty()) return; // empty input resets the reference
86  SeqToList rt_data;
87  bool sorted = getRetentionTimes_(data, rt_data);
88  computeMedians_(rt_data, reference_, sorted);
89  if (reference_.empty())
90  {
91  throw Exception::MissingInformation(__FILE__, __LINE__, OPENMS_PRETTY_FUNCTION, "Could not extract retention time information from the reference file");
92  }
93  }
94 
102  template <typename DataType>
103  void align(std::vector<DataType>& data,
104  std::vector<TransformationDescription>& transformations,
105  Int reference_index = -1)
106  {
107  checkParameters_(data.size());
108  startProgress(0, 3, "aligning maps");
109 
110  reference_index_ = reference_index;
111  // is reference one of the input files?
112  bool use_internal_reference = (reference_index >= 0);
113  if (use_internal_reference)
114  {
115  if (reference_index >= data.size())
116  {
117  throw Exception::IndexOverflow(__FILE__, __LINE__,
118  OPENMS_PRETTY_FUNCTION, reference_index,
119  data.size());
120  }
121  setReference(data[reference_index]);
122  }
123 
124  // one set of RT data for each input map, except reference (if any):
125  std::vector<SeqToList> rt_data(data.size() - use_internal_reference);
126  bool all_sorted = true;
127  for (Size i = 0, j = 0; i < data.size(); ++i)
128  {
129  if ((reference_index >= 0) && (i == Size(reference_index)))
130  {
131  continue; // skip reference map, if any
132  }
133  all_sorted &= getRetentionTimes_(data[i], rt_data[j++]);
134  }
135  setProgress(1);
136 
137  computeTransformations_(rt_data, transformations, all_sorted);
138  setProgress(2);
139 
140  setProgress(3);
141  endProgress();
142  }
143 
144 protected:
145 
147  typedef std::map<String, DoubleList> SeqToList;
148 
150  typedef std::map<String, double> SeqToValue;
151 
154 
157 
160 
170  void computeMedians_(SeqToList& rt_data, SeqToValue& medians,
171  bool sorted = false);
172 
181  bool getRetentionTimes_(std::vector<PeptideIdentification>& peptides,
182  SeqToList& rt_data);
183 
192  bool getRetentionTimes_(PeakMap& experiment, SeqToList& rt_data);
193 
206  template <typename MapType>
207  bool getRetentionTimes_(MapType& features, SeqToList& rt_data)
208  {
209  bool use_feature_rt = param_.getValue("use_feature_rt").toBool();
210  for (typename MapType::Iterator feat_it = features.begin();
211  feat_it != features.end(); ++feat_it)
212  {
213  if (use_feature_rt)
214  {
215  // find the peptide ID closest in RT to the feature centroid:
216  String sequence;
217  double rt_distance = std::numeric_limits<double>::max();
218  bool any_hit = false;
219  for (std::vector<PeptideIdentification>::iterator pep_it =
220  feat_it->getPeptideIdentifications().begin(); pep_it !=
221  feat_it->getPeptideIdentifications().end(); ++pep_it)
222  {
223  if (!pep_it->getHits().empty())
224  {
225  any_hit = true;
226  double current_distance = fabs(pep_it->getRT() -
227  feat_it->getRT());
228  if (current_distance < rt_distance)
229  {
230  pep_it->sort();
231  sequence = pep_it->getHits()[0].getSequence().toString();
232  rt_distance = current_distance;
233  }
234  }
235  }
236 
237  if (any_hit) rt_data[sequence].push_back(feat_it->getRT());
238  }
239  else
240  {
241  getRetentionTimes_(feat_it->getPeptideIdentifications(), rt_data);
242  }
243  }
244 
245  if (!use_feature_rt &&
246  param_.getValue("use_unassigned_peptides").toBool())
247  {
248  getRetentionTimes_(features.getUnassignedPeptideIdentifications(),
249  rt_data);
250  }
251 
252  // remove duplicates (can occur if a peptide ID was assigned to several
253  // features due to overlap or annotation tolerance):
254  for (SeqToList::iterator rt_it = rt_data.begin(); rt_it != rt_data.end();
255  ++rt_it)
256  {
257  DoubleList& rt_values = rt_it->second;
258  sort(rt_values.begin(), rt_values.end());
259  DoubleList::iterator it = unique(rt_values.begin(), rt_values.end());
260  rt_values.resize(it - rt_values.begin());
261  }
262  return true; // RTs were already sorted for duplicate detection
263  }
264 
272  void computeTransformations_(std::vector<SeqToList>& rt_data,
273  std::vector<TransformationDescription>&
274  transforms, bool sorted = false);
275 
283  void checkParameters_(const Size runs);
284 
291  void getReference_();
292 
293 private:
294 
297 
300 
301  };
302 
303 } // namespace OpenMS
304 
A more convenient string class.
Definition: String.h:57
std::vector< double > DoubleList
Vector of double precision real types.
Definition: ListUtils.h:65
std::vector< SpectrumType >::iterator Iterator
Mutable iterator.
Definition: MSExperiment.h:111
Int overflow exception.
Definition: Exception.h:254
Iterator begin()
Definition: MSExperiment.h:157
bool getRetentionTimes_(MapType &features, SeqToList &rt_data)
Collect retention time data ("RT" MetaInfo) from peptide IDs contained in feature maps or consensus m...
Definition: MapAlignmentAlgorithmIdentification.h:207
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:46
Int reference_index_
Index of input file to use as reference (if any)
Definition: MapAlignmentAlgorithmIdentification.h:153
Iterator end()
Definition: MSExperiment.h:167
A map alignment algorithm based on peptide identifications from MS2 spectra.
Definition: MapAlignmentAlgorithmIdentification.h:70
void align(std::vector< DataType > &data, std::vector< TransformationDescription > &transformations, Int reference_index=-1)
Align feature maps, consensus maps, peak maps, or peptide identifications.
Definition: MapAlignmentAlgorithmIdentification.h:103
std::map< String, DoubleList > SeqToList
Type to store retention times given for individual peptide sequences.
Definition: MapAlignmentAlgorithmIdentification.h:147
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:77
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:127
Base class for all classes that want to report their progress.
Definition: ProgressLogger.h:54
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:91
Size min_run_occur_
Minimum number of runs a peptide must occur in.
Definition: MapAlignmentAlgorithmIdentification.h:159
SeqToValue reference_
Reference retention times (per peptide sequence)
Definition: MapAlignmentAlgorithmIdentification.h:156
int Int
Signed integer type.
Definition: Types.h:102
void setReference(DataType &data)
Definition: MapAlignmentAlgorithmIdentification.h:82
std::map< String, double > SeqToValue
Type to store one representative retention time per peptide sequence.
Definition: MapAlignmentAlgorithmIdentification.h:150
Not all required information provided.
Definition: Exception.h:195