libpappsomspp
Library for mass spectrometry
Loading...
Searching...
No Matches
pappso::specpeptidoms::SemiGlobalAlignment Class Reference

#include <semiglobalalignment.h>

Public Member Functions

 SemiGlobalAlignment (const ScoreValues &score_values, const pappso::PrecisionPtr precision_ptr, const AaCode &aaCode)
 ~SemiGlobalAlignment ()
void fastAlign (const SpOMSSpectrum &spectrum, const SpOMSProtein *protein_ptr)
 perform the first alignment search between a protein sequence and a spectrum. The member location heap is filled with the candidates locations.
void preciseAlign (const SpOMSSpectrum &spectrum, const SpOMSProtein *protein_ptr, const std::size_t beginning, const std::size_t length)
 performs the second alignment search between a protein subsequence and a spectrum.
void postProcessingAlign (const SpOMSSpectrum &spectrum, const SpOMSProtein *protein_ptr, std::size_t beginning, std::size_t length, const std::vector< double > &shifts)
 performs the post-processing : generates corrected spectra and align them
LocationSaver getLocationSaver () const
 Returns a copy of m_location_saver.
Scenario getScenario () const
 Returns a copy of m_scenario.
const Alignment & getBestAlignment () const
 Returns a const ref to m_best_alignment.
const std::vector< KeyCell > & getInterestCells () const
 convenient function for degub purpose
const std::vector< KeyCell > & oneAlignStep (const pappso::specpeptidoms::SpOMSProtein &sequence, const std::size_t row_number, const std::vector< AaPosition > &aa_positions, const SpOMSSpectrum &spectrum, const bool fast_align, const pappso::specpeptidoms::SpOMSProtein *protein_ptr)
 function made for testing the fastAlign process, process one line and return the alignment matrix
void initFastAlign (const SpOMSSpectrum &spectrum)
 function made for testing the fastAlign process, initiate the variables for alignment

Static Public Member Functions

static std::vector< double > getPotentialMassErrors (const pappso::AaCode &aa_code, const Alignment &alignment, const QString &protein_seq)
 Returns a list of the potential mass errors corresponding to the provided alignment in the provided protein sequence.
static bool checkSequenceDiversity (const QString &sequence, std::size_t window, std::size_t minimum_aa_diversity)
 check that the sequence has a minimum of amino acid checkSequenceDiversity

Private Member Functions

void saveBestAlignment (const SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, std::size_t offset)
 Stores the best alignment from m_scenario in m_best_alignment.
void correctAlign (const SpOMSProtein &protein_subseq, const SpOMSProtein *protein_ptr, const SpOMSSpectrum &spectrum, std::vector< std::size_t > &peaks_to_remove, std::size_t offset)
 Recursively performs the correction of the alignment.
void updateAlignmentMatrix (const pappso::specpeptidoms::SpOMSProtein &sequence, const std::size_t row_number, const std::vector< AaPosition > &aa_positions, const SpOMSSpectrum &spectrum, const bool fast_align, const pappso::specpeptidoms::SpOMSProtein *protein_ptr)
 updates the scores of the alignment matrix for a given amino acid as well as the location heap/scenario.
bool perfectShiftPossible (const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, const std::size_t origin_row, const std::size_t current_row, const std::size_t l_peak, const std::size_t r_peak) const
 indicates if a perfect shift is possible between the provided positions
std::size_t perfectShiftPossibleFrom0 (const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, const std::size_t current_row, const std::size_t r_peak) const
 indicates if a perfect shift is possible from the spectrum beginning to the provided peak. Returns the perfect shift origin if the shift is possible, otherwise returns the current row.
std::size_t perfectShiftPossibleEnd (const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, std::size_t end_row, std::size_t end_peak) const
 indicates if a perfect shift is possible between the provided positions

Private Attributes

std::vector< KeyCell > m_interest_cells
std::vector< std::pair< std::size_t, KeyCell > > m_updated_cells
const ScoreValues & m_scorevalues
const int min_score = 15
pappso::PrecisionPtr m_precision_ptr
const AaCode & m_aaCode
LocationSaver m_location_saver
Scenario m_scenario
Alignment m_best_alignment
Alignment m_best_corrected_alignment
Alignment m_best_post_processed_alignment

Detailed Description

Definition at line 91 of file semiglobalalignment.h.

Constructor & Destructor Documentation

◆ SemiGlobalAlignment()

pappso::specpeptidoms::SemiGlobalAlignment::SemiGlobalAlignment ( const ScoreValues & score_values,
const pappso::PrecisionPtr precision_ptr,
const AaCode & aaCode )

Default constructor

Definition at line 79 of file semiglobalalignment.cpp.

83 : m_scorevalues(score_values), m_aaCode(aaCode)
84{
85 m_precision_ptr = precision_ptr;
86
87 KeyCell key_cell_init_first;
88 key_cell_init_first.n_row = 0;
89 key_cell_init_first.score = 0;
90 key_cell_init_first.beginning = 0;
91 key_cell_init_first.tree_id = 0;
92 m_interest_cells.push_back(key_cell_init_first);
93}

References pappso::specpeptidoms::KeyCell::beginning, m_aaCode, m_interest_cells, m_precision_ptr, m_scorevalues, pappso::specpeptidoms::KeyCell::n_row, pappso::specpeptidoms::KeyCell::score, and pappso::specpeptidoms::KeyCell::tree_id.

◆ ~SemiGlobalAlignment()

pappso::specpeptidoms::SemiGlobalAlignment::~SemiGlobalAlignment ( )

Destructor

Definition at line 774 of file semiglobalalignment.cpp.

775{
776}

Member Function Documentation

◆ checkSequenceDiversity()

bool pappso::specpeptidoms::SemiGlobalAlignment::checkSequenceDiversity ( const QString & sequence,
std::size_t window,
std::size_t minimum_aa_diversity )
static

check that the sequence has a minimum of amino acid checkSequenceDiversity

Parameters
sequenceprotein sequence
windowthe size of substring to check
minimum_aa_diversityminimum number of different amino acid in this window

Definition at line 1032 of file semiglobalalignment.cpp.

1035{
1036 qDebug() << "sequence=" << sequence << " window=" << window
1037 << " minimum_aa_diversity=" << minimum_aa_diversity;
1038 if(sequence.size() < window)
1039 return false;
1040 auto it_begin = sequence.begin();
1041 auto it_end = sequence.begin() + window;
1042 QString window_copy(sequence.mid(0, window));
1043 while(it_end != sequence.end())
1044 {
1045 std::partial_sort_copy(it_begin, it_end, window_copy.begin(), window_copy.end());
1046
1047 qDebug() << window_copy;
1048 std::size_t uniqueCount =
1049 std::unique(window_copy.begin(), window_copy.end()) - window_copy.begin();
1050
1051 qDebug() << uniqueCount;
1052 if(uniqueCount < minimum_aa_diversity)
1053 return false;
1054 it_begin++;
1055 it_end++;
1056 }
1057 return true;
1058}

Referenced by preciseAlign(), pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment(), and pappso::cbor::psm::PsmSpecPeptidOmsScan::storeAlignment().

◆ correctAlign()

void pappso::specpeptidoms::SemiGlobalAlignment::correctAlign ( const SpOMSProtein & protein_subseq,
const SpOMSProtein * protein_ptr,
const SpOMSSpectrum & spectrum,
std::vector< std::size_t > & peaks_to_remove,
std::size_t offset )
private

Recursively performs the correction of the alignment.

Parameters
protein_seqProtein reversed sequence to align.
protein_ptrProtein pointer on the sequence to align.
spectrumSpectrum to align.
peaks_to_removePeaks to remove from the spectrum.
offsetSize of the protein sequence minus beginning of the alignment. Used to compute the position of the alignment in the protein sequence.

Definition at line 275 of file semiglobalalignment.cpp.

280{
281 std::vector<AaPosition> aa_positions;
282 CorrectionTree correction_tree;
283 std::vector<std::size_t> final_peaks_to_remove;
284
285 KeyCell key_cell_init;
286 key_cell_init.beginning = 0;
287 key_cell_init.n_row = 0;
288 key_cell_init.score = m_scorevalues.get(ScoreType::init);
289 key_cell_init.tree_id = 0;
290
291 std::fill(m_interest_cells.begin(), m_interest_cells.end(), key_cell_init);
292
293 m_interest_cells.at(0).score = 0;
294
295 m_scenario.resetScenario();
296 qDebug();
297 for(qsizetype row_number = 1; row_number <= sequence.size(); row_number++)
298 {
299 qDebug() << row_number - 1 << " " << sequence.size();
300 qDebug() << "sequence[row_number - 1].aa" << (char)sequence[row_number - 1].aa;
301 qDebug();
302 aa_positions = spectrum.getAaPositions(sequence[row_number - 1].code, peaks_to_remove);
303 qDebug();
304 updateAlignmentMatrix(sequence, row_number, aa_positions, spectrum, false, protein_ptr);
305 qDebug();
306 }
307
308 qDebug();
309 // Correction : if complementary peaks are used, corrected spectra without one of the two peaks
310 // are generated and aligned. The best alignment is kept.
311 qDebug() << m_scenario.getBestScore();
312 if(m_scenario.getBestScore() >
313 MIN_ALIGNMENT_SCORE) // We only correct alignments with acceptable scores
314 {
315 qDebug();
316 qDebug() << sequence.getSequence();
317 qDebug() << offset;
318 qDebug() << spectrum.getPrecursorCharge();
319 saveBestAlignment(sequence, spectrum, offset);
320 qDebug();
321 for(std::size_t iter : m_best_alignment.peaks)
322 {
323 qDebug() << "iter:" << iter << "comp:" << spectrum.getComplementaryPeak(iter);
324 if(iter == spectrum.getComplementaryPeak(iter))
325 {
326 continue;
327 }
328 else if(iter > spectrum.getComplementaryPeak(iter))
329 {
330 break;
331 }
332 else if(std::find(m_best_alignment.peaks.begin(),
333 m_best_alignment.peaks.end(),
334 spectrum.getComplementaryPeak(iter)) != m_best_alignment.peaks.end())
335 {
336 correction_tree.addPeaks(iter, spectrum.getComplementaryPeak(iter));
337 }
338 }
339 std::vector<std::vector<std::size_t>> corrections = correction_tree.getPeaks();
340 if(corrections.size() > 0)
341 {
342 for(auto new_peaks_to_remove : corrections)
343 {
344 final_peaks_to_remove = std::vector<std::size_t>(new_peaks_to_remove);
345 final_peaks_to_remove.insert(
346 final_peaks_to_remove.end(), peaks_to_remove.begin(), peaks_to_remove.end());
347 correctAlign(sequence, protein_ptr, spectrum, final_peaks_to_remove, offset);
348 }
349 }
350 else if(m_scenario.getBestScore() > m_best_corrected_alignment.score)
351 {
352 m_best_corrected_alignment = m_best_alignment;
353 }
354 }
355 qDebug();
356}
pappso::specpeptidoms::SemiGlobalAlignment::updateAlignmentMatrix
void updateAlignmentMatrix(const pappso::specpeptidoms::SpOMSProtein &sequence, const std::size_t row_number, const std::vector< AaPosition > &aa_positions, const SpOMSSpectrum &spectrum, const bool fast_align, const pappso::specpeptidoms::SpOMSProtein *protein_ptr)
updates the scores of the alignment matrix for a given amino acid as well as the location heap/scenar...
Definition semiglobalalignment.cpp:385
pappso::specpeptidoms::SemiGlobalAlignment::correctAlign
void correctAlign(const SpOMSProtein &protein_subseq, const SpOMSProtein *protein_ptr, const SpOMSSpectrum &spectrum, std::vector< std::size_t > &peaks_to_remove, std::size_t offset)
Recursively performs the correction of the alignment.
Definition semiglobalalignment.cpp:275
pappso::specpeptidoms::SemiGlobalAlignment::saveBestAlignment
void saveBestAlignment(const SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, std::size_t offset)
Stores the best alignment from m_scenario in m_best_alignment.
Definition semiglobalalignment.cpp:791
pappso::specpeptidoms::MIN_ALIGNMENT_SCORE
const int MIN_ALIGNMENT_SCORE(15)

References pappso::specpeptidoms::CorrectionTree::addPeaks(), pappso::specpeptidoms::KeyCell::beginning, correctAlign(), pappso::specpeptidoms::SpOMSSpectrum::getAaPositions(), pappso::specpeptidoms::SpOMSSpectrum::getComplementaryPeak(), pappso::specpeptidoms::CorrectionTree::getPeaks(), pappso::specpeptidoms::SpOMSSpectrum::getPrecursorCharge(), pappso::specpeptidoms::SpOMSProtein::getSequence(), pappso::specpeptidoms::init, m_best_alignment, m_best_corrected_alignment, m_interest_cells, m_scenario, m_scorevalues, pappso::specpeptidoms::MIN_ALIGNMENT_SCORE(), pappso::specpeptidoms::KeyCell::n_row, saveBestAlignment(), pappso::specpeptidoms::KeyCell::score, pappso::specpeptidoms::KeyCell::tree_id, and updateAlignmentMatrix().

Referenced by correctAlign(), and preciseAlign().

◆ fastAlign()

void pappso::specpeptidoms::SemiGlobalAlignment::fastAlign ( const SpOMSSpectrum & spectrum,
const SpOMSProtein * protein_ptr )

perform the first alignment search between a protein sequence and a spectrum. The member location heap is filled with the candidates locations.

Parameters
spectrumSpectrum to align
protein_ptrProtein pointer on the sequence to align.

Definition at line 102 of file semiglobalalignment.cpp.

104{
105 std::size_t sequence_length = protein_ptr->size();
106
107 initFastAlign(spectrum);
108
109 for(std::size_t row_number = 1; row_number <= sequence_length; row_number++)
110 {
111 updateAlignmentMatrix(*protein_ptr,
112 row_number,
113 spectrum.getAaPositions(protein_ptr->at(row_number - 1).code),
114 spectrum,
115 true,
116 protein_ptr);
117 }
118}
pappso::specpeptidoms::SemiGlobalAlignment::initFastAlign
void initFastAlign(const SpOMSSpectrum &spectrum)
function made for testing the fastAlign process, initiate the variables for alignment
Definition semiglobalalignment.cpp:121

References pappso::specpeptidoms::SpOMSSpectrum::getAaPositions(), initFastAlign(), and updateAlignmentMatrix().

Referenced by pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ getBestAlignment()

const pappso::specpeptidoms::Alignment & pappso::specpeptidoms::SemiGlobalAlignment::getBestAlignment ( ) const

Returns a const ref to m_best_alignment.

Definition at line 981 of file semiglobalalignment.cpp.

982{
983 return m_best_alignment;
984}

References m_best_alignment.

Referenced by pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ getInterestCells()

const std::vector< pappso::specpeptidoms::KeyCell > & pappso::specpeptidoms::SemiGlobalAlignment::getInterestCells ( ) const

convenient function for degub purpose

Definition at line 96 of file semiglobalalignment.cpp.

97{
98 return m_interest_cells;
99}

References m_interest_cells.

Referenced by oneAlignStep().

◆ getLocationSaver()

pappso::specpeptidoms::LocationSaver pappso::specpeptidoms::SemiGlobalAlignment::getLocationSaver ( ) const

Returns a copy of m_location_saver.

Definition at line 779 of file semiglobalalignment.cpp.

780{
781 return m_location_saver;
782}

References m_location_saver.

Referenced by pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ getPotentialMassErrors()

std::vector< double > pappso::specpeptidoms::SemiGlobalAlignment::getPotentialMassErrors ( const pappso::AaCode & aa_code,
const Alignment & alignment,
const QString & protein_seq )
static

Returns a list of the potential mass errors corresponding to the provided alignment in the provided protein sequence.

Parameters
aa_codethe amino acid code of reference to get aminon acid masses
alignmentAlignment for which to get the potential mass errors.
protein_seqProtein sequence corresponding to the provided alignment.

Definition at line 987 of file semiglobalalignment.cpp.

990{
991 // qDebug() << protein_seq;
992 if(alignment.end > (std::size_t)protein_seq.size())
993 {
994 throw pappso::ExceptionOutOfRange(QString("alignment.end > protein_seq.size() %1 %2")
995 .arg(alignment.end)
996 .arg(protein_seq.size()));
997 }
998 std::vector<double> potential_mass_errors(alignment.shifts);
999 double shift = alignment.end_shift;
1000 std::size_t index;
1001 if(alignment.beginning > 0)
1002 { // -1 on unsigned int makes it wrong
1003 index = alignment.beginning - 1;
1004 while(shift > 0 && index > 0)
1005 {
1006 potential_mass_errors.push_back(shift);
1007 // qDebug() << " shift=" << shift << " index=" << index
1008 // << " letter=" << protein_seq.at(index).toLatin1();
1009 shift -= aa_code.getMass(
1010 protein_seq.at(index).toLatin1()); // Aa(protein_seq.at(index).unicode()).getMass();
1011 index--;
1012 }
1013 }
1014
1015 // qDebug() << "second";
1016 shift = alignment.begin_shift;
1017 index = alignment.end + 1;
1018 while(shift > 0 && index < (std::size_t)protein_seq.size())
1019 {
1020 potential_mass_errors.push_back(shift);
1021 qDebug() << " shift=" << shift << " index=" << index
1022 << " letter=" << protein_seq.at(index).toLatin1();
1023 shift -= aa_code.getMass(
1024 protein_seq.at(index).toLatin1()); // Aa(protein_seq.at(index).unicode()).getMass();
1025 index++;
1026 }
1027 // qDebug();
1028 return potential_mass_errors;
1029}
pappso::AaCode::getMass
double getMass(uint8_t aa_code) const
get the mass of the amino acid given its integer code the amino acid can bear some modification (if a...
Definition aacode.cpp:224

References pappso::specpeptidoms::Alignment::begin_shift, pappso::specpeptidoms::Alignment::beginning, pappso::specpeptidoms::Alignment::end, pappso::specpeptidoms::Alignment::end_shift, pappso::AaCode::getMass(), pappso::specpeptidoms::shift, and pappso::specpeptidoms::Alignment::shifts.

Referenced by pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ getScenario()

pappso::specpeptidoms::Scenario pappso::specpeptidoms::SemiGlobalAlignment::getScenario ( ) const

Returns a copy of m_scenario.

Definition at line 785 of file semiglobalalignment.cpp.

786{
787 return m_scenario;
788}

References m_scenario.

◆ initFastAlign()

void pappso::specpeptidoms::SemiGlobalAlignment::initFastAlign ( const SpOMSSpectrum & spectrum)

function made for testing the fastAlign process, initiate the variables for alignment

Definition at line 121 of file semiglobalalignment.cpp.

122{
123 // m_scenario.clear();
124 // TODO don't forget to reset any important variable
125 // m_best_alignment.reset();
126 // m_best_corrected_alignment.reset();
127 // m_best_post_processed_alignment.reset();
128
129 KeyCell key_cell_init;
130 key_cell_init.n_row = 0;
131 key_cell_init.score = m_scorevalues.get(ScoreType::init);
132 key_cell_init.beginning = 0;
133 key_cell_init.tree_id = 0;
134
135 m_interest_cells.resize(spectrum.size());
136 std::fill(m_interest_cells.begin(), m_interest_cells.end(), key_cell_init);
137
138 m_interest_cells.at(0).score = 0;
139
140 // m_location_saver.resetLocationSaver();
141 m_updated_cells.clear();
142}
pappso::specpeptidoms::SemiGlobalAlignment::m_updated_cells
std::vector< std::pair< std::size_t, KeyCell > > m_updated_cells
Definition semiglobalalignment.h:291

References pappso::specpeptidoms::KeyCell::beginning, pappso::specpeptidoms::init, m_interest_cells, m_scorevalues, m_updated_cells, pappso::specpeptidoms::KeyCell::n_row, pappso::specpeptidoms::KeyCell::score, and pappso::specpeptidoms::KeyCell::tree_id.

Referenced by fastAlign().

◆ oneAlignStep()

const std::vector< pappso::specpeptidoms::KeyCell > & pappso::specpeptidoms::SemiGlobalAlignment::oneAlignStep ( const pappso::specpeptidoms::SpOMSProtein & sequence,
const std::size_t row_number,
const std::vector< AaPosition > & aa_positions,
const SpOMSSpectrum & spectrum,
const bool fast_align,
const pappso::specpeptidoms::SpOMSProtein * protein_ptr )

function made for testing the fastAlign process, process one line and return the alignment matrix

Definition at line 1061 of file semiglobalalignment.cpp.

1068{
1069 updateAlignmentMatrix(sequence, row_number, aa_positions, spectrum, fast_align, protein_ptr);
1070 return getInterestCells();
1071}
pappso::specpeptidoms::SemiGlobalAlignment::getInterestCells
const std::vector< KeyCell > & getInterestCells() const
convenient function for degub purpose
Definition semiglobalalignment.cpp:96

References getInterestCells(), and updateAlignmentMatrix().

◆ perfectShiftPossible()

bool pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossible ( const pappso::specpeptidoms::SpOMSProtein & sequence,
const SpOMSSpectrum & spectrum,
const std::size_t origin_row,
const std::size_t current_row,
const std::size_t l_peak,
const std::size_t r_peak ) const
private

indicates if a perfect shift is possible between the provided positions

Parameters
sequenceReversed sequence of the protein being aligned
spectrumSpectrum being aligned
origin_rowbeginning row of the aa gap to verify (== index of the first missing aa in sequence)
current_rowrow being processed (== index of the current AaPosition in sequence)
l_peakleft peak index of the mz gap to verify
r_peakright peak index of the mz gap to verify

Definition at line 672 of file semiglobalalignment.cpp.

679{
680 try
681 {
682 double missing_mass = 0;
683 auto it_end = sequence.begin() + current_row;
684 for(auto iter = sequence.begin() + origin_row; (iter != it_end) && (iter != sequence.end());
685 iter++)
686 {
687 missing_mass += iter->mass; // Aa(iter->unicode()).getMass();
688 }
689 if(missing_mass > 0)
690 {
691 pappso::MzRange mz_range(missing_mass, m_precision_ptr);
692 return mz_range.contains(spectrum.getMZShift(l_peak, r_peak));
693 }
694 else
695 {
696 return false;
697 }
698 }
699 catch(const std::exception &error)
700 {
701 throw pappso::PappsoException(
702 QObject::tr("perfectShiftPossible failed std exception:\n%1").arg(error.what()));
703 }
704 catch(const pappso::PappsoException &err)
705 {
706 throw pappso::PappsoException(
707 QObject::tr("perfectShiftPossible failed :\n%1").arg(err.qwhat()));
708 }
709}
pappso::PappsoException::qwhat
virtual const QString & qwhat() const
Definition pappsoexception.h:66

References pappso::MzRange::contains(), pappso::specpeptidoms::SpOMSSpectrum::getMZShift(), m_precision_ptr, and pappso::PappsoException::qwhat().

Referenced by updateAlignmentMatrix().

◆ perfectShiftPossibleEnd()

std::size_t pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossibleEnd ( const pappso::specpeptidoms::SpOMSProtein & sequence,
const SpOMSSpectrum & spectrum,
std::size_t end_row,
std::size_t end_peak ) const
private

indicates if a perfect shift is possible between the provided positions

Parameters
sequenceReversed sequence of the protein being aligned
spectrumSpectrum being aligned
end_rowIndex of the last aligned row.
end_peakIndex of the last aligned peak.

Definition at line 739 of file semiglobalalignment.cpp.

744{
745 try
746 {
747 std::size_t perfect_shift_end = end_row + 1;
748 double missing_mass = spectrum.getMissingMass(end_peak);
749 pappso::MzRange mz_range(missing_mass, m_precision_ptr);
750 double aa_mass = 0;
751 while(aa_mass < missing_mass && perfect_shift_end < (std::size_t)sequence.size() &&
752 !mz_range.contains(aa_mass))
753 {
754 aa_mass += sequence.at(perfect_shift_end - 1)
755 .mass; // Aa(sequence.at(perfect_shift_end - 1).unicode()).getMass();
756 perfect_shift_end++;
757 }
758 if(mz_range.contains(aa_mass))
759 {
760 return perfect_shift_end - 1;
761 }
762 else
763 {
764 return end_row;
765 }
766 }
767 catch(const pappso::PappsoException &err)
768 {
769 throw pappso::PappsoException(
770 QObject::tr("perfectShiftPossibleEnd failed :\n%1").arg(err.qwhat()));
771 }
772}

References pappso::MzRange::contains(), pappso::specpeptidoms::SpOMSSpectrum::getMissingMass(), m_precision_ptr, and pappso::PappsoException::qwhat().

Referenced by saveBestAlignment().

◆ perfectShiftPossibleFrom0()

std::size_t pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossibleFrom0 ( const pappso::specpeptidoms::SpOMSProtein & sequence,
const SpOMSSpectrum & spectrum,
const std::size_t current_row,
const std::size_t r_peak ) const
private

indicates if a perfect shift is possible from the spectrum beginning to the provided peak. Returns the perfect shift origin if the shift is possible, otherwise returns the current row.

Parameters
sequenceReversed sequence of the protein being aligned
spectrumSpectrum being aligned
current_rowrow being processed (== index of the current AaPosition in sequence)
r_peakright peak index of the mz gap to verify

Definition at line 712 of file semiglobalalignment.cpp.

717{
718 std::size_t perfect_shift_origin = current_row;
719 double missing_mass = spectrum.getMZShift(0, r_peak);
720 pappso::MzRange mz_range(missing_mass, m_precision_ptr);
721 double aa_mass = 0;
722 while(aa_mass < missing_mass && perfect_shift_origin > 0 && !mz_range.contains(aa_mass))
723 {
724 aa_mass += sequence.at(perfect_shift_origin - 1)
725 .mass; // Aa(sequence.at(perfect_shift_origin - 1).unicode()).getMass();
726 perfect_shift_origin--;
727 }
728 if(mz_range.contains(aa_mass))
729 {
730 return perfect_shift_origin;
731 }
732 else
733 {
734 return current_row;
735 }
736}

References pappso::MzRange::contains(), pappso::specpeptidoms::SpOMSSpectrum::getMZShift(), and m_precision_ptr.

Referenced by updateAlignmentMatrix().

◆ postProcessingAlign()

void pappso::specpeptidoms::SemiGlobalAlignment::postProcessingAlign ( const SpOMSSpectrum & spectrum,
const SpOMSProtein * protein_ptr,
std::size_t beginning,
std::size_t length,
const std::vector< double > & shifts )

performs the post-processing : generates corrected spectra and align them

Parameters
spectrumSpectrum to align
protein_ptrProtein pointer on the sequence to align.
beginningIndex of the beginning of the subsequence to align.
lengthLength of the subsequence to align.
shiftsList of potential precursor mass errors to test.

Definition at line 359 of file semiglobalalignment.cpp.

364{
365 std::size_t current_SPC = m_best_alignment.SPC;
366 int current_best_score = m_best_alignment.score;
367 m_best_post_processed_alignment.score = 0;
368 for(double precursor_mass_error : shifts)
369 {
370 SpOMSSpectrum corrected_spectrum(spectrum, precursor_mass_error);
371 preciseAlign(corrected_spectrum, protein_ptr, beginning, length);
372 if(m_best_alignment.score >= m_best_post_processed_alignment.score)
373 {
374 m_best_post_processed_alignment = m_best_alignment;
375 }
376 }
377 if(m_best_post_processed_alignment.SPC > current_SPC &&
378 m_best_post_processed_alignment.score >= current_best_score)
379 {
380 m_best_alignment = m_best_post_processed_alignment;
381 }
382}
pappso::specpeptidoms::SemiGlobalAlignment::preciseAlign
void preciseAlign(const SpOMSSpectrum &spectrum, const SpOMSProtein *protein_ptr, const std::size_t beginning, const std::size_t length)
performs the second alignment search between a protein subsequence and a spectrum.
Definition semiglobalalignment.cpp:145

References m_best_alignment, m_best_post_processed_alignment, and preciseAlign().

Referenced by pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ preciseAlign()

void pappso::specpeptidoms::SemiGlobalAlignment::preciseAlign ( const SpOMSSpectrum & spectrum,
const SpOMSProtein * protein_ptr,
const std::size_t beginning,
const std::size_t length )

performs the second alignment search between a protein subsequence and a spectrum.

Parameters
spectrumSpectrum to align
protein_ptrProtein pointer on the sequence to align.
beginningIndex of the beginning of the subsequence to align.
lengthLength of the subsequence to align.

Definition at line 145 of file semiglobalalignment.cpp.

149{
150 try
151 {
152 qDebug();
153 const QString &protein_seq = protein_ptr->getSequence();
154 std::size_t length2;
155 if((qsizetype)(beginning + length) <= protein_seq.size())
156 {
157 length2 = length;
158 }
159 else
160 {
161 length2 = protein_seq.size() - beginning;
162 }
163
164 qDebug();
165 QString sequence_str = protein_seq.sliced(protein_seq.size() - beginning - length2, length2);
166
167 SpOMSProtein sequence("sub_sequence", sequence_str, m_aaCode);
168
169 // std::reverse(sequence.begin(), sequence.end());
170 std::vector<AaPosition> aa_positions;
171 CorrectionTree correction_tree;
172
173 qDebug();
174 m_scenario.reserve(length2 + 1, spectrum.size());
175 m_interest_cells.reserve(spectrum.size());
176 m_interest_cells.at(0).n_row = 0;
177 m_interest_cells.at(0).score = 0;
178 m_interest_cells.at(0).beginning = 0;
179 m_interest_cells.at(0).tree_id = 0;
180 for(std::size_t i = 1; i < m_interest_cells.size(); i++)
181 {
182 m_interest_cells.at(i).n_row = 0;
183 m_interest_cells.at(i).score = m_scorevalues.get(ScoreType::init);
184 m_interest_cells.at(i).beginning = 0;
185 m_interest_cells.at(i).tree_id = 0;
186 }
187 qDebug();
188 for(std::size_t iter = m_interest_cells.size(); iter < spectrum.size(); iter++)
189 {
190 m_interest_cells.push_back({0, m_scorevalues.get(ScoreType::init), 0, 0});
191 }
192 qDebug();
193 m_scenario.resetScenario();
194 qDebug();
195 for(std::size_t row_number = 1; row_number <= length2; row_number++)
196 {
197
198 qDebug() << "row_number - 1=" << row_number - 1 << " sequence.size()=" << sequence.size();
199 // aa = Aa(sequence[row_number - 1].unicode());
200 updateAlignmentMatrix(sequence,
201 row_number,
202 spectrum.getAaPositions(sequence[row_number - 1].code),
203 spectrum,
204 false,
205 protein_ptr);
206 }
207 qDebug();
208 saveBestAlignment(sequence, spectrum, protein_seq.size() - beginning);
209
210 qDebug() << m_scenario.getBestScore() << " " << MIN_ALIGNMENT_SCORE;
211 // Correction : if complementary peaks are used, corrected spectra without one of the two
212 // peaks are generated and aligned. The best alignment is kept.
213 if(m_scenario.getBestScore() >
214 MIN_ALIGNMENT_SCORE) // We only correct alignments with acceptable scores
215 {
216 // we only correct alignment if the sequence has a minimum amino acid diversity
217 if(checkSequenceDiversity(sequence.getSequence(), 5, 2))
218 {
219
220 qDebug();
221 m_best_corrected_alignment.score = 0;
222 for(std::size_t iter : m_best_alignment.peaks)
223 {
224 if(iter > spectrum.getComplementaryPeak(iter))
225 {
226 break;
227 }
228 else if(std::find(m_best_alignment.peaks.begin(),
229 m_best_alignment.peaks.end(),
230 spectrum.getComplementaryPeak(iter)) !=
231 m_best_alignment.peaks.end())
232 {
233 correction_tree.addPeaks(iter, spectrum.getComplementaryPeak(iter));
234 }
235 }
236 qDebug();
237 std::vector<std::vector<std::size_t>> corrections = correction_tree.getPeaks();
238 if(corrections.size() > 0)
239 {
240 m_best_alignment.score =
241 0; // Reset the best alignment score (we dont want to keep
242 // the original alignment if corrections are needed)
243 qDebug();
244 for(auto peaks_to_remove : corrections)
245 {
246 qDebug();
247 correctAlign(sequence,
248 protein_ptr,
249 spectrum,
250 peaks_to_remove,
251 protein_seq.size() - beginning);
252 qDebug();
253 }
254 qDebug();
255 m_best_alignment = m_best_corrected_alignment;
256 }
257 }
258 }
259 else
260 {
261 // this sequence has too much redundancy
262 // we have to lower the score
263 m_best_alignment.score = 0;
264 }
265 qDebug();
266 }
267 catch(const pappso::PappsoException &err)
268 {
269 throw pappso::PappsoException(
270 QObject::tr("SemiGlobalAlignment::preciseAlign failed :\n%1").arg(err.qwhat()));
271 }
272}
pappso::specpeptidoms::SemiGlobalAlignment::checkSequenceDiversity
static bool checkSequenceDiversity(const QString &sequence, std::size_t window, std::size_t minimum_aa_diversity)
check that the sequence has a minimum of amino acid checkSequenceDiversity
Definition semiglobalalignment.cpp:1032

References pappso::specpeptidoms::CorrectionTree::addPeaks(), checkSequenceDiversity(), correctAlign(), pappso::specpeptidoms::SpOMSSpectrum::getAaPositions(), pappso::specpeptidoms::SpOMSSpectrum::getComplementaryPeak(), pappso::specpeptidoms::CorrectionTree::getPeaks(), pappso::specpeptidoms::SpOMSProtein::getSequence(), pappso::specpeptidoms::init, m_aaCode, m_best_alignment, m_best_corrected_alignment, m_interest_cells, m_scenario, m_scorevalues, pappso::specpeptidoms::MIN_ALIGNMENT_SCORE(), pappso::PappsoException::qwhat(), saveBestAlignment(), and updateAlignmentMatrix().

Referenced by postProcessingAlign(), and pappso::cbor::psm::PsmSpecPeptidOmsScan::sequenceAlignment().

◆ saveBestAlignment()

void pappso::specpeptidoms::SemiGlobalAlignment::saveBestAlignment ( const SpOMSProtein & sequence,
const SpOMSSpectrum & spectrum,
std::size_t offset )
private

Stores the best alignment from m_scenario in m_best_alignment.

Parameters
sequencereversed sequence of the current alignment.
spectrumSpectrum currently being aligned.
offsetSize of the protein sequence minus beginning of the alignment. Used to compute the position of the alignment in the protein sequence.

Definition at line 791 of file semiglobalalignment.cpp.

795{
796 qDebug();
797 m_best_alignment.m_peptideModel.reset();
798 m_best_alignment.peaks.clear();
799 m_best_alignment.shifts.clear();
800 std::size_t previous_row; // FIXME : may be used uninitialised
801 std::size_t previous_column = 0;
802 std::size_t perfect_shift_end;
803 std::pair<std::vector<ScenarioCell>, int> best_alignment = m_scenario.getBestAlignment();
804 m_best_alignment.score = best_alignment.second;
805 std::vector<SpOMSAa> skipped_aa;
806 double skipped_mass;
807 // Retrieving beginning and end
808 if(best_alignment.first.front().previous_row > offset)
809 {
810 throw pappso::PappsoException(
811 QString("best_alignment.first.front().previous_row > offset %1 %2")
812 .arg(offset)
813 .arg(best_alignment.first.front().previous_row));
814 }
815 if(best_alignment.first.back().previous_row > offset)
816 {
817 throw pappso::PappsoException(
818 QString("best_alignment.first.back().previous_row > offset %1 %2")
819 .arg(offset)
820 .arg(best_alignment.first.back().previous_row));
821 }
822 m_best_alignment.beginning = offset - best_alignment.first.front().previous_row;
823 m_best_alignment.end = offset - best_alignment.first.back().previous_row - 1;
824
825 qDebug();
826 AminoAcidModel aa_model;
827 aa_model.m_massDifference = 0;
828 // Filling temp_interpretation and peaks vectors
829 for(auto cell : best_alignment.first)
830 {
831 switch(cell.alignment_type)
832 {
833 case AlignType::found:
834 aa_model.m_aminoAcid = sequence.at(previous_row - 1).aa;
835 aa_model.m_massDifference = 0;
836 aa_model.m_skipped = false;
837 m_best_alignment.m_peptideModel.push_back(aa_model);
838 if(previous_row > cell.previous_row + 1)
839 {
840 skipped_mass = sequence.at(previous_row - 1)
841 .mass; // Aa(sequence.at(previous_row - 1).unicode()).getMass();
842 skipped_aa =
843 sequence.sliced(cell.previous_row, previous_row - cell.previous_row - 1);
844 aa_model.m_massDifference = 0;
845 aa_model.m_skipped = true;
846 for(auto aa : skipped_aa)
847 {
848 aa_model.m_aminoAcid = aa.aa;
849 m_best_alignment.m_peptideModel.push_back(aa_model);
850 skipped_mass += aa.mass; // Aa(aa.unicode()).getMass();
851 }
852 m_best_alignment.m_peptideModel.back().m_massDifference =
853 spectrum.getMZShift(cell.previous_column, previous_column) - skipped_mass;
854 }
855 m_best_alignment.peaks.push_back(cell.previous_column);
856 break;
857 case AlignType::notFound:
858 aa_model.m_aminoAcid = sequence.at(previous_row - 1).aa;
859 aa_model.m_massDifference = 0;
860 aa_model.m_skipped = true;
861 m_best_alignment.m_peptideModel.push_back(aa_model);
862 break;
863 case AlignType::shift:
864
865 aa_model.m_aminoAcid = sequence.at(previous_row - 1).aa;
866 aa_model.m_massDifference = spectrum.getMZShift(cell.previous_column, previous_column) -
867 aa_model.m_aminoAcid.getMass();
868 aa_model.m_skipped = false;
869 m_best_alignment.m_peptideModel.push_back(aa_model);
870 m_best_alignment.peaks.push_back(cell.previous_column);
871 m_best_alignment.shifts.push_back(
872 spectrum.getMZShift(cell.previous_column, previous_column) -
873 sequence.at(previous_row - 1).mass);
874 break;
875 case AlignType::perfectShift:
876 m_best_alignment.peaks.push_back(cell.previous_column);
877 skipped_aa = sequence.sliced(cell.previous_row, previous_row - cell.previous_row);
878 std::reverse(skipped_aa.begin(), skipped_aa.end());
879 aa_model.m_massDifference = 0;
880 aa_model.m_skipped = false;
881 for(auto aa : skipped_aa)
882 {
883 aa_model.m_aminoAcid = aa.aa;
884 m_best_alignment.m_peptideModel.push_back(aa_model);
885 }
886 break;
887 case AlignType::init:
888 previous_row = cell.previous_row;
889 previous_column = cell.previous_column;
890 m_best_alignment.peaks.push_back(cell.previous_column);
891 break;
892 }
893 previous_row = cell.previous_row;
894 previous_column = cell.previous_column;
895 }
896 std::reverse(m_best_alignment.peaks.begin(), m_best_alignment.peaks.end());
897 std::reverse(m_best_alignment.m_peptideModel.begin(), m_best_alignment.m_peptideModel.end());
898
899 qDebug();
900 // Compute begin_shift and end_shift
901 MzRange zero(0, m_precision_ptr);
902 m_best_alignment.begin_shift = spectrum.getMZShift(0, m_best_alignment.peaks.front());
903 m_best_alignment.end_shift = spectrum.getMissingMass(m_best_alignment.peaks.back());
904 if(zero.contains(m_best_alignment.end_shift))
905 {
906 m_best_alignment.end_shift = 0;
907 }
908
909 qDebug();
910 // Computing SPC
911 m_best_alignment.SPC = 0;
912 for(auto peak : m_best_alignment.peaks)
913 {
914 switch(spectrum.at(peak).type)
915 {
916 case specglob::ExperimentalSpectrumDataPointType::native:
917 qDebug() << peak << "native";
918 m_best_alignment.SPC += 1;
919 break;
920 case specglob::ExperimentalSpectrumDataPointType::both:
921 qDebug() << peak << "both";
922 m_best_alignment.SPC += 2;
923 break;
924 case specglob::ExperimentalSpectrumDataPointType::synthetic:
925 qDebug() << peak << "synthetic";
926 break;
927 case specglob::ExperimentalSpectrumDataPointType::symmetric:
928 qDebug() << peak << "symmetric";
929 m_best_alignment.SPC += 1;
930 break;
931 }
932 }
933
934 qDebug();
935 // Final check of the end shift
936 if(m_best_alignment.end_shift > 0)
937 {
938 perfect_shift_end = perfectShiftPossibleEnd(sequence,
939 spectrum,
940 best_alignment.first.front().previous_row,
941 m_best_alignment.peaks.back());
942 if(perfect_shift_end != best_alignment.first.front().previous_row)
943 {
944 skipped_aa =
945 sequence.sliced(best_alignment.first.front().previous_row,
946 perfect_shift_end - best_alignment.first.front().previous_row);
947 aa_model.m_massDifference = 0;
948 aa_model.m_skipped = true;
949 for(auto aa = skipped_aa.begin(); aa != skipped_aa.end(); aa++)
950 {
951 aa_model.m_aminoAcid = aa->aa;
952 m_best_alignment.m_peptideModel.push_back(aa_model);
953 }
954 m_best_alignment.beginning = offset - perfect_shift_end;
955 m_best_alignment.end_shift = 0;
956 }
957 else
958 {
959 m_best_alignment.score += m_scorevalues.get(ScoreType::foundShift);
960 }
961 }
962
963 qDebug();
964 // Writing final interpretation
965 if(m_best_alignment.end_shift > 0)
966 {
967 m_best_alignment.m_peptideModel.setNterShift(m_best_alignment.end_shift);
968 }
969
970 std::reverse(m_best_alignment.m_peptideModel.begin(), m_best_alignment.m_peptideModel.end());
971 if(m_best_alignment.begin_shift > 0)
972 {
973 m_best_alignment.m_peptideModel.setCterShift(m_best_alignment.begin_shift);
974 }
975
976 m_best_alignment.m_peptideModel.setPrecursorMass(spectrum.getPrecursorMass());
977 qDebug();
978}
pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossibleEnd
std::size_t perfectShiftPossibleEnd(const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, std::size_t end_row, std::size_t end_peak) const
indicates if a perfect shift is possible between the provided positions
Definition semiglobalalignment.cpp:739
pappso::masschroq::PeakQualityCategory::aa
@ aa
best possible : more than one direct MS2 fragmentation in same MSRUN
Definition types.h:45
pappso::specglob::ExperimentalSpectrumDataPointType::synthetic
@ synthetic
does not correspond to existing peak, for computational purpose
Definition types.h:85
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

References pappso::specglob::both, pappso::MzRange::contains(), pappso::specpeptidoms::found, pappso::specpeptidoms::foundShift, pappso::Aa::getMass(), pappso::specpeptidoms::SpOMSSpectrum::getMissingMass(), pappso::specpeptidoms::SpOMSSpectrum::getMZShift(), pappso::specpeptidoms::SpOMSSpectrum::getPrecursorMass(), pappso::specpeptidoms::init, pappso::specpeptidoms::AminoAcidModel::m_aminoAcid, m_best_alignment, pappso::specpeptidoms::AminoAcidModel::m_massDifference, m_precision_ptr, m_scenario, m_scorevalues, pappso::specpeptidoms::AminoAcidModel::m_skipped, pappso::specglob::native, pappso::specpeptidoms::notFound, pappso::specpeptidoms::perfectShift, perfectShiftPossibleEnd(), pappso::specpeptidoms::shift, pappso::specpeptidoms::SpOMSProtein::sliced(), pappso::specglob::symmetric, and pappso::specglob::synthetic.

Referenced by correctAlign(), and preciseAlign().

◆ updateAlignmentMatrix()

void pappso::specpeptidoms::SemiGlobalAlignment::updateAlignmentMatrix ( const pappso::specpeptidoms::SpOMSProtein & sequence,
const std::size_t row_number,
const std::vector< AaPosition > & aa_positions,
const SpOMSSpectrum & spectrum,
const bool fast_align,
const pappso::specpeptidoms::SpOMSProtein * protein_ptr )
private

updates the scores of the alignment matrix for a given amino acid as well as the location heap/scenario.

Parameters
sequenceReversed sequence of the protein being aligned /!\ In case of a subsequence alignment, the provided protein must only contain the sequence to align.
row_numbernumber of the row to update (== index in sequence of the amino acid being aligned)
aa_positionslist of the AaPositions of the current amino acid
spectrumSpectrum being aligned
fast_alignWhether to use the fast version of the algorithm (for 1st alignemnt step)
protein_ptrProtein pointer on the sequence to align. /!\ In case of a subsequence alignment, the provided protein must be the protein of origin (i.e. complete sequence).

Definition at line 385 of file semiglobalalignment.cpp.

392{
393 int where = 0;
394 try
395 {
396 int score_found, score_shift, best_score, alt_score, tree_id;
397 uint32_t condition; // FIXME : may be used uninitialised
398 std::size_t best_column, shift, beginning, missing_aas, length, perfect_shift_origin;
399 KeyCell *current_cell_ptr, *tested_cell_ptr;
400 AlignType alignment_type, temp_align_type;
401
402 double smallest_aa_mass = m_aaCode.getMass((std::uint8_t)1);
403
404 m_updated_cells.reserve(aa_positions.size());
405 where = 1;
406 // Computation of the threePeaks condition, see spomsspectrum.h for more details.
407 if(fast_align)
408 {
409 condition = 3;
410 if(row_number > 1)
411 {
412 qDebug() << (char)sequence.at(row_number - 2).aa;
413 qDebug() << "condition" << condition;
414 condition += 2 << sequence.at(row_number - 2).code;
415 qDebug();
416 qDebug() << "condition" << condition;
417 }
418 }
419 where = 2;
420 for(std::vector<AaPosition>::const_iterator aa_position = aa_positions.begin();
421 aa_position != aa_positions.end();
422 aa_position++)
423 {
424
425 where = 3;
426 if(((condition & aa_position->condition) != 0) ||
427 !fast_align) // Verification of the threePeaks condition (only during first alignment).
428 {
429 current_cell_ptr = &m_interest_cells.at(aa_position->r_peak);
430 if(spectrum.peakType(aa_position->r_peak) ==
431 specglob::ExperimentalSpectrumDataPointType::both)
432 {
433 score_found = m_scorevalues.get(ScoreType::foundDouble);
434 score_shift = m_scorevalues.get(ScoreType::foundShiftDouble);
435 }
436 else
437 {
438 score_found = m_scorevalues.get(ScoreType::found);
439 score_shift = m_scorevalues.get(ScoreType::foundShift);
440 }
441
442 // not found case (always computed)
443 best_column = aa_position->r_peak;
444 best_score = current_cell_ptr->score + (row_number - current_cell_ptr->n_row) *
445 m_scorevalues.get(ScoreType::notFound);
446 beginning = current_cell_ptr->beginning;
447 tree_id = current_cell_ptr->tree_id;
448 alignment_type = AlignType::notFound;
449
450 // found case (Can only happen if the left peak is supported)
451 if(aa_position->l_support)
452 {
453 tested_cell_ptr = &m_interest_cells.at(aa_position->l_peak);
454 if(aa_position->l_peak == 0)
455 {
456 alt_score = tested_cell_ptr->score + score_found;
457 }
458 else
459 {
460 if(tested_cell_ptr->n_row == row_number - 1)
461 {
462 alt_score = tested_cell_ptr->score +
463 (row_number - tested_cell_ptr->n_row - 1) *
464 m_scorevalues.get(ScoreType::notFound) +
465 score_found;
466 }
467 else
468 {
469 alt_score = tested_cell_ptr->score +
470 (row_number - tested_cell_ptr->n_row - 1) *
471 m_scorevalues.get(ScoreType::notFound) +
472 score_shift;
473 }
474 }
475 if(alt_score >= best_score)
476 {
477 alignment_type = AlignType::found;
478 best_score = alt_score;
479 best_column = aa_position->l_peak;
480 if(best_column == 0)
481 {
482 if(row_number < ALIGNMENT_SURPLUS)
483 {
484 beginning = 0;
485 }
486 else
487 {
488 beginning = std::max((std::size_t)(row_number - ALIGNMENT_SURPLUS),
489 (std::size_t)0);
490 }
491 if(fast_align)
492 {
493 tree_id = m_location_saver.getNextTree();
494 }
495 }
496 else
497 {
498 beginning = tested_cell_ptr->beginning;
499 tree_id = tested_cell_ptr->tree_id;
500 }
501 }
502 }
503
504 where = 4;
505 // generic shift case (all shifts are tested)
506 if(aa_position->l_support)
507 {
508 shift = 1;
509 }
510 else
511 {
512 shift = 0;
513 }
514 while(shift < aa_position->l_peak)
515 {
516 tested_cell_ptr = &m_interest_cells.at(aa_position->l_peak - shift);
517 // verification saut parfait
518 if(perfectShiftPossible(sequence,
519 spectrum,
520 tested_cell_ptr->n_row,
521 row_number,
522 aa_position->l_peak - shift,
523 aa_position->r_peak) &&
524 shift <= TOL_PEAKS_MISSING)
525 {
526 alt_score = tested_cell_ptr->score +
527 (row_number - tested_cell_ptr->n_row - 1) *
528 m_scorevalues.get(ScoreType::notFound) +
529 score_found;
530 temp_align_type = AlignType::perfectShift;
531 }
532 else
533 {
534 alt_score = tested_cell_ptr->score +
535 (row_number - tested_cell_ptr->n_row - 1) *
536 m_scorevalues.get(ScoreType::notFound) +
537 score_shift;
538 temp_align_type = AlignType::shift;
539 }
540 if(alt_score > best_score)
541 {
542 alignment_type = temp_align_type;
543 best_score = alt_score;
544 best_column = aa_position->l_peak - shift;
545 beginning = tested_cell_ptr->beginning;
546 tree_id = tested_cell_ptr->tree_id;
547 }
548 shift++;
549 }
550
551 where = 5;
552 // case shift from column 0 (no penalties if all precedent amino acids are missed)
553 tested_cell_ptr = &m_interest_cells.at(0);
554 // verification saut parfait
555 if(aa_position->r_peak <= TOL_PEAKS_MISSING_FIRST_COLUMN)
556 {
557 perfect_shift_origin =
558 perfectShiftPossibleFrom0(sequence, spectrum, row_number, aa_position->r_peak);
559 }
560 else
561 {
562 perfect_shift_origin = row_number;
563 }
564
565 if(perfect_shift_origin != row_number)
566 {
567 alt_score = tested_cell_ptr->score + score_found;
568 temp_align_type = AlignType::perfectShift;
569 }
570 else
571 {
572 alt_score = tested_cell_ptr->score + score_shift;
573 temp_align_type = AlignType::shift;
574 }
575
576 where = 6;
577 if(alt_score > best_score)
578 {
579 alignment_type = temp_align_type;
580 best_score = alt_score;
581 best_column = 0;
582 missing_aas =
583 std::floor(spectrum.getMZShift(0, aa_position->l_peak) / smallest_aa_mass);
584 if(row_number < ALIGNMENT_SURPLUS + missing_aas)
585 {
586 beginning = 0;
587 }
588 else
589 {
590 beginning =
591 std::max((std::size_t)(row_number - missing_aas - ALIGNMENT_SURPLUS),
592 (std::size_t)0);
593 }
594 where = 7;
595 if(fast_align)
596 {
597 tree_id = m_location_saver.getNextTree();
598 }
599 }
600
601 where = 8;
602 if(best_column != aa_position->r_peak)
603 {
604 m_updated_cells.push_back(
605 {aa_position->r_peak, {row_number, best_score, beginning, tree_id}});
606 }
607
608 where = 9;
609 if(best_score > m_location_saver.getMinScore(tree_id) && fast_align)
610 {
611 length =
612 row_number - beginning + 1 +
613 std::ceil(spectrum.getMissingMass(aa_position->r_peak) / smallest_aa_mass) +
614 ALIGNMENT_SURPLUS;
615 where = 10;
616 m_location_saver.addLocation(beginning, length, tree_id, best_score, protein_ptr);
617 }
618 else if(!fast_align)
619 {
620
621 where = 11;
622 if(alignment_type == AlignType::perfectShift && best_column == 0)
623 {
624 m_scenario.saveOrigin(row_number,
625 aa_position->r_peak,
626 perfect_shift_origin,
627 0,
628 best_score,
629 AlignType::perfectShift);
630 }
631 else
632 {
633 m_scenario.saveOrigin(row_number,
634 aa_position->r_peak,
635 m_interest_cells.at(best_column).n_row,
636 best_column,
637 best_score,
638 alignment_type);
639 }
640 }
641 }
642 }
643
644 where = 30;
645 // Update row number in column 0
646 m_updated_cells.push_back({0, {row_number, 0, 0, 0}});
647
648 // Save updated key cells in the matrix
649 while(m_updated_cells.size() > 0)
650 {
651 qDebug() << m_interest_cells.size() << " " << m_updated_cells.back().first;
652 m_interest_cells.at(m_updated_cells.back().first) = m_updated_cells.back().second;
653 m_updated_cells.pop_back();
654 }
655 where++;
656 }
657 catch(const std::exception &error)
658 {
659 throw pappso::PappsoException(
660 QObject::tr("updateAlignmentMatrix failed std::exception :\n%1 %2")
661 .arg(where)
662 .arg(error.what()));
663 }
664 catch(const pappso::PappsoException &err)
665 {
666 throw pappso::PappsoException(
667 QObject::tr("updateAlignmentMatrix failed :\n%1").arg(err.qwhat()));
668 }
669}
pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossible
bool perfectShiftPossible(const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, const std::size_t origin_row, const std::size_t current_row, const std::size_t l_peak, const std::size_t r_peak) const
indicates if a perfect shift is possible between the provided positions
Definition semiglobalalignment.cpp:672
pappso::specpeptidoms::SemiGlobalAlignment::perfectShiftPossibleFrom0
std::size_t perfectShiftPossibleFrom0(const pappso::specpeptidoms::SpOMSProtein &sequence, const SpOMSSpectrum &spectrum, const std::size_t current_row, const std::size_t r_peak) const
indicates if a perfect shift is possible from the spectrum beginning to the provided peak....
Definition semiglobalalignment.cpp:712
pappso::specpeptidoms::ALIGNMENT_SURPLUS
const uint ALIGNMENT_SURPLUS(5)
pappso::specpeptidoms::TOL_PEAKS_MISSING
const uint TOL_PEAKS_MISSING(4)
pappso::specpeptidoms::TOL_PEAKS_MISSING_FIRST_COLUMN
const uint TOL_PEAKS_MISSING_FIRST_COLUMN(5)

References pappso::specpeptidoms::ALIGNMENT_SURPLUS(), pappso::specpeptidoms::KeyCell::beginning, pappso::specglob::both, pappso::specpeptidoms::found, pappso::specpeptidoms::foundDouble, pappso::specpeptidoms::foundShift, pappso::specpeptidoms::foundShiftDouble, pappso::specpeptidoms::SpOMSSpectrum::getMissingMass(), pappso::specpeptidoms::SpOMSSpectrum::getMZShift(), m_aaCode, m_interest_cells, m_location_saver, m_scenario, m_scorevalues, m_updated_cells, pappso::specpeptidoms::KeyCell::n_row, pappso::specpeptidoms::notFound, pappso::specpeptidoms::SpOMSSpectrum::peakType(), pappso::specpeptidoms::perfectShift, perfectShiftPossible(), perfectShiftPossibleFrom0(), pappso::PappsoException::qwhat(), pappso::specpeptidoms::KeyCell::score, pappso::specpeptidoms::shift, pappso::specpeptidoms::TOL_PEAKS_MISSING(), pappso::specpeptidoms::TOL_PEAKS_MISSING_FIRST_COLUMN(), and pappso::specpeptidoms::KeyCell::tree_id.

Referenced by correctAlign(), fastAlign(), oneAlignStep(), and preciseAlign().

Member Data Documentation

◆ m_aaCode

const AaCode& pappso::specpeptidoms::SemiGlobalAlignment::m_aaCode
private

Definition at line 295 of file semiglobalalignment.h.

Referenced by SemiGlobalAlignment(), preciseAlign(), and updateAlignmentMatrix().

◆ m_best_alignment

Alignment pappso::specpeptidoms::SemiGlobalAlignment::m_best_alignment
private

Definition at line 298 of file semiglobalalignment.h.

Referenced by correctAlign(), getBestAlignment(), postProcessingAlign(), preciseAlign(), and saveBestAlignment().

◆ m_best_corrected_alignment

Alignment pappso::specpeptidoms::SemiGlobalAlignment::m_best_corrected_alignment
private

Definition at line 299 of file semiglobalalignment.h.

Referenced by correctAlign(), and preciseAlign().

◆ m_best_post_processed_alignment

Alignment pappso::specpeptidoms::SemiGlobalAlignment::m_best_post_processed_alignment
private

Definition at line 300 of file semiglobalalignment.h.

Referenced by postProcessingAlign().

◆ m_interest_cells

std::vector<KeyCell> pappso::specpeptidoms::SemiGlobalAlignment::m_interest_cells
private

Definition at line 290 of file semiglobalalignment.h.

Referenced by SemiGlobalAlignment(), correctAlign(), getInterestCells(), initFastAlign(), preciseAlign(), and updateAlignmentMatrix().

◆ m_location_saver

LocationSaver pappso::specpeptidoms::SemiGlobalAlignment::m_location_saver
private

Definition at line 296 of file semiglobalalignment.h.

Referenced by getLocationSaver(), and updateAlignmentMatrix().

◆ m_precision_ptr

pappso::PrecisionPtr pappso::specpeptidoms::SemiGlobalAlignment::m_precision_ptr
private

Definition at line 294 of file semiglobalalignment.h.

Referenced by SemiGlobalAlignment(), perfectShiftPossible(), perfectShiftPossibleEnd(), perfectShiftPossibleFrom0(), and saveBestAlignment().

◆ m_scenario

Scenario pappso::specpeptidoms::SemiGlobalAlignment::m_scenario
private

Definition at line 297 of file semiglobalalignment.h.

Referenced by correctAlign(), getScenario(), preciseAlign(), saveBestAlignment(), and updateAlignmentMatrix().

◆ m_scorevalues

const ScoreValues& pappso::specpeptidoms::SemiGlobalAlignment::m_scorevalues
private

Definition at line 292 of file semiglobalalignment.h.

Referenced by SemiGlobalAlignment(), correctAlign(), initFastAlign(), preciseAlign(), saveBestAlignment(), and updateAlignmentMatrix().

◆ m_updated_cells

std::vector<std::pair<std::size_t, KeyCell> > pappso::specpeptidoms::SemiGlobalAlignment::m_updated_cells
private

Definition at line 291 of file semiglobalalignment.h.

Referenced by initFastAlign(), and updateAlignmentMatrix().

◆ min_score

const int pappso::specpeptidoms::SemiGlobalAlignment::min_score = 15
private

Definition at line 293 of file semiglobalalignment.h.

The documentation for this class was generated from the following files: