Repository: matchms/matchms
Pull request #329 (2022-02-24 20:34:27)

add NeutralLossesCosine score
- added neutral losses score (see discussion #328 ).
- refactored `test_modified_cosine.py`
- added `test_neutral_losses_cosine.py`

@hechth just fyi to avoid double-work: I here also refactored part of the tests for `ModifiedCosine` according to what you did before for `CosineGreedy`.


Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 0c7b996c..54deaa76 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -7,6 +7,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Added
+
+- new spectral similarity score: `NeutralLossesCosine` which is based on matches between neutral losses of two spectra [#329](https://github.com/matchms/matchms/pull/329)
+
 ## [0.14.0] - 2022-02-18
 
 This is the first of a few releases to work our way towards matchms 1.0.0, which also means that a few things in the API will likely change. Here the main change is that `Spectrum.metadata` is no longer a simple Python dictionary but became a `Metadata` object. In this context metadata field-names/keys will now be harmonized by default (e.g. "Precursor Mass" will become "precursor_mz). For list of conversions see [matchms key conversion table](https://github.com/matchms/matchms/blob/development/matchms/data/known_key_conversions.csv).
diff --git a/matchms/similarity/ModifiedCosine.py b/matchms/similarity/ModifiedCosine.py
index 5b572469..ff0a42dd 100644
--- a/matchms/similarity/ModifiedCosine.py
+++ b/matchms/similarity/ModifiedCosine.py
@@ -1,6 +1,6 @@
 import logging
 from typing import Tuple
-import numpy
+import numpy as np
 from matchms.filtering.add_precursor_mz import _convert_precursor_mz
 from matchms.typing import SpectrumType
 from .BaseSimilarity import BaseSimilarity
@@ -55,7 +55,7 @@ class ModifiedCosine(BaseSimilarity):
     # Set key characteristics as class attributes
     is_commutative = True
     # Set output data type, e.g. ("score", "float") or [("score", "float"), ("matches", "int")]
-    score_datatype = [("score", numpy.float64), ("matches", "int")]
+    score_datatype = [("score", np.float64), ("matches", "int")]
 
     def __init__(self, tolerance: float = 0.1, mz_power: float = 0.0,
                  intensity_power: float = 1.0):
@@ -120,19 +120,19 @@ class ModifiedCosine(BaseSimilarity):
                                                intensity_power=self.intensity_power)
 
             if zero_pairs is None:
-                zero_pairs = numpy.zeros((0, 3))
+                zero_pairs = np.zeros((0, 3))
             if nonzero_pairs is None:
-                nonzero_pairs = numpy.zeros((0, 3))
-            matching_pairs = numpy.concatenate((zero_pairs, nonzero_pairs), axis=0)
+                nonzero_pairs = np.zeros((0, 3))
+            matching_pairs = np.concatenate((zero_pairs, nonzero_pairs), axis=0)
             if matching_pairs.shape[0] > 0:
-                matching_pairs = matching_pairs[numpy.argsort(matching_pairs[:, 2])[::-1], :]
+                matching_pairs = matching_pairs[np.argsort(matching_pairs[:, 2])[::-1], :]
             return matching_pairs
 
         spec1 = reference.peaks.to_numpy
         spec2 = query.peaks.to_numpy
         matching_pairs = get_matching_pairs()
         if matching_pairs.shape[0] == 0:
-            return numpy.asarray((float(0), 0), dtype=self.score_datatype)
+            return np.asarray((float(0), 0), dtype=self.score_datatype)
         score = score_best_matches(matching_pairs, spec1, spec2,
                                    self.mz_power, self.intensity_power)
-        return numpy.asarray(score, dtype=self.score_datatype)
+        return np.asarray(score, dtype=self.score_datatype)
diff --git a/matchms/similarity/NeutralLossesCosine.py b/matchms/similarity/NeutralLossesCosine.py
new file mode 100644
index 00000000..be26a668
--- /dev/null
+++ b/matchms/similarity/NeutralLossesCosine.py
@@ -0,0 +1,108 @@
+import logging
+from typing import Tuple
+import numpy as np
+from matchms.filtering.add_precursor_mz import _convert_precursor_mz
+from matchms.typing import SpectrumType
+from .BaseSimilarity import BaseSimilarity
+from .spectrum_similarity_functions import (collect_peak_pairs,
+                                            score_best_matches)
+
+
+logger = logging.getLogger("matchms")
+
+
+class NeutralLossesCosine(BaseSimilarity):
+    """Calculate 'neutral losses cosine score' between mass spectra.
+
+    The neutral losses cosine score aims at quantifying the similarity between two
+    mass spectra. The score is calculated by finding best possible matches between
+    peaks of two spectra. Two peaks are considered a potential match if their
+    m/z ratios lie within the given 'tolerance' once a mass-shift is applied.
+    The mass shift is the difference in precursor-m/z between the two spectra.
+
+    """
+    # Set key characteristics as class attributes
+    is_commutative = True
+    # Set output data type, e.g. ("score", "float") or [("score", "float"), ("matches", "int")]
+    score_datatype = [("score", np.float64), ("matches", "int")]
+
+    def __init__(self, tolerance: float = 0.1, mz_power: float = 0.0,
+                 intensity_power: float = 1.0, ignore_peaks_above_precursor: bool = True):
+        """
+        Parameters
+        ----------
+        tolerance:
+            Peaks will be considered a match when <= tolerance apart. Default is 0.1.
+        mz_power:
+            The power to raise mz to in the cosine function. The default is 0, in which
+            case the peak intensity products will not depend on the m/z ratios.
+        intensity_power:
+            The power to raise intensity to in the cosine function. The default is 1.
+        ignore_peaks_above_precursor:
+            By default this is set to True, meaning that peaks with m/z values larger
+            than the precursor-m/z will be ignored (since those would correspond to negative
+            "neutral losses").
+        """
+        self.tolerance = tolerance
+        self.mz_power = mz_power
+        self.intensity_power = intensity_power
+        self.ignore_peaks_above_precursor = ignore_peaks_above_precursor
+
+    def pair(self, reference: SpectrumType, query: SpectrumType) -> Tuple[float, int]:
+        """Calculate neutral losses cosine score between two spectra.
+
+        Parameters
+        ----------
+        reference
+            Single reference spectrum.
+        query
+            Single query spectrum.
+
+        Returns
+        -------
+
+        Tuple with cosine score and number of matched peaks.
+        """
+        def get_valid_precursor_mz(spectrum):
+            """Extract valid precursor_mz from spectrum if possible. If not raise exception."""
+            message_precursor_missing = \
+                "Precursor_mz missing. Apply 'add_precursor_mz' filter first."
+            message_precursor_no_number = \
+                "Precursor_mz must be of type int or float. Apply 'add_precursor_mz' filter first."
+            message_precursor_below_0 = "Expect precursor to be positive number." \
+                                        "Apply 'require_precursor_mz' first"
+
+            precursor_mz = spectrum.get("precursor_mz", None)
+            if not isinstance(precursor_mz, (int, float)):
+                logger.warning(message_precursor_no_number)
+            precursor_mz = _convert_precursor_mz(precursor_mz)
+            assert precursor_mz is not None, message_precursor_missing
+            assert precursor_mz > 0, message_precursor_below_0
+            return precursor_mz
+
+        def get_matching_pairs():
+            """Find all pairs of peaks that match within the given tolerance."""
+            matching_pairs = collect_peak_pairs(spec1, spec2, self.tolerance, shift=mass_shift,
+                                                mz_power=self.mz_power,
+                                                intensity_power=self.intensity_power)
+            if matching_pairs is None:
+                return None
+            if matching_pairs.shape[0] > 0:
+                matching_pairs = matching_pairs[np.argsort(matching_pairs[:, 2])[::-1], :]
+            return matching_pairs
+
+        precursor_mz_ref = get_valid_precursor_mz(reference)
+        precursor_mz_query = get_valid_precursor_mz(query)
+        mass_shift = precursor_mz_ref - precursor_mz_query
+
+        spec1 = reference.peaks.to_numpy
+        spec2 = query.peaks.to_numpy
+        if self.ignore_peaks_above_precursor:
+            spec1 = spec1[np.where(spec1[:, 0] < precursor_mz_ref)]
+            spec2 = spec2[np.where(spec2[:, 0] < precursor_mz_query)]
+        matching_pairs = get_matching_pairs()
+        if matching_pairs is None:
+            return np.asarray((float(0), 0), dtype=self.score_datatype)
+        score = score_best_matches(matching_pairs, spec1, spec2,
+                                   self.mz_power, self.intensity_power)
+        return np.asarray(score, dtype=self.score_datatype)
diff --git a/matchms/similarity/__init__.py b/matchms/similarity/__init__.py
index 427d16a7..1310b359 100644
--- a/matchms/similarity/__init__.py
+++ b/matchms/similarity/__init__.py
@@ -24,6 +24,7 @@ from .FingerprintSimilarity import FingerprintSimilarity
 from .IntersectMz import IntersectMz
 from .MetadataMatch import MetadataMatch
 from .ModifiedCosine import ModifiedCosine
+from .NeutralLossesCosine import NeutralLossesCosine
 from .ParentMassMatch import ParentMassMatch
 from .PrecursorMzMatch import PrecursorMzMatch
 
@@ -35,6 +36,7 @@ __all__ = [
     "IntersectMz",
     "MetadataMatch",
     "ModifiedCosine",
+    "NeutralLossesCosine",
     "ParentMassMatch",
     "PrecursorMzMatch",
 ]


Test patch:
diff --git a/tests/test_modified_cosine.py b/tests/test_modified_cosine.py
index 4b9704e7..46439272 100644
--- a/tests/test_modified_cosine.py
+++ b/tests/test_modified_cosine.py
@@ -1,17 +1,27 @@
-import numpy
+import numpy as np
 import pytest
 from matchms import Spectrum
 from matchms.filtering import normalize_intensities
 from matchms.similarity import ModifiedCosine
+from .builder_Spectrum import SpectrumBuilder
+
+
+def compute_expected_score(spectrum_1, spectrum_2, matches):
+    spec1 = spectrum_1.peaks.intensities
+    spec2 = spectrum_2.peaks.intensities
+    peak_pairs_multiplied = 0
+    for match in matches:
+        peak_pairs_multiplied += spec1[match[0]] * spec2[match[1]]
+    return peak_pairs_multiplied / np.sqrt(np.sum(spec1 ** 2) * np.sum(spec2 ** 2))
 
 
 def test_modified_cosine_without_precursor_mz():
     """Test without precursor-m/z. Should raise assertion error."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"))
-
-    spectrum_2 = Spectrum(mz=numpy.array([100, 140, 190, 300, 490, 510, 1090], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"))
+    mz = np.array([100, 150, 200, 300, 500, 510, 1100], dtype="float")
+    intensities = np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float")
+    builder = SpectrumBuilder()
+    spectrum_1 = builder.with_mz(mz).with_intensities(intensities).build()
+    spectrum_2 = builder.with_mz(mz).with_intensities(intensities).build()
 
     norm_spectrum_1 = normalize_intensities(spectrum_1)
     norm_spectrum_2 = normalize_intensities(spectrum_2)
@@ -24,74 +34,60 @@ def test_modified_cosine_without_precursor_mz():
     assert str(msg.value) == expected_message
 
 
-def test_modified_cosine_with_mass_shift_5():
-    """Test modified cosine on two spectra with mass set shift."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
-                          metadata={"precursor_mz": 1000.0})
-
-    spectrum_2 = Spectrum(mz=numpy.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
-                          metadata={"precursor_mz": 1005.0})
+@pytest.mark.parametrize("peaks, tolerance, masses, expected_matches", [
+    [
+        [
+            [np.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"), np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float")],
+            [np.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"), np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float")]
+        ],
+        None, (1000.0, 1005.0), [(0, 1), (2, 2)]
+    ], [
+        [
+            [np.array([100, 299, 300, 301, 500, 510], dtype="float"), np.array([10, 500, 100, 200, 20, 100], dtype="float")],
+            [np.array([105, 305, 306, 505, 517], dtype="float"), np.array([10, 500, 100, 20, 100], dtype="float")],
+        ],
+        2.0, (1000.0, 1005.0), [(0, 0), (1, 1), (3, 2), (4, 3), (5, 4)]
+    ], [
+        [
+            [np.array([100, 110, 200, 300, 400, 500, 600], dtype="float"), np.array([100, 50, 1, 80, 1, 1, 50], dtype="float")],
+            [np.array([110, 200, 300, 310, 700, 800], dtype="float"), np.array([100, 1, 90, 90, 1, 100], dtype="float")],
+        ],
+        None, (1000.0, 1010.0), [(0, 0), (2, 1), (3, 3)]
+    ], [
+        [
+            [np.array([100, 200, 300], dtype="float"), np.array([10, 10, 500], dtype="float")],
+            [np.array([120, 220, 320], dtype="float"), np.array([10, 10, 500], dtype="float")],
+        ],
+        None, (1000.0, 1010.0), []
+    ]
+])
+def test_modified_cosine_with_mass_shift(peaks, tolerance, masses, expected_matches):
+    """Test modified cosine on two spectra with mass shift."""
+    builder = SpectrumBuilder()
+    spectrum_1 = builder.with_mz(peaks[0][0]).with_intensities(peaks[0][1]).with_metadata(metadata={"precursor_mz": masses[0]}).build()
+    spectrum_2 = builder.with_mz(peaks[1][0]).with_intensities(peaks[1][1]).with_metadata(metadata={"precursor_mz": masses[1]}).build()
 
     norm_spectrum_1 = normalize_intensities(spectrum_1)
     norm_spectrum_2 = normalize_intensities(spectrum_2)
-    modified_cosine = ModifiedCosine()
-    score = modified_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+    if tolerance is None:
+        modified_cosine = ModifiedCosine()
+    else:
+        modified_cosine = ModifiedCosine(tolerance=tolerance)
 
-    assert score["score"] == pytest.approx(0.081966, 0.0001), "Expected different cosine score."
-    assert score["matches"] == 2, "Expected 2 matching peaks."
-
-
-def test_modified_cosine_with_mass_shift_5_tolerance_2():
-    """Test modified cosine on two spectra with mass set shift and tolerance."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 200, 299, 300, 301, 500, 510], dtype="float"),
-                          intensities=numpy.array([10, 10, 500, 100, 200, 20, 100], dtype="float"),
-                          metadata={"precursor_mz": 1000.0})
-
-    spectrum_2 = Spectrum(mz=numpy.array([105, 205, 305, 306, 505, 517], dtype="float"),
-                          intensities=numpy.array([10, 10, 500, 100, 20, 100], dtype="float"),
-                          metadata={"precursor_mz": 1005})
-
-    norm_spectrum_1 = normalize_intensities(spectrum_1)
-    norm_spectrum_2 = normalize_intensities(spectrum_2)
-    modified_cosine = ModifiedCosine(tolerance=2.0)
     score = modified_cosine.pair(norm_spectrum_1, norm_spectrum_2)
-
-    assert score["score"] == pytest.approx(0.96788, 0.0001), "Expected different modified cosine score."
-    assert score["matches"] == 6, "Expected 6 matching peaks."
-
-
-def test_modified_cosine_with_mass_shifted_and_unshifted_matches():
-    """Test modified cosine on two spectra with mass set shift.
-    In this example 5 peak pairs are possible, but only 3 should be selected (every peak
-    can only be counted once!)"""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 110, 200, 300, 400, 500, 600], dtype="float"),
-                          intensities=numpy.array([100, 50, 1, 80, 1, 1, 50], dtype="float"),
-                          metadata={"precursor_mz": 1000.0})
-
-    spectrum_2 = Spectrum(mz=numpy.array([110, 200, 300, 310, 700, 800], dtype="float"),
-                          intensities=numpy.array([100, 1, 90, 90, 1, 100], dtype="float"),
-                          metadata={"precursor_mz": 1010.0})
-
-    modified_cosine = ModifiedCosine()
-    score = modified_cosine.pair(spectrum_1, spectrum_2)
-    spec1 = spectrum_1.peaks.intensities
-    spec2 = spectrum_2.peaks.intensities
-    peak_pairs_multiplied = spec1[0] * spec2[0] + spec1[3] * spec2[3] + spec1[2] * spec2[1]
-    expected_score = peak_pairs_multiplied / numpy.sqrt(numpy.sum(spec1 ** 2) * numpy.sum(spec2 ** 2))
-    assert score["score"] == pytest.approx(expected_score, 0.00001), "Expected different cosine score."
-    assert score["matches"] == 3, "Expected 3 matching peaks."
+    expected_score = compute_expected_score(norm_spectrum_1, norm_spectrum_2, expected_matches)
+    assert score["score"] == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
+    assert score["matches"] == len(expected_matches), "Expected differnt number of matching peaks."
 
 
 def test_modified_cosine_order_of_input_spectrums():
     """Test modified cosine on two spectra in changing order."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
+    spectrum_1 = Spectrum(mz=np.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"),
+                          intensities=np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
                           metadata={"precursor_mz": 1000.0})
 
-    spectrum_2 = Spectrum(mz=numpy.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"),
-                          intensities=numpy.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
+    spectrum_2 = Spectrum(mz=np.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"),
+                          intensities=np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
                           metadata={"precursor_mz": 1005.0})
 
     norm_spectrum_1 = normalize_intensities(spectrum_1)
@@ -104,33 +100,14 @@ def test_modified_cosine_order_of_input_spectrums():
     assert score_1_2 == score_2_1, "Expected that the order of the arguments would not matter."
 
 
-def test_modified_cosine_with_mass_shift_5_no_matches_expected():
-    """Test modified cosine on two spectra with no expected matches."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 200, 300], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
-                          metadata={"precursor_mz": 1000.0})
-
-    spectrum_2 = Spectrum(mz=numpy.array([120, 220, 320], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
-                          metadata={"precursor_mz": 1005})
-
-    norm_spectrum_1 = normalize_intensities(spectrum_1)
-    norm_spectrum_2 = normalize_intensities(spectrum_2)
-    modified_cosine = ModifiedCosine(tolerance=1.0)
-    score = modified_cosine.pair(norm_spectrum_1, norm_spectrum_2)
-
-    assert score["score"] == pytest.approx(0.0, 1e-5), "Expected different modified cosine score."
-    assert score["matches"] == 0, "Expected 0 matching peaks."
-
-
 def test_modified_cosine_precursor_mz_as_invalid_string():
     """Test modified cosine on two spectra with precursor_mz given as string."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 200, 300], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
+    spectrum_1 = Spectrum(mz=np.array([100, 200, 300], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
                           metadata={"precursor_mz": 1000.0})
 
-    spectrum_2 = Spectrum(mz=numpy.array([120, 220, 320], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
+    spectrum_2 = Spectrum(mz=np.array([120, 220, 320], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
                           metadata={"precursor_mz": "mz 1005.0"})
 
     norm_spectrum_1 = normalize_intensities(spectrum_1)
@@ -145,13 +122,13 @@ def test_modified_cosine_precursor_mz_as_invalid_string():
 
 def test_modified_cosine_precursor_mz_as_string(caplog):
     """Test modified cosine on two spectra with precursor_mz given as string."""
-    spectrum_1 = Spectrum(mz=numpy.array([100, 200, 300], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
+    spectrum_1 = Spectrum(mz=np.array([100, 200, 300], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
                           metadata={"precursor_mz": 1000.0},
                           metadata_harmonization=False)
 
-    spectrum_2 = Spectrum(mz=numpy.array([120, 220, 320], dtype="float"),
-                          intensities=numpy.array([10, 10, 500], dtype="float"),
+    spectrum_2 = Spectrum(mz=np.array([120, 220, 320], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
                           metadata={"precursor_mz": "1005.0"},
                           metadata_harmonization=False)
 
diff --git a/tests/test_neutral_losses_cosine.py b/tests/test_neutral_losses_cosine.py
new file mode 100644
index 00000000..d1e865cd
--- /dev/null
+++ b/tests/test_neutral_losses_cosine.py
@@ -0,0 +1,151 @@
+import numpy as np
+import pytest
+from matchms import Spectrum
+from matchms.filtering import normalize_intensities
+from matchms.similarity import NeutralLossesCosine
+from .builder_Spectrum import SpectrumBuilder
+
+
+def compute_expected_score(spectrum_1, spectrum_2, matches):
+    mass1 = spectrum_1.get("precursor_mz")
+    mass2 = spectrum_2.get("precursor_mz")
+    spec1 = spectrum_1.peaks.intensities[np.where(spectrum_1.peaks.mz < mass1)]
+    spec2 = spectrum_2.peaks.intensities[np.where(spectrum_2.peaks.mz < mass2)]
+    peak_pairs_multiplied = 0
+    for match in matches:
+        peak_pairs_multiplied += spec1[match[0]] * spec2[match[1]]
+    return peak_pairs_multiplied / np.sqrt(np.sum(spec1 ** 2) * np.sum(spec2 ** 2))
+
+
+def test_neutral_losses_cosine_without_precursor_mz():
+    """Test without precursor-m/z. Should raise assertion error."""
+    mz = np.array([100, 150, 200], dtype="float")
+    intensities = np.array([700, 200, 100], dtype="float")
+    builder = SpectrumBuilder()
+    spectrum_1 = builder.with_mz(mz).with_intensities(intensities).build()
+    spectrum_2 = builder.with_mz(mz).with_intensities(intensities).build()
+
+    norm_spectrum_1 = normalize_intensities(spectrum_1)
+    norm_spectrum_2 = normalize_intensities(spectrum_2)
+    neutral_losses_cosine = NeutralLossesCosine()
+
+    with pytest.raises(AssertionError) as msg:
+        neutral_losses_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+
+    expected_message = "Precursor_mz missing. Apply 'add_precursor_mz' filter first."
+    assert str(msg.value) == expected_message
+
+
+@pytest.mark.parametrize("peaks, tolerance, masses, expected_matches", [
+    [
+        [
+            [np.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"), np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float")],
+            [np.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"), np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float")]
+        ],
+        None, (1000.0, 1005.0), [(0, 1), (2, 2)]
+    ], [
+        [
+            [np.array([100, 299, 300, 301, 500, 510], dtype="float"), np.array([10, 500, 100, 200, 20, 100], dtype="float")],
+            [np.array([105, 305, 306, 505, 517], dtype="float"), np.array([10, 500, 100, 20, 100], dtype="float")],
+        ],
+        2.0, (1000.0, 1005.0), [(0, 0), (1, 1), (3, 2), (4, 3), (5, 4)]
+    ], [
+        [
+            [np.array([100, 110, 200, 300, 400, 500], dtype="float"), np.array([100, 50, 1, 80, 1, 1], dtype="float")],
+            [np.array([110, 200, 300, 310, 700], dtype="float"), np.array([100, 1, 90, 50, 1], dtype="float")],
+        ],
+        None, (1000.0, 1010.0), [(0, 0), (3, 3)]
+    ], [
+        [
+            [np.array([100, 200, 300], dtype="float"), np.array([10, 10, 500], dtype="float")],
+            [np.array([120, 220, 320], dtype="float"), np.array([10, 10, 500], dtype="float")],
+        ],
+        8.0, (1000.0, 1010.0), []
+    ], [
+        [
+            [np.array([100, 200, 300], dtype="float"), np.array([0.1, 1.0, 0.5], dtype="float")],
+            [np.array([120, 220, 320], dtype="float"), np.array([10, 10, 500], dtype="float")],
+        ],
+        2.0, (270.0, 289.5), [(0, 0), (1, 1)]
+    ]
+])
+def test_neutral_losses_cosine_with_mass_shift(peaks, tolerance, masses, expected_matches):
+    """Test neutral losses cosine on two spectra with mass shift."""
+    builder = SpectrumBuilder()
+    spectrum_1 = builder.with_mz(peaks[0][0]).with_intensities(peaks[0][1]).with_metadata(metadata={"precursor_mz": masses[0]}).build()
+    spectrum_2 = builder.with_mz(peaks[1][0]).with_intensities(peaks[1][1]).with_metadata(metadata={"precursor_mz": masses[1]}).build()
+
+    norm_spectrum_1 = normalize_intensities(spectrum_1)
+    norm_spectrum_2 = normalize_intensities(spectrum_2)
+    if tolerance is None:
+        neutral_losses_cosine = NeutralLossesCosine()
+    else:
+        neutral_losses_cosine = NeutralLossesCosine(tolerance=tolerance)
+
+    score = neutral_losses_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+    expected_score = compute_expected_score(norm_spectrum_1, norm_spectrum_2, expected_matches)
+    assert score["score"] == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
+    assert score["matches"] == len(expected_matches), "Expected differnt number of matching peaks."
+
+
+def test_neutral_losses_cosine_order_of_input_spectrums():
+    """Test neutral losses cosine on two spectra in changing order."""
+    spectrum_1 = Spectrum(mz=np.array([100, 150, 200, 300, 500, 510, 1100], dtype="float"),
+                          intensities=np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
+                          metadata={"precursor_mz": 1000.0})
+
+    spectrum_2 = Spectrum(mz=np.array([55, 105, 205, 304.5, 494.5, 515.5, 1045], dtype="float"),
+                          intensities=np.array([700, 200, 100, 1000, 200, 5, 500], dtype="float"),
+                          metadata={"precursor_mz": 1005.0})
+
+    norm_spectrum_1 = normalize_intensities(spectrum_1)
+    norm_spectrum_2 = normalize_intensities(spectrum_2)
+    neutral_losses_cosine = NeutralLossesCosine(tolerance=2.0)
+    score_1_2 = neutral_losses_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+    score_2_1 = neutral_losses_cosine.pair(norm_spectrum_2, norm_spectrum_1)
+
+    assert score_1_2["score"] == score_2_1["score"], "Expected that the order of the arguments would not matter."
+    assert score_1_2 == score_2_1, "Expected that the order of the arguments would not matter."
+
+
+def test_neutral_losses_cosine_precursor_mz_as_invalid_string():
+    """Test neutral losses cosine on two spectra with precursor_mz given as string."""
+    spectrum_1 = Spectrum(mz=np.array([100, 200, 300], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
+                          metadata={"precursor_mz": 1000.0})
+
+    spectrum_2 = Spectrum(mz=np.array([120, 220, 320], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
+                          metadata={"precursor_mz": "mz 1005.0"})
+
+    norm_spectrum_1 = normalize_intensities(spectrum_1)
+    norm_spectrum_2 = normalize_intensities(spectrum_2)
+    neutral_losses_cosine = NeutralLossesCosine(tolerance=1.0)
+    with pytest.raises(AssertionError) as msg:
+        _ = neutral_losses_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+
+    expected_message = "Precursor_mz missing. Apply 'add_precursor_mz' filter first."
+    assert str(msg.value) == expected_message
+
+
+def test_neutral_losses_cosine_precursor_mz_as_string(caplog):
+    """Test neutral losses cosine on two spectra with precursor_mz given as string."""
+    spectrum_1 = Spectrum(mz=np.array([100, 200, 300], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
+                          metadata={"precursor_mz": 1000.0},
+                          metadata_harmonization=False)
+
+    spectrum_2 = Spectrum(mz=np.array([120, 220, 320], dtype="float"),
+                          intensities=np.array([10, 10, 500], dtype="float"),
+                          metadata={"precursor_mz": "1005.0"},
+                          metadata_harmonization=False)
+
+    norm_spectrum_1 = normalize_intensities(spectrum_1)
+    norm_spectrum_2 = normalize_intensities(spectrum_2)
+    neutral_losses_cosine = NeutralLossesCosine(tolerance=1.0)
+    score = neutral_losses_cosine.pair(norm_spectrum_1, norm_spectrum_2)
+
+    assert score["score"] == pytest.approx(0.0, 1e-5), "Expected different neutral losses cosine score."
+    assert score["matches"] == 0, "Expected 0 matching peaks."
+    expected_msg = "Precursor_mz must be of type int or float. Apply 'add_precursor_mz' filter first."
+    assert expected_msg in caplog.text, "Expected different log message"