FLASHDeconv algorithm: ultrafast mass deconvolution algorithm for top down mass spectrometry dataset From MSSpectrum, this class outputs DeconvolvedSpectrum. Deconvolution takes three steps: i) decharging and select candidate masses - speed up via binning ii) collecting isotopes from the candidate masses and deisotope - peak groups are defined here iii) scoring and filter out low scoring masses (i.e., peak groups)
More...
#include <OpenMS/ANALYSIS/TOPDOWN/SpectralDeconvolution.h>
|
| static int | getNominalMass (double mass) |
| | convert double to nominal mass More...
|
| |
| static float | getCosine (const std::vector< float > &a, int a_start, int a_end, const IsotopeDistribution &b, int offset, int min_iso_len) |
| |
| static float | getIsotopeCosineAndIsoOffset (double mono_mass, const std::vector< float > &per_isotope_intensities, int &offset, const PrecalculatedAveragine &avg, int iso_int_shift=0, int window_width=-1, int allowed_isotope_error=0, PeakGroup::TargetDecoyType target_decoy_type=PeakGroup::TargetDecoyType::target) |
| | Examine intensity distribution over isotope indices. Also determines the most plausible isotope index or, monoisotopic mono_mass. More...
|
| |
Static Public Member Functions inherited from DefaultParamHandler |
| static void | writeParametersToMetaValues (const Param &write_this, MetaInfoInterface &write_here, const String &key_prefix="") |
| | Writes all parameters to meta values. More...
|
| |
|
| void | updateLogMzPeaks_ () |
| | generate log mz peaks from the input spectrum More...
|
| |
| void | binLogMzPeaks_ (const Size bin_number, std::vector< float > &binned_log_mz_peak_intensities) |
| | generate mz bins and intensity per mz bin from log mz peaks More...
|
| |
| double | getMassFromMassBin_ (Size mass_bin, double bin_mul_factor) const |
| | get mass value for input mass bin More...
|
| |
| double | getMzFromMzBin_ (Size mass_bin, double bin_mul_factor) const |
| | get mz value for input mz bin More...
|
| |
| void | generatePeakGroupsFromSpectrum_ () |
| | Generate peak groups from the input spectrum. More...
|
| |
| void | removeOverlappingPeakGroups_ (DeconvolvedSpectrum &dspec, double tol, PeakGroup::TargetDecoyType target_decoy_type=PeakGroup::TargetDecoyType::target) |
| | filter out overlapping masses More...
|
| |
| Matrix< int > | updateMassBins_ (const std::vector< float > &mz_intensities) |
| | Update binned_log_masses_. It select candidate mass bins using the universal pattern, eliminate possible harmonic masses. This function does not perform deisotoping. More...
|
| |
| Matrix< int > | filterMassBins_ (const std::vector< float > &mass_intensities) |
| | Subfunction of updateMassBins_. More...
|
| |
| void | updateCandidateMassBins_ (std::vector< float > &mass_intensities, const std::vector< float > &mz_intensities) |
| | Subfunction of updateMassBins_. It select candidate masses and update binned_log_masses_ using the universal pattern, eliminate possible harmonic masses. More...
|
| |
| void | getCandidatePeakGroups_ (const Matrix< int > &per_mass_abs_charge_ranges) |
| | For selected masses in binned_log_masses_, select the peaks from the original spectrum. Also isotopic peaks are clustered in this function. More...
|
| |
| void | setFilters_ () |
| | Make the universal pattern. More...
|
| |
| void | scoreAndFilterPeakGroups_ () |
| | function for peak group scoring and filtering More...
|
| |
| void | removeChargeErrorPeakGroups_ (DeconvolvedSpectrum &dspec, const PeakGroup::TargetDecoyType &target_decoy_type) const |
| | filter out charge error masses More...
|
| |
| void | removeExcludedMasses_ (DeconvolvedSpectrum &dspec, std::vector< double > excluded_masses) const |
| | filter out excluded masses More...
|
| |
| void | setTargetPrecursorCharge_ () |
| |
|
| static double | getBinValue_ (Size bin, double min_value, double bin_mul_factor) |
| | static function that converts bin to value More...
|
| |
| static Size | getBinNumber_ (double value, double min_value, double bin_mul_factor) |
| | static function that converts value to bin More...
|
| |
FLASHDeconv algorithm: ultrafast mass deconvolution algorithm for top down mass spectrometry dataset From MSSpectrum, this class outputs DeconvolvedSpectrum. Deconvolution takes three steps: i) decharging and select candidate masses - speed up via binning ii) collecting isotopes from the candidate masses and deisotope - peak groups are defined here iii) scoring and filter out low scoring masses (i.e., peak groups)
◆ LogMzPeak
◆ PrecalculatedAveragine
◆ SpectralDeconvolution() [1/3]
◆ SpectralDeconvolution() [2/3]
◆ SpectralDeconvolution() [3/3]
◆ ~SpectralDeconvolution()
◆ binLogMzPeaks_()
| void binLogMzPeaks_ |
( |
const Size |
bin_number, |
|
|
std::vector< float > & |
binned_log_mz_peak_intensities |
|
) |
| |
|
private |
generate mz bins and intensity per mz bin from log mz peaks
- Parameters
-
| bin_number | number of mz bins |
| binned_log_mz_peak_intensities | intensity per mz bin |
◆ calculateAveragine()
| void calculateAveragine |
( |
bool |
use_RNA_averagine, |
|
|
const bool |
is_centroid = true |
|
) |
| |
precalculate averagine (for predefined mass bins) to speed up averagine generation
- Parameters
-
| use_RNA_averagine | if set, averagine for RNA (nucleotides) is calculated |
| is_centroid | this is for noise averagine calculation. For noise, centroid and profile averagines are different. |
◆ filterMassBins_()
| Matrix<int> filterMassBins_ |
( |
const std::vector< float > & |
mass_intensities | ) |
|
|
private |
Subfunction of updateMassBins_.
- Parameters
-
| mass_intensities | per mass bin intensity |
- Returns
- a matrix containing charge ranges for all found masses
◆ generatePeakGroupsFromSpectrum_()
| void generatePeakGroupsFromSpectrum_ |
( |
| ) |
|
|
private |
Generate peak groups from the input spectrum.
◆ getAveragine()
get calculated averagine. Call after calculateAveragine is called.
◆ getBinNumber_()
| static Size getBinNumber_ |
( |
double |
value, |
|
|
double |
min_value, |
|
|
double |
bin_mul_factor |
|
) |
| |
|
staticprivate |
static function that converts value to bin
- Parameters
-
| value | value |
| min_value | minimum value (corresponding to bin number = 0) |
| bin_mul_factor | bin multiplication factor: bin_number = (bin_value * bin_mul_factors_ - min_value) |
- Returns
- bin corresponding to value
◆ getBinValue_()
| static double getBinValue_ |
( |
Size |
bin, |
|
|
double |
min_value, |
|
|
double |
bin_mul_factor |
|
) |
| |
|
staticprivate |
static function that converts bin to value
- Parameters
-
| bin | bin number |
| min_value | minimum value (corresponding to bin number = 0) |
| bin_mul_factor | bin multiplication factor: bin_value = (min_value + bin_number/ bin_mul_factors_) |
- Returns
- value corresponding to bin
◆ getCandidatePeakGroups_()
| void getCandidatePeakGroups_ |
( |
const Matrix< int > & |
per_mass_abs_charge_ranges | ) |
|
|
private |
For selected masses in binned_log_masses_, select the peaks from the original spectrum. Also isotopic peaks are clustered in this function.
- Parameters
-
| per_mass_abs_charge_ranges | charge range per mass |
◆ getCosine()
| static float getCosine |
( |
const std::vector< float > & |
a, |
|
|
int |
a_start, |
|
|
int |
a_end, |
|
|
const IsotopeDistribution & |
b, |
|
|
int |
offset, |
|
|
int |
min_iso_len |
|
) |
| |
|
static |
calculate cosine between two vectors a and b with additional parameters for fast calculation
- Parameters
-
| a | vector a |
| a_start | non zero start index of a |
| a_end | non zero end index of a (exclusive) |
| b | vector b |
| offset | element index offset between a and b |
| min_iso_len | minimum isotope size. If isotope size is less than this, return 0 |
◆ getDeconvolvedSpectrum()
return deconvolved spectrum
◆ getIsotopeCosineAndIsoOffset()
| static float getIsotopeCosineAndIsoOffset |
( |
double |
mono_mass, |
|
|
const std::vector< float > & |
per_isotope_intensities, |
|
|
int & |
offset, |
|
|
const PrecalculatedAveragine & |
avg, |
|
|
int |
iso_int_shift = 0, |
|
|
int |
window_width = -1, |
|
|
int |
allowed_isotope_error = 0, |
|
|
PeakGroup::TargetDecoyType |
target_decoy_type = PeakGroup::TargetDecoyType::target |
|
) |
| |
|
static |
Examine intensity distribution over isotope indices. Also determines the most plausible isotope index or, monoisotopic mono_mass.
- Parameters
-
| mono_mass | monoisotopic mass |
| per_isotope_intensities | vector of intensities associated with each isotope - aggregated through charges |
| offset | output offset between input monoisotopic mono_mass and determined monoisotopic mono_mass |
| avg | precalculated averagine |
| iso_int_shift | isotope shift in per_isotope_intensities. |
| window_width | isotope offset value range. If -1, set automatically. |
| allowed_isotope_error | allowed isotope error to calculate the second best cos. If target_decoy_type is not PeakGroup::TargetDecoyType::target, the second best cosine and its corresponding offset will be output |
| target_decoy_type | This target_decoy_type specifies if a PeakGroup is a target (0), charge dummy (1), noise dummy (2), or isotope dummy (3) |
- Returns
- calculated cosine similarity score
◆ getMassFromMassBin_()
| double getMassFromMassBin_ |
( |
Size |
mass_bin, |
|
|
double |
bin_mul_factor |
|
) |
| const |
|
private |
get mass value for input mass bin
◆ getMzFromMzBin_()
| double getMzFromMzBin_ |
( |
Size |
mass_bin, |
|
|
double |
bin_mul_factor |
|
) |
| const |
|
private |
get mz value for input mz bin
◆ getNominalMass()
| static int getNominalMass |
( |
double |
mass | ) |
|
|
static |
convert double to nominal mass
◆ operator=() [1/2]
◆ operator=() [2/2]
◆ performSpectrumDeconvolution()
| void performSpectrumDeconvolution |
( |
const MSSpectrum & |
spec, |
|
|
int |
scan_number, |
|
|
const PeakGroup & |
precursor_peak_group |
|
) |
| |
main deconvolution function that generates the deconvolved target and dummy spectrum based on the original spectrum.
- Parameters
-
| spec | the original spectrum |
| scan_number | scan number from input spectrum. |
| precursor_peak_group | precursor peak group |
◆ removeChargeErrorPeakGroups_()
filter out charge error masses
◆ removeExcludedMasses_()
| void removeExcludedMasses_ |
( |
DeconvolvedSpectrum & |
dspec, |
|
|
std::vector< double > |
excluded_masses |
|
) |
| const |
|
private |
filter out excluded masses
◆ removeOverlappingPeakGroups_()
filter out overlapping masses
◆ scoreAndFilterPeakGroups_()
| void scoreAndFilterPeakGroups_ |
( |
| ) |
|
|
private |
function for peak group scoring and filtering
◆ setAveragine()
◆ setFilters_()
Make the universal pattern.
◆ setTargetDecoyType()
set target dummy type for the SpectralDeconvolution run. All masses from the target SpectralDeconvolution run will have the target_decoy_type_.
- Parameters
-
| target_decoy_type | This target_decoy_type_ specifies if a PeakGroup is a target (0), charge dummy (1), noise dummy (2), or isotope dummy (3) |
| target_dspec_for_decoy_calcualtion | target masses from normal deconvolution |
◆ setTargetMasses()
| void setTargetMasses |
( |
const std::vector< double > & |
masses, |
|
|
bool |
exclude = false |
|
) |
| |
set targeted or excluded masses for targeted deconvolution. Masses are targeted or excluded in all ms levels.
- Parameters
-
| masses | target masses to set |
| exclude | if set, masses are excluded. |
◆ setTargetPrecursorCharge_()
| void setTargetPrecursorCharge_ |
( |
| ) |
|
|
private |
◆ setToleranceEstimation()
| void setToleranceEstimation |
( |
| ) |
|
|
inline |
when estimating tolerance, max_mass_dalton_tolerance_ should be large
◆ updateCandidateMassBins_()
| void updateCandidateMassBins_ |
( |
std::vector< float > & |
mass_intensities, |
|
|
const std::vector< float > & |
mz_intensities |
|
) |
| |
|
private |
Subfunction of updateMassBins_. It select candidate masses and update binned_log_masses_ using the universal pattern, eliminate possible harmonic masses.
- Parameters
-
| mass_intensities | mass bin intensities which are updated in this function |
| mz_intensities | mz bin intensities |
◆ updateLogMzPeaks_()
| void updateLogMzPeaks_ |
( |
| ) |
|
|
private |
generate log mz peaks from the input spectrum
◆ updateMassBins_()
| Matrix<int> updateMassBins_ |
( |
const std::vector< float > & |
mz_intensities | ) |
|
|
private |
Update binned_log_masses_. It select candidate mass bins using the universal pattern, eliminate possible harmonic masses. This function does not perform deisotoping.
- Parameters
-
| mz_intensities | per mz bin intensity |
- Returns
- a matrix containing charge ranges for all found masses
◆ updateMembers_()
This method is used to update extra member variables at the end of the setParameters() method.
Also call it at the end of the derived classes' copy constructor and assignment operator.
The default implementation is empty.
Reimplemented from DefaultParamHandler.
◆ allowed_iso_error_
| int allowed_iso_error_ = 1 |
|
private |
FLASHDeconv parameters.
allowed isotope error in deconvolved mass to calculate qvalue
◆ avg_
precalculated averagine distributions for fast averagine generation
◆ bin_mul_factors_
bin multiplication factor (log mz * bin_mul_factors_ = bin number) - for fast convolution, binning is used
◆ binned_harmonic_patterns
| Matrix<int> binned_harmonic_patterns |
|
private |
This stores the patterns for harmonic reduction in binned dimension.
◆ binned_log_masses_
| boost::dynamic_bitset binned_log_masses_ |
|
private |
binned_log_masses_ stores the selected bins for this spectrum + overlapped spectrum (previous a few spectra).
◆ binned_log_mz_peaks_
| boost::dynamic_bitset binned_log_mz_peaks_ |
|
private |
binned_log_mz_peaks_ stores the binned log mz peaks
◆ binned_universal_pattern_
| std::vector<int> binned_universal_pattern_ |
|
private |
This stores the "universal pattern" in binned dimension.
◆ current_max_charge_
current_max_charge_: controlled by precursor charge for MSn n>1; otherwise just max_abs_charge_
◆ current_max_mass_
max mass is controlled by precursor mass for MSn n>1; otherwise just max_mass
◆ current_min_mass_
max mass is max_mass for MS1 and 50 for MS2
◆ deconvolved_spectrum_
selected_peak_groups_ stores the deconvolved mass peak groups
◆ excluded_masses_
| std::vector<double> excluded_masses_ |
|
private |
mass bins that are excluded for FLASHIda global targeting mode
◆ harmonic_pattern_matrix_
| Matrix<double> harmonic_pattern_matrix_ |
|
private |
This stores the patterns for harmonic reduction.
◆ is_positive_
◆ iso_da_distance_
◆ log_mz_peaks_
◆ mass_bin_min_value_
| double mass_bin_min_value_ |
|
private |
minimum mass and mz values representing the first bin of massBin and mzBin, respectively: to save memory space
◆ max_abs_charge_
◆ max_mass_
◆ max_mass_dalton_tolerance_
| double max_mass_dalton_tolerance_ = .16 |
|
private |
this is additional mass tolerance in Da to get more high signal-to-ratio peaks in this candidate peakgroup finding
◆ min_abs_charge_
min charge and max charge subject to analysis, set by users
◆ min_iso_size
| const int min_iso_size = 2 |
|
static |
minimum isotopologue count in a peak group
◆ min_isotope_cosine_
Isotope cosine threshold for each MS level.
◆ min_mass_
mass ranges of deconvolution, set by users
◆ min_snr_
SNR threshold for each MS level.
◆ min_support_peak_count_
| const int min_support_peak_count_ = 2 |
|
staticprivate |
minimum number of peaks supporting a mass minus one
◆ ms_level_
◆ mz_bin_min_value_
◆ noise_iso_delta_
| double noise_iso_delta_ = .9444 |
|
private |
isotope distance for noise decoy
◆ previously_deconved_mass_bins_for_decoy_
| boost::dynamic_bitset previously_deconved_mass_bins_for_decoy_ |
|
private |
mass bins that are previously deconvolved and excluded for dummy mass generation
◆ previously_deconved_peak_masses_for_decoy_
| std::vector<double> previously_deconved_peak_masses_for_decoy_ |
|
private |
◆ target_decoy_type_
◆ target_dspec_for_decoy_calculation_
the peak group vector from normal run. This is used when dummy masses are generated.
◆ target_mass_bins_
| boost::dynamic_bitset target_mass_bins_ |
|
private |
mass bins that are targeted for FLASHIda global targeting mode
◆ target_mono_masses_
| std::vector<double> target_mono_masses_ |
|
private |
◆ target_precursor_charge_
| int target_precursor_charge_ = 0 |
|
private |
◆ target_precursor_mz_
| double target_precursor_mz_ = 0 |
|
private |
◆ tolerance_
tolerance in ppm for each MS level
◆ universal_pattern_
| std::vector<double> universal_pattern_ |
|
private |
This stores the "universal pattern".