Loading...
Searching...
No Matches
spectree.cpp
Go to the documentation of this file.
/**
* \file pappsomspp/processing/spectree/spectree.cpp
* \date 11/12/2023
* \author Olivier Langella
* \brief Matthieu David's SpecTree structure
*
* C++ implementation of algorithm already described in :
* 1. David, M., Fertin, G., Rogniaux, H. & Tessier, D. SpecOMS: A Full Open
* Modification Search Method Performing All-to-All Spectra Comparisons within
* Minutes. J. Proteome Res. 16, 3030–3038 (2017).
*
* https://www.theses.fr/2019NANT4092
*/
/*
* SpecTree
* Copyright (C) 2023 Olivier Langella
* <olivier.langella@universite-paris-saclay.fr>
*
* This program is free software: you can redistribute ipetide to spectrum
* alignmentt 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 "spectree.h"
#include <QDebug>
#include <QObject>
#include "../../exception/exceptioninterrupted.h"
#include "../../exception/exceptionoutofrange.h"
namespace pappso
{
namespace spectree
{
SpecTree::SpecTree(const BucketClustering &bucket_clustering)
{
std::vector<Bucket> bins = bucket_clustering.asSortedList();
SpecTreeNode node;
node.value = bins[0].getCartList().at(0);
node.count = 1;
qDebug() << node.value;
// Specific processing for the first bucket : the first node and all its
// children are added in SpecTrees and the transversal accessor is updated
// accordingly
std::vector<std::size_t> spectrumLastNodeIndex;
spectrumLastNodeIndex.resize(bucket_clustering.getItemCartCount(),
index_not_defined);
m_spectrumFirstNodeIndex.resize(bucket_clustering.getItemCartCount(),
index_not_defined);
m_spectrumFirstNodeIndex[node.value] = m_nodeList.size() - 1;
m_nodeList.push_back(node);
manageSideAccess(spectrumLastNodeIndex);
for(std::size_t cpt = 1; cpt < bins[0].size(); ++cpt)
{
node.parentIndex = m_nodeList.size() - 1;
node.nextIndex = index_not_defined;
node.value = bins[0].getCartList().at(cpt);
node.count = 1;
m_nodeList.push_back(node);
manageSideAccess(spectrumLastNodeIndex);
}
qDebug();
// For the rest of the buckets, we increment the node counter in prefix
// common parts without creating new nodes. Once out of the prefix common
// part, new nodes are created similarily as they were for the first bucket.
for(std::size_t bins_index = 1; bins_index < bins.size(); ++bins_index)
{
qDebug() << "bins_index=" << bins_index;
std::size_t pos = 0;
bool common = false;
std::size_t parentNode = index_not_defined;
std::size_t previous_bin_size = bins[bins_index - 1].size();
// detection of a common first element
if(pos < previous_bin_size)
{
if(bins[bins_index].getCartList().at(pos) ==
bins[bins_index - 1].getCartList().at(pos))
{
common = true;
}
}
// insertion procedure for a common prefix part insertion : we go
// through existing nodes and increment their counter value by one
while(common)
{
qDebug() << "common pos=" << pos;
std::size_t spectrum_index = bins[bins_index].getCartList().at(pos);
qDebug() << "common spectrum_index=" << spectrum_index;
std::size_t tempNodeIndex = spectrumLastNodeIndex[spectrum_index];
qDebug() << "common tempNodeIndex=" << tempNodeIndex;
m_nodeList[tempNodeIndex].count++;
// checking on the prefix common part persistence conditions ; if
// not valid, we are not in a common prefix part anymore for the
// currently inserted bucket
if((previous_bin_size - 1 == pos) ||
(bins[bins_index].size() - 1 == pos))
{
qDebug() << "last one";
common = false;
}
else
{
qDebug() << "bins[bins_index - 1].getSpectrumIndex(pos + 1)="
<< bins[bins_index - 1].getCartList().at(pos + 1);
qDebug() << "bins[bins_index].getSpectrumIndex(pos + 1)="
<< bins[bins_index].getCartList().at(pos + 1);
if(bins[bins_index - 1].getCartList().at(pos + 1) !=
bins[bins_index].getCartList().at(pos + 1))
{
common = false;
}
}
parentNode = tempNodeIndex;
++pos;
qDebug() << "common=" << common;
}
qDebug();
// insertion procedure for all the remaining nodes not in the prefix
// common part : we create new nodes with a counter value of 1 and
// update the transversal accessor and last inserted node list
while(pos < bins[bins_index].size())
{
std::size_t spectrum_index = bins[bins_index].getCartList().at(pos);
node.parentIndex = parentNode;
node.nextIndex = index_not_defined;
node.value = spectrum_index;
node.count = 1;
m_nodeList.push_back(node);
manageSideAccess(spectrumLastNodeIndex);
parentNode = m_nodeList.size() - 1;
++pos;
}
}
}
SpecTree::~SpecTree()
{
}
void
SpecTree::addNewNode(const SpecTree::SpecTreeNode &node)
{
m_nodeList.push_back(node);
}
void
SpecTree::manageSideAccess(std::vector<std::size_t> &spectrumLastNodeIndex)
{
auto spectrum_index = m_nodeList.back().value;
auto node_position = m_nodeList.size() - 1;
if(m_spectrumFirstNodeIndex[spectrum_index] == index_not_defined)
{
// first occurence of this spectrum
m_spectrumFirstNodeIndex[spectrum_index] = node_position;
spectrumLastNodeIndex[spectrum_index] = node_position;
}
else
{
if(spectrumLastNodeIndex[spectrum_index] == node_position)
{
}
else
{
m_nodeList[spectrumLastNodeIndex[spectrum_index]].nextIndex =
node_position;
spectrumLastNodeIndex[spectrum_index] = node_position;
}
}
}
QString
SpecTree::toString() const
{
QString node_str("node|spectrum|parent|next|count|\n");
std::size_t i = 0;
for(auto node : m_nodeList)
{
int parent = -1;
if(node.parentIndex < index_not_defined)
parent = node.parentIndex;
int next = -1;
if(node.nextIndex < index_not_defined)
next = node.nextIndex;
node_str.append(QString("node_%1|%2|%3|%4|%5end\n")
.arg(i)
.arg(node.value)
.arg(parent)
.arg(next)
.arg(node.count));
i++;
}
return node_str;
}
const std::vector<std::size_t> &
SpecTree::getSpectrumFirstNodeIndex() const
{
return m_spectrumFirstNodeIndex;
}
void
SpecTree::xtract(UiMonitorInterface &monitor,
SpecXtractInterface &spec_xtract,
std::size_t minimum_count,
std::size_t cart_id_range_max,
std::size_t cart_id_range_min,
std::size_t target_cart_id_max,
std::size_t target_cart_id_min) const
{
if(cart_id_range_max < cart_id_range_min)
{
std::swap(cart_id_range_max, cart_id_range_min);
}
if(cart_id_range_max < m_spectrumFirstNodeIndex.size())
{
monitor.setStatus(QObject::tr("starting SpecXtract operation"));
MapSimilarityCount map_count;
map_count.map_id_count.resize(cart_id_range_max);
// initialisation of the structure to store the compared pairs
// int[][] results = new int[m_count][3];
// int count = 0;
// iteration on the transversal accessor starting from the deepest spetrum
// and climbing up to limit
std::size_t spectra1 = cart_id_range_max;
monitor.setTotalSteps(cart_id_range_max - cart_id_range_min);
while(true)
{
if(monitor.shouldIstop())
{
throw pappso::ExceptionInterrupted(
QObject::tr("SpecXtract job stopped"));
}
qDebug() << "spectra1=" << spectra1;
spec_xtract.beginItemCartExtraction(spectra1);
monitor.count();
// System.out.println(start);
map_count.keys.clear();
map_count.aboveThreshold.clear();
extractSpectrumSimilarityCount(map_count,
minimum_count,
spectra1,
target_cart_id_max,
target_cart_id_min);
qDebug() << "spectra1=" << spectra1
<< " map_count.aboveThreshold.size()="
<< map_count.aboveThreshold.size()
<< " map_count.keys.size()=" << map_count.keys.size();
for(std::size_t spectra2 : map_count.aboveThreshold)
{
// std::size_t spectra2 = *it;
qDebug() << "spectra2=" << spectra2;
// std::size_t total_extracted_count =
// map_count.map_id_count.at(spectra2);
spec_xtract.reportSimilarity(
spectra1, spectra2, map_count.map_id_count.at(spectra2).count);
}
spec_xtract.endItemCartExtraction(spectra1);
// monitor.count();
if(spectra1 == cart_id_range_min)
break;
spectra1--;
}
}
else
{
throw pappso::ExceptionOutOfRange(
QObject::tr(
"item cart index %1 out of range (item cart id max size=%2)")
.arg(cart_id_range_max)
.arg(m_spectrumFirstNodeIndex.size()));
}
qDebug();
}
void
SpecTree::walkBackInBranchFromNode(const SpecTree::SpecTreeNode &start_node,
MapSimilarityCount &map_count,
std::size_t minimum_count,
std::size_t target_cart_id_max,
std::size_t target_cart_id_min) const
{
std::size_t parent_id = start_node.parentIndex;
qDebug() << " start node spectrum=" << start_node.value
<< " parent_id=" << parent_id << " init_count=" << start_node.count
<< " minimum_count=" << minimum_count
<< " target_cart_id_min=" << target_cart_id_min;
while(parent_id != index_not_defined)
{
// next_node :
const SpecTree::SpecTreeNode &current_node = m_nodeList[parent_id];
if(current_node.value < target_cart_id_min)
break;
if(current_node.value <= target_cart_id_max)
{
auto &map_id_count = map_count.map_id_count[current_node.value];
if(map_id_count.lastWitness == start_node.value)
{
map_id_count.count += start_node.count;
// if at one time, a similarity exceeds the threshold value for
// the first time, it is written into the proper stack
if((map_id_count.count >= minimum_count) &&
(map_id_count.aboveThreshold == false))
{
map_id_count.aboveThreshold = true;
map_count.aboveThreshold.push_back(current_node.value);
}
}
else
{
map_id_count.lastWitness = start_node.value;
map_id_count.count = start_node.count;
map_id_count.aboveThreshold = false;
// if at one time, a similarity exceeds the threshold value for
// the first time, it is written into the proper stack
if(map_id_count.count >= minimum_count)
{
map_id_count.aboveThreshold = true;
map_count.aboveThreshold.push_back(current_node.value);
}
map_count.keys.push_back(current_node.value);
}
}
parent_id = current_node.parentIndex;
}
}
std::size_t
SpecTree::walkBackInBranchFromNodeToTarget(
const SpecTree::SpecTreeNode &start_node,
std::size_t spectrum_index_target) const
{
std::size_t parent_id = start_node.parentIndex;
while(parent_id != index_not_defined)
{
// next_node :
const SpecTree::SpecTreeNode &current_node = m_nodeList[parent_id];
if(current_node.value < spectrum_index_target)
return 0;
if(current_node.value == spectrum_index_target)
return start_node.count;
parent_id = current_node.parentIndex;
}
return 0;
}
/**
* Extraction of the similarities above a given threshold between a given
* spectrum and all other spectra higher located in SpecTrees.
* @param tree_indices SpecTrees the data has to be extracted from
* @param spectrum_index The spectrum with ated in SpecTrees.
* @param tree_indices SpecTrees the data has to be extracted from
* @param spectrum_index The spectrum with which we want to extract similarities
* @param minimum_count The threshold to use for the extraction
* @since 0.1
*/
void
SpecTree::extractSpectrumSimilarityCount(MapSimilarityCount &map_count,
std::size_t minimum_count,
std::size_t spectrum_index,
std::size_t target_cart_id_max,
std::size_t target_cart_id_min) const
{
qDebug() << " spectrum_index=" << spectrum_index;
// sOne;
// we position the current node at the first occurence of the spectrum in
// the transversal accessor and initialise variablethresholds accordingly
std::size_t node_index = m_spectrumFirstNodeIndex[spectrum_index];
while(node_index != index_not_defined)
{
qDebug() << "node_index=" << node_index;
const SpecTreeNode &current_node = m_nodeList[node_index];
// go back in the tree branch from the current_node
walkBackInBranchFromNode(current_node,
map_count,
minimum_count,
target_cart_id_max,
target_cart_id_min);
node_index = current_node.nextIndex;
}
qDebug();
}
std::size_t
SpecTree::extractSpectrumPairSimilarityCount(std::size_t spectrum_a_index,
std::size_t spectrum_b_index) const
{
if(spectrum_a_index > spectrum_b_index)
std::swap(spectrum_a_index, spectrum_b_index);
std::size_t count_target = 0;
if(spectrum_b_index < m_spectrumFirstNodeIndex.size())
{
std::size_t node_index = m_spectrumFirstNodeIndex[spectrum_b_index];
while(node_index != index_not_defined)
{
qDebug() << "node_index=" << node_index;
const SpecTreeNode &current_node = m_nodeList[node_index];
// go back in the tree branch from the current_node
count_target +=
walkBackInBranchFromNodeToTarget(current_node, spectrum_a_index);
node_index = current_node.nextIndex;
}
}
else
{
throw pappso::ExceptionOutOfRange(
QObject::tr(
"spectrum_index %1 out of range (spectrum_index max size=%2)")
.arg(spectrum_b_index)
.arg(m_spectrumFirstNodeIndex.size()));
}
return count_target;
}
} // namespace spectree
} // namespace pappso
pappso::ExceptionInterrupted
Definition exceptioninterrupted.h:36
pappso::ExceptionOutOfRange
Definition exceptionoutofrange.h:32
pappso::UiMonitorInterface
Definition uimonitorinterface.h:50
pappso::UiMonitorInterface::setStatus
virtual void setStatus(const QString &status)=0
current status of the process
pappso::UiMonitorInterface::setTotalSteps
virtual void setTotalSteps(std::size_t total_number_of_steps)
use it if the number of steps is known in an algorithm the total number of steps is usefull to report...
Definition uimonitorinterface.h:66
pappso::UiMonitorInterface::shouldIstop
virtual bool shouldIstop()=0
should the procces be stopped ? If true, then cancel process Use this function at strategic point of ...
pappso::UiMonitorInterface::count
virtual void count()=0
count steps report when a step is computed in an algorithm
pappso::spectree::BucketClustering
Definition bucketclustering.h:56
pappso::spectree::BucketClustering::asSortedList
std::vector< Bucket > asSortedList() const
Definition bucketclustering.cpp:102
pappso::spectree::BucketClustering::getItemCartCount
std::size_t getItemCartCount() const
Definition bucketclustering.cpp:72
pappso::spectree::SpecTree::m_nodeList
std::vector< SpecTreeNode > m_nodeList
Definition spectree.h:172
pappso::spectree::SpecTree::walkBackInBranchFromNodeToTarget
std::size_t walkBackInBranchFromNodeToTarget(const SpecTree::SpecTreeNode &start_node, std::size_t spectrum_index_target) const
Definition spectree.cpp:374
pappso::spectree::SpecTree::getSpectrumFirstNodeIndex
const std::vector< std::size_t > & getSpectrumFirstNodeIndex() const
get the adress map of sepctrum index and their first node index
Definition spectree.cpp:224
pappso::spectree::SpecTree::~SpecTree
virtual ~SpecTree()
Definition spectree.cpp:163
pappso::spectree::SpecTree::index_not_defined
static constexpr std::size_t index_not_defined
Definition spectree.h:118
pappso::spectree::SpecTree::walkBackInBranchFromNode
void walkBackInBranchFromNode(const SpecTree::SpecTreeNode &start_node, MapSimilarityCount &map_count, std::size_t minimum_count, std::size_t target_cart_id_max, std::size_t target_cart_id_min) const
Definition spectree.cpp:318
pappso::spectree::SpecTree::m_spectrumFirstNodeIndex
std::vector< std::size_t > m_spectrumFirstNodeIndex
Definition spectree.h:173
pappso::spectree::SpecTree::toString
QString toString() const
Definition spectree.cpp:199
pappso::spectree::SpecTree::extractSpectrumPairSimilarityCount
std::size_t extractSpectrumPairSimilarityCount(std::size_t spectrum_a_index, std::size_t spectrum_b_index) const
get the number of common component for a pair of spectrum
Definition spectree.cpp:435
pappso::spectree::SpecTree::manageSideAccess
void manageSideAccess(std::vector< std::size_t > &spectrumLastNodeIndex)
Definition spectree.cpp:174
pappso::spectree::SpecTree::addNewNode
void addNewNode(const SpecTreeNode &node)
Definition spectree.cpp:168
pappso::spectree::SpecTree::extractSpectrumSimilarityCount
void extractSpectrumSimilarityCount(MapSimilarityCount &map_count, std::size_t minimum_count, std::size_t spectrum_index, std::size_t target_cart_id_max, std::size_t target_cart_id_min) const
get a map of similarities for a given spectrum index
Definition spectree.cpp:407
pappso::spectree::SpecTree::SpecTree
SpecTree(const BucketClustering &bucket_clustering)
Definition spectree.cpp:48
pappso::spectree::SpecTree::xtract
void xtract(UiMonitorInterface &monitor, SpecXtractInterface &spec_xtract, std::size_t minimum_count, std::size_t cart_id_range_max, std::size_t cart_id_range_min, std::size_t target_cart_id_max, std::size_t target_cart_id_min) const
Definition spectree.cpp:230
pappso::spectree::SpecXtractInterface
yield similarities between pairs of ItemCart
Definition specxtractinterface.h:53
pappso::spectree::SpecXtractInterface::beginItemCartExtraction
virtual void beginItemCartExtraction(std::size_t cart_id_a)
Definition specxtractinterface.cpp:59
pappso::spectree::SpecXtractInterface::reportSimilarity
virtual void reportSimilarity(std::size_t cart_id_a, std::size_t cart_id_b, std::size_t similarity)=0
pappso::spectree::SpecXtractInterface::endItemCartExtraction
virtual void endItemCartExtraction(std::size_t cart_id_a)
Definition specxtractinterface.cpp:54
pappso
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition aa.cpp:39
spectree.h
Matthieu David's SpecTree structure.
pappso::spectree::SpecTree::MapSimilarityCount
Definition spectree.h:135
pappso::spectree::SpecTree::MapSimilarityCount::map_id_count
std::vector< MapSimilarityCountElement > map_id_count
Definition spectree.h:138
pappso::spectree::SpecTree::MapSimilarityCount::aboveThreshold
std::vector< std::size_t > aboveThreshold
Definition spectree.h:137
pappso::spectree::SpecTree::MapSimilarityCount::keys
std::vector< std::size_t > keys
Definition spectree.h:136
pappso::spectree::SpecTree::SpecTreeNode
Definition spectree.h:121
pappso::spectree::SpecTree::SpecTreeNode::count
std::size_t count
Definition spectree.h:125
pappso::spectree::SpecTree::SpecTreeNode::parentIndex
std::size_t parentIndex
Definition spectree.h:122
pappso::spectree::SpecTree::SpecTreeNode::nextIndex
std::size_t nextIndex
Definition spectree.h:123
pappso::spectree::SpecTree::SpecTreeNode::value
std::size_t value
Definition spectree.h:124