cosinesimilarity.cpp

Go to the documentation of this file.

/**
 * \file pappsomspp/processing/compartraces/cosinesimilarity.cpp
 * \date 12/06/2023
 * \author Olivier Langella
 * \brief computes cosine similarity of 2 traces vector
 * https://en.wikipedia.org/wiki/Cosine_similarity
 */
 
/*******************************************************************************
 * Copyright (c) 2023 Olivier Langella <Olivier.Langella@u-psud.fr>.
 *
 * This file is part of PAPPSOms-tools.
 *
 *     PAPPSOms-tools 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.
 *
 *     PAPPSOms-tools 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 PAPPSOms-tools.  If not, see <http://www.gnu.org/licenses/>.
 *
 ******************************************************************************/
 
 
#include "cosinesimilarity.h"
#include <cmath>
 
 
using namespace pappso;
 
CosineSimilarity::CosineSimilarity(PrecisionPtr precision)
{
  mp_precision = precision;
 
  msp_filterExclusion = std::make_shared<FilterMzExclusion>(precision);
}
 
CosineSimilarity::CosineSimilarity(const CosineSimilarity &other)
{
  mp_precision        = other.mp_precision;
  msp_filterExclusion = other.msp_filterExclusion;
}
 
CosineSimilarity::~CosineSimilarity()
{
}
 
double
pappso::CosineSimilarity::similarity(pappso::Trace trace_a,
                                     pappso::Trace trace_b) const
{
  // get inner peaks
  // we need a filter to get pairs of common peaks between a & b
  //  or at least compute quickly intensity square of common pairs
 
  msp_filterExclusion.get()->filter(trace_a);
  msp_filterExclusion.get()->filter(trace_b);
 
  auto itb                       = trace_b.begin();
  double inner_intensity_product = 0;
  for(const auto &peak_a : trace_a)
    {
      MzRange range(peak_a.x, mp_precision);
      double low = range.lower();
      double up  = range.upper();
 
      while((itb != trace_b.end()) && (itb->x < low))
        {
          itb++;
        }
      if(itb->x < up)
        {
          inner_intensity_product += peak_a.y * itb->y;
        }
    }
 
  double tracea_product = 0;
  // a intensity sum of intensity square
  for(const auto &peak_a : trace_a)
    {
      tracea_product += peak_a.y * peak_a.y;
    }
 
 
  double traceb_product = 0;
  // a intensity sum of intensity square
  for(const auto &peak_a : trace_b)
    {
      traceb_product += peak_a.y * peak_a.y;
    }
 
  return (inner_intensity_product /
          (sqrt(tracea_product) * sqrt(traceb_product)));
}