Repository: matchms/matchms
Pull request #112 (2020-08-03 14:13:57)

Add intensity and mz power options to Cosine Hungarian
Update ``CosineHungarian`` to be analogue to ``CosineGreedy`` after PR #98 

--> After merging this it seems time for a __new release__!

florian-huber: @jspaaks I slightly modified the test function. I added some spectra that I used to figure out if it indeed gave me the same results than the "original" score, which contains some non-trival peak mixture that tends to go wrong with too simple cosine score implementation. But feel free to change.
florian-huber: Just realized that your function does return the score, but not the number of matching peaks. MS people often select scores based on how many matching peaks this is based on (``min_match``). 
I now added it and hope I didn't mess up your code too much (feel free to implement it differently).
jspaaks: Hi Florian,
just a note on process: it is customary to not append commits to somebody else's pull request --instead you should leave your suggestions as comments, or possibly as issues to work on later. 

Could you roll back your commits on ``58-cosine-greedy``? The easiest way to do so is probably

1. make a branch ``cosine-greedy-tuple-output`` off of this branch (including your changes)

    ```
    git checkout 58-cosine-greedy
    git checkout -b cosine-greedy-tuple-output
    git checkout 58-cosine-greedy
    ```
1. Revert the commits by adding a new commit:

    ```
    git revert 788ca1d
    ```
1. Resolve the conflicts by picking the HEAD part and deleting the other part.
1. git add and commit your resolved file, then push the resolve commit to this branch.

If I would be doing the review, I would add one or more comments, something along the lines of "this unit test is not discriminatory enough", "units test should use the same spectrum combination", and/or "add a unit test that tests whether the argument order matters" with the permalink to the code that I'm talking about, or by leaving it as a comment in the 'margin' of the diff.

Secondly, I would make a new issue that says "CosineGreedy should return score as well as n_peaks", and in the body describe why you think that is needed. That could then be the start of a discussion there, resulting in a plan on how to implement it, which then would be done in a branch, which could then become its own PR.

As to your question https://github.com/matchms/matchms/pull/112#pullrequestreview-393537745 about performance: I have not compared one against the other, but I would assume that vectorized matrix operations are faster than looping, even for spectrums with many peaks. 

Thanks for reviewing and let me know if I can help you in some way,
-Jurriaan

jspaaks: Hi Florian, thanks for rolling back. I made some changes, and suggested that we make the "return multiple values" issue a new issue, such that we can discuss the implications there, without blocking further developments based on this PR. Let me know what you think,
-Jurriaan

florian-huber: Hi Juriaan, looks good. Thanks!
I made a new issue #118 for later modifying the cosine score returns.


Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 0d9c0cf3..8afbeda4 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -10,7 +10,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Added
 
 - Read msp files to create Spectrum objects from it [#102](https://github.com/matchms/matchms/pull/102)
-- Peak weighting option for cosine score and modified cosine score [#96](https://github.com/matchms/matchms/issues/96)
+- Peak weighting option for CosineGreedy and ModifiedCosine score [#96](https://github.com/matchms/matchms/issues/96)
+- Peak weighting option for CosineHungarian score [#112](https://github.com/matchms/matchms/pull/112)
 - Similarity score based on comparing parent masses [#79](https://github.com/matchms/matchms/pull/79)
 - Method for instantiating a spectrum from the metabolomics USI [#93](https://github.com/matchms/matchms/pull/93)
 
diff --git a/matchms/similarity/CosineHungarian.py b/matchms/similarity/CosineHungarian.py
index 414a7cc2..b4364e1b 100644
--- a/matchms/similarity/CosineHungarian.py
+++ b/matchms/similarity/CosineHungarian.py
@@ -18,15 +18,22 @@ class CosineHungarian:
     does represent a mathematically proper solution to the problem.
 
     """
-    def __init__(self, tolerance=0.1):
+    def __init__(self, tolerance: float = 0.1, mz_power: float = 0.0,
+                 intensity_power: float = 1.0):
         """
-
-        Args:
-        ----
-        tolerance: float
-            Peaks will be considered a match when <= tolerance appart. Default is 0.1.
+        Parameters
+        ----------
+        tolerance:
+            Peaks will be considered a match when <= tolerance apart. Default is 0.1.
+        mz_power:
+            The power to raise m/z 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.
         """
         self.tolerance = tolerance
+        self.mz_power = mz_power
+        self.intensity_power = intensity_power
 
     def __call__(self, spectrum1: SpectrumType, spectrum2: SpectrumType) -> Tuple[float, int]:
         """Calculate cosine score between two spectra.
@@ -37,10 +44,17 @@ class CosineHungarian:
             Input spectrum 1.
         spectrum2: SpectrumType
             Input spectrum 2.
+
+        Returns:
+        --------
+
+        Tuple with cosine score and number of matched peaks.
         """
         def get_matching_pairs():
             """Get pairs of peaks that match within the given tolerance."""
-            matching_pairs = collect_peak_pairs(spec1, spec2, self.tolerance, shift=0.0)
+            matching_pairs = collect_peak_pairs(spec1, spec2, self.tolerance, shift=0.0,
+                                                mz_power=self.mz_power,
+                                                intensity_power=self.intensity_power)
             return sorted(matching_pairs, key=lambda x: x[2], reverse=True)
 
         def get_matching_pairs_matrix():
@@ -77,7 +91,11 @@ class CosineHungarian:
             if matching_pairs_matrix is not None:
                 score, used_matches = solve_hungarian()
                 # Normalize score:
-                score = score/max(numpy.sum(spec1[:, 1]**2), numpy.sum(spec2[:, 1]**2))
+                spec1_power = numpy.power(spec1[:, 0], self.mz_power) \
+                    * numpy.power(spec1[:, 1], self.intensity_power)
+                spec2_power = numpy.power(spec2[:, 0], self.mz_power) \
+                    * numpy.power(spec2[:, 1], self.intensity_power)
+                score = score/(numpy.sqrt(numpy.sum(spec1_power**2)) * numpy.sqrt(numpy.sum(spec2_power**2)))
                 return score, len(used_matches)
             return 0.0, 0
 


Test patch:
diff --git a/tests/test_cosine_hungarian.py b/tests/test_cosine_hungarian.py
index 569df705..27f2d883 100644
--- a/tests/test_cosine_hungarian.py
+++ b/tests/test_cosine_hungarian.py
@@ -18,7 +18,9 @@ def test_cosine_hungarian_without_parameters():
     expected_matches = [0, 1, 4]  # Those peaks have matching mz values (within given tolerance)
     multiply_matching_intensities = spectrum_1.peaks.intensities[expected_matches] \
         * spectrum_2.peaks.intensities[expected_matches]
-    expected_score = multiply_matching_intensities.sum() / (spectrum_1.peaks.intensities ** 2).sum()
+    denominator = numpy.sqrt((spectrum_1.peaks.intensities ** 2).sum()) \
+        * numpy.sqrt((spectrum_2.peaks.intensities ** 2).sum())
+    expected_score = multiply_matching_intensities.sum() / denominator
 
     assert score == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
     assert n_matches == len(expected_matches), "Expected different number of matching peaks."
@@ -38,7 +40,9 @@ def test_cosine_hungarian_with_tolerance_0_2():
     expected_matches = [[0, 2, 3], [0, 1, 2]]  # Those peaks have matching mz values (within given tolerance)
     multiply_matching_intensities = spectrum_1.peaks.intensities[expected_matches[0]] \
         * spectrum_2.peaks.intensities[expected_matches[1]]
-    expected_score = multiply_matching_intensities.sum() / (spectrum_1.peaks.intensities ** 2).sum()
+    denominator = numpy.sqrt((spectrum_1.peaks.intensities ** 2).sum()) \
+        * numpy.sqrt((spectrum_2.peaks.intensities ** 2).sum())
+    expected_score = multiply_matching_intensities.sum() / denominator
 
     assert score == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
     assert n_matches == len(expected_matches[0]), "Expected different number of matching peaks."
@@ -58,7 +62,9 @@ def test_cosine_hungarian_with_tolerance_2_0():
     expected_matches = [[0, 1, 3, 4], [0, 1, 2, 3]]  # Those peaks have matching mz values (within given tolerance)
     multiply_matching_intensities = spectrum_1.peaks.intensities[expected_matches[0]] \
         * spectrum_2.peaks.intensities[expected_matches[1]]
-    expected_score = multiply_matching_intensities.sum() / (spectrum_1.peaks.intensities ** 2).sum()
+    denominator = numpy.sqrt((spectrum_1.peaks.intensities ** 2).sum()) \
+        * numpy.sqrt((spectrum_2.peaks.intensities ** 2).sum())
+    expected_score = multiply_matching_intensities.sum() / denominator
 
     assert score == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
     assert n_matches == len(expected_matches[0]), "Expected different number of matching peaks."
@@ -112,3 +118,33 @@ def test_cosine_hungarian_case_without_matches():
     score, n_matches = cosine_hungarian(spectrum_1, spectrum_2)
     assert score == 0.0, "Expected different cosine score."
     assert n_matches == 0, "Expected different number of matching peaks."
+
+
+def test_cosine_hungarian_with_peak_powers():
+    """Compare output cosine score with own calculation on simple dummy spectrums.
+    Here testing the options to raise peak intensities to given powers.
+    """
+    mz_power = 0.5
+    intensity_power = 2.0
+    spectrum_1 = Spectrum(mz=numpy.array([100, 200, 300, 500, 510], dtype="float"),
+                          intensities=numpy.array([0.1, 0.2, 1.0, 0.3, 0.4], dtype="float"))
+
+    spectrum_2 = Spectrum(mz=numpy.array([100, 200, 290, 490, 510], dtype="float"),
+                          intensities=numpy.array([0.1, 0.2, 1.0, 0.3, 0.4], dtype="float"))
+    cosine_hungarian = CosineHungarian(tolerance=1.0, mz_power=mz_power, intensity_power=intensity_power)
+    score, n_matches = cosine_hungarian(spectrum_1, spectrum_2)
+
+    # Derive expected cosine score
+    matches = [0, 1, 4]  # Those peaks have matching mz values (within given tolerance)
+    intensity1 = spectrum_1.peaks.intensities
+    mz1 = spectrum_1.peaks.mz
+    intensity2 = spectrum_2.peaks.intensities
+    mz2 = spectrum_2.peaks.mz
+    multiply_matching_intensities = (mz1[matches] ** mz_power) * (intensity1[matches] ** intensity_power) \
+        * (mz2[matches] ** mz_power) * (intensity2[matches] ** intensity_power)
+    denominator = numpy.sqrt((((mz1 ** mz_power) * (intensity1 ** intensity_power)) ** 2).sum()) \
+        * numpy.sqrt((((mz2 ** mz_power) * (intensity2 ** intensity_power)) ** 2).sum())
+    expected_score = multiply_matching_intensities.sum() / denominator
+
+    assert score == pytest.approx(expected_score, 0.0001), "Expected different cosine score."
+    assert n_matches == len(matches), "Expected different number of matching peaks."