| OpenMS 2.5.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-2020.
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
37 #include <OpenMS/KERNEL/StandardTypes.h>
38 #include <OpenMS/ANALYSIS/MAPMATCHING/TransformationDescription.h>
39 #include <OpenMS/CONCEPT/ProgressLogger.h>
40 #include <OpenMS/DATASTRUCTURES/DefaultParamHandler.h>
41 #include <OpenMS/DATASTRUCTURES/ListUtils.h>
42 #include <OpenMS/KERNEL/ConsensusMap.h>
43 #include <OpenMS/KERNEL/FeatureMap.h>
44 #include <OpenMS/KERNEL/MSExperiment.h>
45 #include <OpenMS/METADATA/PeptideIdentification.h>
46
47 #include <cmath> // for "abs"
48 #include <limits> // for "max"
49 #include <map>
50
51 namespace OpenMS
52 {
71 class OPENMS_DLLAPI MapAlignmentAlgorithmIdentification :
72 public DefaultParamHandler,
73 public ProgressLogger
74 {
75 public:
77 MapAlignmentAlgorithmIdentification();
78
80 ~MapAlignmentAlgorithmIdentification() override;
81
82 // Set a reference for the alignment
83 template <typename DataType> void setReference(DataType& data)
84 {
85 reference_.clear();
86 if (data.empty()) return; // empty input resets the reference
87 use_feature_rt_ = param_.getValue("use_feature_rt").toBool();
88 SeqToList rt_data;
89 bool sorted = getRetentionTimes_(data, rt_data);
90 computeMedians_(rt_data, reference_, sorted);
91 if (reference_.empty())
92 {
93 throw Exception::MissingInformation(__FILE__, __LINE__, OPENMS_PRETTY_FUNCTION, "Could not extract retention time information from the reference file");
94 }
95 }
96
104 template <typename DataType>
105 void align(std::vector<DataType>& data,
106 std::vector<TransformationDescription>& transformations,
107 Int reference_index = -1)
108 {
109 checkParameters_(data.size());
110 startProgress(0, 3, "aligning maps");
111
112 reference_index_ = reference_index;
113 // is reference one of the input files?
114 bool use_internal_reference = (reference_index >= 0);
115 if (use_internal_reference)
116 {
117 if (reference_index >= static_cast<Int>(data.size()))
118 {
119 throw Exception::IndexOverflow(__FILE__, __LINE__,
120 OPENMS_PRETTY_FUNCTION, reference_index,
121 data.size());
122 }
123 setReference(data[reference_index]);
124 }
125
126 // one set of RT data for each input map, except reference (if any):
127 std::vector<SeqToList> rt_data(data.size() - use_internal_reference);
128 bool all_sorted = true;
129 for (Size i = 0, j = 0; i < data.size(); ++i)
130 {
131 if ((reference_index >= 0) && (i == Size(reference_index)))
132 {
133 continue; // skip reference map, if any
134 }
135 all_sorted &= getRetentionTimes_(data[i], rt_data[j++]);
136 }
137 setProgress(1);
138
139 computeTransformations_(rt_data, transformations, all_sorted);
140 setProgress(2);
141
142 setProgress(3);
143 endProgress();
144 }
145
146 protected:
147
149 typedef std::map<String, DoubleList> SeqToList;
150
152 typedef std::map<String, double> SeqToValue;
153
155 Int reference_index_;
156
158 SeqToValue reference_;
159
161 Size min_run_occur_;
162
164 bool use_feature_rt_;
165
167 double min_score_;
168
170 bool score_cutoff_;
171
173 bool (*better_) (double,double) = [](double, double) {return true;};
174
184 void computeMedians_(SeqToList& rt_data, SeqToValue& medians,
185 bool sorted = false);
186
195 bool getRetentionTimes_(std::vector<PeptideIdentification>& peptides,
196 SeqToList& rt_data);
197
206 bool getRetentionTimes_(PeakMap& experiment, SeqToList& rt_data);
207
222 template <typename MapType>
223 bool getRetentionTimes_(MapType& features, SeqToList& rt_data)
224 {
225 if (!score_cutoff_)
226 {
227 better_ = [](double, double)
228 {return true;};
229 }
230 else if (features[0].getPeptideIdentifications()[0].isHigherScoreBetter())
231 {
232 better_ = [](double a, double b)
233 { return a >= b; };
234 }
235 else
236 {
237 better_ = [](double a, double b)
238 { return a <= b; };
239 }
240
241 for (typename MapType::Iterator feat_it = features.begin();
242 feat_it != features.end(); ++feat_it)
243 {
244 if (use_feature_rt_)
245 {
246 // find the peptide ID closest in RT to the feature centroid:
247 String sequence;
248 double rt_distance = std::numeric_limits<double>::max();
249 bool any_hit = false;
250 for (std::vector<PeptideIdentification>::iterator pep_it =
251 feat_it->getPeptideIdentifications().begin(); pep_it !=
252 feat_it->getPeptideIdentifications().end(); ++pep_it)
253 {
254 if (!pep_it->getHits().empty())
255 {
256 any_hit = true;
257 double current_distance = fabs(pep_it->getRT() -
258 feat_it->getRT());
259 if (current_distance < rt_distance)
260 {
261 pep_it->sort();
262 if (better_(pep_it->getHits()[0].getScore(), min_score_))
263 {
264 sequence = pep_it->getHits()[0].getSequence().toString();
265 rt_distance = current_distance;
266 }
267 }
268 }
269 }
270
271 if (any_hit) rt_data[sequence].push_back(feat_it->getRT());
272 }
273 else
274 {
275 getRetentionTimes_(feat_it->getPeptideIdentifications(), rt_data);
276 }
277 }
278
279 if (!use_feature_rt_ &&
280 param_.getValue("use_unassigned_peptides").toBool())
281 {
282 getRetentionTimes_(features.getUnassignedPeptideIdentifications(),
283 rt_data);
284 }
285
286 // remove duplicates (can occur if a peptide ID was assigned to several
287 // features due to overlap or annotation tolerance):
288 for (SeqToList::iterator rt_it = rt_data.begin(); rt_it != rt_data.end();
289 ++rt_it)
290 {
291 DoubleList& rt_values = rt_it->second;
292 sort(rt_values.begin(), rt_values.end());
293 DoubleList::iterator it = unique(rt_values.begin(), rt_values.end());
294 rt_values.resize(it - rt_values.begin());
295 }
296 return true; // RTs were already sorted for duplicate detection
297 }
298
306 void computeTransformations_(std::vector<SeqToList>& rt_data,
307 std::vector<TransformationDescription>&
308 transforms, bool sorted = false);
309
317 void checkParameters_(const Size runs);
318
325 void getReference_();
326
327 private:
328
330 MapAlignmentAlgorithmIdentification(const MapAlignmentAlgorithmIdentification&);
331
333 MapAlignmentAlgorithmIdentification& operator=(const MapAlignmentAlgorithmIdentification&);
334
335 };
336
337 } // namespace OpenMS
338
OpenMS::MapAlignmentAlgorithmIdentification::score_cutoff_
bool score_cutoff_
Actually use the above defined score_cutoff? Needed since it is hard to define a non-cutting score fo...
Definition: MapAlignmentAlgorithmIdentification.h:170
OpenMS::String
A more convenient string class.
Definition: String.h:58
OpenMS::DefaultParamHandler
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:91
OpenMS::MapAlignmentAlgorithmIdentification::align
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:105
OpenMS::MapAlignmentAlgorithmIdentification::reference_
SeqToValue reference_
Reference retention times (per peptide sequence)
Definition: MapAlignmentAlgorithmIdentification.h:158
OpenMS::MapAlignmentAlgorithmIdentification::getRetentionTimes_
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:223
OpenMS::MSExperiment
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:77
OpenMS::MapAlignmentAlgorithmIdentification::min_run_occur_
Size min_run_occur_
Minimum number of runs a peptide must occur in.
Definition: MapAlignmentAlgorithmIdentification.h:161
StandardTypes.h
OpenMS::MSExperiment::Iterator
std::vector< SpectrumType >::iterator Iterator
Mutable iterator.
Definition: MSExperiment.h:111
OpenMS::MapAlignmentAlgorithmIdentification::reference_index_
Int reference_index_
Index of input file to use as reference (if any)
Definition: MapAlignmentAlgorithmIdentification.h:155
ConsensusMap.h
FeatureMap.h
OpenMS::MapAlignmentAlgorithmIdentification::use_feature_rt_
bool use_feature_rt_
Use feature RT instead of RT from best peptide ID in the feature.
Definition: MapAlignmentAlgorithmIdentification.h:164
ProgressLogger.h
OpenMS::DoubleList
std::vector< double > DoubleList
Vector of double precision real types.
Definition: ListUtils.h:62
PeptideIdentification.h
OpenMS
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:46
OpenMS::MapAlignmentAlgorithmIdentification::SeqToValue
std::map< String, double > SeqToValue
Type to store one representative retention time per peptide sequence.
Definition: MapAlignmentAlgorithmIdentification.h:152
OpenMS::Exception::IndexOverflow
Int overflow exception.
Definition: Exception.h:254
MSExperiment.h
OpenMS::ProgressLogger
Base class for all classes that want to report their progress.
Definition: ProgressLogger.h:54
OpenMS::Exception::MissingInformation
Not all required information provided.
Definition: Exception.h:195
TransformationDescription.h
DefaultParamHandler.h
double
OpenMS::MapAlignmentAlgorithmIdentification
A map alignment algorithm based on peptide identifications from MS2 spectra.
Definition: MapAlignmentAlgorithmIdentification.h:71
int
OpenMS::MSExperiment::begin
Iterator begin()
Definition: MSExperiment.h:157
OpenMS::MapAlignmentAlgorithmIdentification::min_score_
double min_score_
Minimum score to reach for a peptide to be considered.
Definition: MapAlignmentAlgorithmIdentification.h:167
ListUtils.h
OpenMS::MSExperiment::end
Iterator end()
Definition: MSExperiment.h:167
OpenMS::MapAlignmentAlgorithmIdentification::setReference
void setReference(DataType &data)
Definition: MapAlignmentAlgorithmIdentification.h:83
OpenMS::Size
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:127
OpenMS::MapAlignmentAlgorithmIdentification::SeqToList
std::map< String, DoubleList > SeqToList
Type to store retention times given for individual peptide sequences.
Definition: MapAlignmentAlgorithmIdentification.h:149