Loading...
Searching...
No Matches
experimentalspectrum.cpp
Go to the documentation of this file.
/**
* \file pappsomspp/processing/specglob/experimentalspectrum.cpp
* \date 08/11/2023
* \author Olivier Langella
* \brief transform a spectrum to SpecGlob spectra
*
* C++ implementation of the SpecGlob algorithm described in :
* 1. Prunier, G. et al. Fast alignment of mass spectra in large proteomics
* datasets, capturing dissimilarities arising from multiple complex
* modifications of peptides. BMC Bioinformatics 24, 421 (2023).
*
* HAL Id : hal-04296170 , version 1
* Mot de passe : hxo20cl
* DOI : 10.1186/s12859-023-05555-y
*/
/*
* SpecGlobTool, Spectra to peptide alignment tool
* Copyright (C) 2023 Olivier Langella
* <olivier.langella@universite-paris-saclay.fr>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*
*/
#include "experimentalspectrum.h"
#include <QDebug>
#include "pappsomspp/core/pappsoexception.h"
#include <QObject>
namespace pappso
{
namespace specglob
{
ExperimentalSpectrum::ExperimentalSpectrum(const pappso::QualifiedMassSpectrum &qmass_spectrum,
pappso::PrecisionPtr precision_ptr)
: std::vector<ExperimentalSpectrumDataPoint>()
{
m_qualifiedMassSpectrum = qmass_spectrum;
m_precisionPtr = precision_ptr;
createSymetricPeakList();
}
ExperimentalSpectrum::ExperimentalSpectrum(const pappso::QualifiedMassSpectrum &qmass_spectrum,
pappso::PrecisionPtr precision_ptr,
double precursor_mass_error)
: std::vector<ExperimentalSpectrumDataPoint>()
{
m_qualifiedMassSpectrum = qmass_spectrum;
m_precisionPtr = precision_ptr;
createSymetricPeakList(precursor_mass_error);
}
ExperimentalSpectrum::ExperimentalSpectrum(const ExperimentalSpectrum &other)
: std::vector<ExperimentalSpectrumDataPoint>(other)
{
m_targetMzSum = other.m_targetMzSum;
m_qualifiedMassSpectrum = other.m_qualifiedMassSpectrum;
m_precisionPtr = other.m_precisionPtr;
}
ExperimentalSpectrum::~ExperimentalSpectrum()
{
}
void
ExperimentalSpectrum::createSymetricPeakList()
{
createSymetricPeakList(0);
}
void
ExperimentalSpectrum::createSymetricPeakList(double precursor_mass_error)
{
bool ok;
auto charge = m_qualifiedMassSpectrum.getPrecursorCharge(&ok);
if(!ok)
{
throw pappso::PappsoException(
QObject::tr("precursor charge is not defined in spectrum %1")
.arg(m_qualifiedMassSpectrum.getMassSpectrumId().getNativeId()));
}
double mz_prec = m_qualifiedMassSpectrum.getPrecursorMz(&ok);
if(!ok)
{
throw pappso::PappsoException(
QObject::tr("precursor m/z is not defined in spectrum %1")
.arg(m_qualifiedMassSpectrum.getMassSpectrumId().getNativeId()));
}
// compute precursor mass given the charge state
m_precursorMass = mz_prec * (double)charge - precursor_mass_error;
m_precursorMass -= pappso::MHPLUS * (double)charge;
m_targetMzSum = m_precursorMass + pappso::MHPLUS + pappso::MHPLUS;
std::vector<double> mz_list = m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->xValues();
std::size_t mz_current_indice = 0;
auto itend = m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->end();
for(std::vector<pappso::DataPoint>::const_iterator it =
m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->begin();
it != itend;
it++)
{
double current_mz = it->x;
double symmetric_current_mz = getSymetricMz(current_mz);
auto itpair_symmetric = findMz(symmetric_current_mz);
if(itpair_symmetric != itend)
{
// there is a counterpart
push_back({ExperimentalSpectrumDataPointType::both, it->x, 0});
qDebug() << "current_mz=" << current_mz << " both";
}
else
{
// there is no counterpart
push_back({ExperimentalSpectrumDataPointType::native, it->x, 0});
push_back({ExperimentalSpectrumDataPointType::symmetric, symmetric_current_mz, 0});
qDebug() << "current_mz=" << current_mz << " symmetrics";
}
mz_current_indice++;
}
// we add a peak with NT mass (1.0078) if it is not detected, to give better
// chance to align the first amino acid in b if it is present
if(findMz(pappso::MPROTIUM) == itend)
{
push_back({ExperimentalSpectrumDataPointType::synthetic, pappso::MPROTIUM, 0});
}
// We add the B peak corresponding to the precursor (complete peptide)
double precusorBion = m_targetMzSum - pappso::MASSH2O - pappso::MPROTIUM;
if(findMz(precusorBion) == itend)
{
push_back({ExperimentalSpectrumDataPointType::synthetic, precusorBion, 0});
}
std::sort(begin(),
end(),
[](const ExperimentalSpectrumDataPoint &a, const ExperimentalSpectrumDataPoint &b) {
return (a.peak_mz < b.peak_mz);
});
std::size_t i = 0;
for(auto &data_point : *this)
{
data_point.indice = i;
i++;
}
}
std::vector<pappso::DataPoint>::const_iterator
ExperimentalSpectrum::findMz(double mz)
{
pappso::MzRange mz_range(mz, m_precisionPtr);
auto itend = m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->end();
auto it =
findFirstEqualOrGreaterX(m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->begin(),
m_qualifiedMassSpectrum.getMassSpectrumCstSPtr().get()->end(),
mz_range.lower());
if(it != itend)
{
if(it->x <= mz_range.upper())
{
return it;
}
}
return itend;
}
double
ExperimentalSpectrum::getTargetMzSum() const
{
return m_targetMzSum;
}
std::vector<double>
ExperimentalSpectrum::getMassList() const
{
std::vector<double> mass_list;
for(const ExperimentalSpectrumDataPoint &n : *this)
{
mass_list.push_back(n.peak_mz);
}
return mass_list;
}
double
ExperimentalSpectrum::getSymetricMz(double mz) const
{
return m_targetMzSum - mz;
}
std::vector<double>
ExperimentalSpectrum::getMassList(ExperimentalSpectrumDataPointType type) const
{
std::vector<double> mass_list;
for(const ExperimentalSpectrumDataPoint &n : *this)
{
if(n.type == type)
mass_list.push_back(n.peak_mz);
}
return mass_list;
}
QString
ExperimentalSpectrum::toString() const
{
QStringList all_element;
for(const ExperimentalSpectrumDataPoint &n : *this)
{
all_element << QString("%1 %2").arg(n.peak_mz).arg((std::uint8_t)n.type);
}
return QString("[%1]").arg(all_element.join("] ["));
}
double
ExperimentalSpectrum::getPrecursorMass() const
{
return m_precursorMass;
}
std::vector<ExperimentalSpectrumDataPoint>::const_reverse_iterator
ExperimentalSpectrum::reverseFindDiffMz(std::size_t start_position,
const pappso::MzRange &aaTheoMzRange) const
{
qDebug() << "start_position" << start_position << " lookfor=" << aaTheoMzRange.getMz();
std::vector<ExperimentalSpectrumDataPoint>::const_reverse_iterator itrbegin =
rbegin() + (size() - 1 - start_position);
qDebug() << itrbegin->indice << "mz=" << itrbegin->peak_mz;
double eperimentalMzReference = itrbegin->peak_mz;
auto itrend = this->rend();
// We check all row from j to 0
qDebug();
for(auto itr = itrbegin + 1; itr != itrend; ++itr)
{
qDebug() << itr->indice;
double experimentalMzDifference = eperimentalMzReference - itr->peak_mz;
if(experimentalMzDifference > aaTheoMzRange.upper())
{
// if we pass the mass of the theoretical amino acid, we stop to not
// over
// calculate
qDebug() << experimentalMzDifference << ">" << aaTheoMzRange.upper();
return itrend;
}
else if(experimentalMzDifference < aaTheoMzRange.lower())
{
continue;
}
qDebug() << itr->indice << " diff=" << experimentalMzDifference;
// if we found that j-k give an amino acid, we keep k value
return itr;
}
qDebug() << "rend";
return itrend; // a value of -1 to show that
// there is no k value
}
const pappso::QualifiedMassSpectrum &
ExperimentalSpectrum::getQualifiedMassSpectrum() const
{
return m_qualifiedMassSpectrum;
}
} // namespace specglob
} // namespace pappso
pappso::MzRange
Definition mzrange.h:46
pappso::MzRange::getMz
pappso_double getMz() const
Definition mzrange.cpp:109
pappso::MzRange::lower
pappso_double lower() const
Definition mzrange.h:71
pappso::MzRange::upper
pappso_double upper() const
Definition mzrange.h:77
pappso::PappsoException
Definition pappsoexception.h:41
pappso::QualifiedMassSpectrum
Class representing a fully specified mass spectrum.
Definition qualifiedmassspectrum.h:99
pappso::specglob::ExperimentalSpectrum::getSymetricMz
double getSymetricMz(double mz) const
compute the symmetric mass for debuggin purpose
Definition experimentalspectrum.cpp:207
pappso::specglob::ExperimentalSpectrum::getMassList
std::vector< double > getMassList() const
Definition experimentalspectrum.cpp:195
pappso::specglob::ExperimentalSpectrum::~ExperimentalSpectrum
virtual ~ExperimentalSpectrum()
Definition experimentalspectrum.cpp:72
pappso::specglob::ExperimentalSpectrum::getPrecursorMass
double getPrecursorMass() const
Definition experimentalspectrum.cpp:237
pappso::specglob::ExperimentalSpectrum::getQualifiedMassSpectrum
const pappso::QualifiedMassSpectrum & getQualifiedMassSpectrum() const
Definition experimentalspectrum.cpp:284
pappso::specglob::ExperimentalSpectrum::getTargetMzSum
double getTargetMzSum() const
Definition experimentalspectrum.cpp:189
pappso::specglob::ExperimentalSpectrum::m_precisionPtr
pappso::PrecisionPtr m_precisionPtr
Definition experimentalspectrum.h:141
pappso::specglob::ExperimentalSpectrum::m_precursorMass
double m_precursorMass
Definition experimentalspectrum.h:140
pappso::specglob::ExperimentalSpectrum::reverseFindDiffMz
std::vector< ExperimentalSpectrumDataPoint >::const_reverse_iterator reverseFindDiffMz(std::size_t start_position, const pappso::MzRange &targeted_mass_range) const
find the peak for wich mass difference from rbegin corresponds to aaTheoMass Find if a peak back in t...
Definition experimentalspectrum.cpp:243
pappso::specglob::ExperimentalSpectrum::toString
QString toString() const
Definition experimentalspectrum.cpp:226
pappso::specglob::ExperimentalSpectrum::ExperimentalSpectrum
ExperimentalSpectrum(const pappso::QualifiedMassSpectrum &qmass_spectrum, pappso::PrecisionPtr precision_ptr)
Definition experimentalspectrum.cpp:45
pappso::specglob::ExperimentalSpectrum::createSymetricPeakList
void createSymetricPeakList()
add symmetric peaks to the spectrum Create a SymetricPeakList that contain symmetric peaks and the in...
Definition experimentalspectrum.cpp:77
pappso::specglob::ExperimentalSpectrum::findMz
std::vector< pappso::DataPoint >::const_iterator findMz(double mz)
find the correspondin mz in the mass spectrum (given the precision)
Definition experimentalspectrum.cpp:169
pappso::specglob::ExperimentalSpectrum::m_targetMzSum
double m_targetMzSum
Definition experimentalspectrum.h:139
pappso::specglob::ExperimentalSpectrum::m_qualifiedMassSpectrum
pappso::QualifiedMassSpectrum m_qualifiedMassSpectrum
Definition experimentalspectrum.h:138
experimentalspectrum.h
pappso::specglob
Definition experimentalspectrum.cpp:44
pappso::specglob::ExperimentalSpectrumDataPointType
ExperimentalSpectrumDataPointType
Definition types.h:78
pappso::specglob::ExperimentalSpectrumDataPointType::synthetic
@ synthetic
does not correspond to existing peak, for computational purpose
Definition types.h:85
pappso::specglob::ExperimentalSpectrumDataPointType::native
@ native
Definition types.h:79
pappso::specglob::ExperimentalSpectrumDataPointType::both
@ both
both, the ion and the complement exists in the original spectrum
Definition types.h:83
pappso::specglob::ExperimentalSpectrumDataPointType::symmetric
@ symmetric
new peak : computed symmetric mass from a corresponding native peak
Definition types.h:81
pappso
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition aa.cpp:39
pappso::findFirstEqualOrGreaterX
std::vector< DataPoint >::iterator findFirstEqualOrGreaterX(std::vector< DataPoint >::iterator begin, std::vector< DataPoint >::iterator end, const double &value)
find the first element in which X is equal or greater than the value searched important : it implies ...
Definition trace.cpp:65
pappso::PeptideIonNter::a
@ a
Definition peptide.h:77
pappso::PeptideIonNter::b
@ b
Definition peptide.h:74
pappso::MHPLUS
const pappso_double MHPLUS(1.007276466879)
pappso::MPROTIUM
const pappso_double MPROTIUM(1.007825032241)
pappso::MASSH2O
const pappso_double MASSH2O((MPROTIUM *2)+MASSOXYGEN)
pappso::PrecisionPtr
const PrecisionBase * PrecisionPtr
Definition precision.h:122
pappsoexception.h
pappso::specglob::ExperimentalSpectrumDataPoint
Definition experimentalspectrum.h:51
pappso::specglob::ExperimentalSpectrumDataPoint::peak_mz
double peak_mz
Definition experimentalspectrum.h:53
pappso::specglob::ExperimentalSpectrumDataPoint::type
ExperimentalSpectrumDataPointType type
Definition experimentalspectrum.h:52