Merge branch '286_ph_inten_min_abs' into 'dev'

Add an absolute intensity floor for LC-MS peak picking

See merge request mass-spectrometry/corems!272

Author: corilo

corems/encapsulation/factory/processingSetting.py

@@ -323,6 +323,11 @@ class LiquidChromatographSetting:
         Should be greater to or equal to ph_inten_min_rel.
         Called within the PH_Calculations.find_mass_features() method.
         Default is 0.001.
+    ph_inten_min_abs : float, optional
+        Absolute minimum intensity for points that enter persistent homology.
+        Applied with ph_inten_min_rel on profile and centroided data.
+        The relative floor is still a fraction of the unfiltered maximum.
+        0 disables the absolute floor. Default is 0.0.
     mass_feature_cluster_mz_tolerance_rel : float, optional
         Relative m/z tolerance to use for clustering mass features.
         Used for both "persistent homology" and "centroided_persistent_homology"
@@ -438,6 +443,7 @@ class LiquidChromatographSetting:
     ph_smooth_radius_scan = 1
     ph_inten_min_rel = 0.001
     ph_persis_min_rel = 0.001
+    ph_inten_min_abs: float = 0.0
 
     # Parameters used to cluster mass features
     mass_feature_cluster_mz_tolerance_rel: float = 5e-6

corems/mass_spectra/calc/lc_calc.py

@@ -2141,7 +2141,10 @@ def find_mass_features_ph(self, ms_level=1, grid=True, targeted_search=False, ta
                 self.mass_features = {}
                 return
         else:
-            threshold = self.parameters.lc_ms.ph_inten_min_rel * data.intensity.max()
+            threshold = max(
+                self.parameters.lc_ms.ph_inten_min_rel * data.intensity.max(),
+                self.parameters.lc_ms.ph_inten_min_abs,
+            )
             persistence_threshold = (
                 self.parameters.lc_ms.ph_persis_min_rel * data.intensity.max()
             )
@@ -2469,7 +2472,10 @@ def find_mass_features_ph_centroid(self, ms_level=1, targeted_search=False, targ
         else:
             # Normal mode with threshold
             max_intensity = data_with_time["intensity"].max()
-            threshold = self.parameters.lc_ms.ph_inten_min_rel * max_intensity
+            threshold = max(
+                self.parameters.lc_ms.ph_inten_min_rel * max_intensity,
+                self.parameters.lc_ms.ph_inten_min_abs,
+            )
             valid_mask = data_with_time["intensity"].notna() & (data_with_time["intensity"] > threshold)
             required_cols = ["mz", "intensity", "scan", "scan_time"]
             data_thres = data_with_time.loc[valid_mask, required_cols].copy()

tests/test_ph_inten_min_abs.py

@@ -0,0 +1,123 @@
+"""Absolute intensity floor before LC-MS persistent-homology peak picking."""
+
+import pandas as pd
+
+from corems.encapsulation.factory.parameters import LCMSParameters
+from corems.encapsulation.factory.processingSetting import (
+    LiquidChromatographSetting,
+    settings_to_export_dict,
+)
+from corems.mass_spectra.calc.lc_calc import PHCalculations
+
+# Max intensity is 1000. A relative floor of 0.1 keeps intensity > 100.
+_INTENSITIES = (50.0, 200.0, 1000.0)
+
+
+def _frame():
+    data = pd.DataFrame(
+        {
+            "mz": [200.0, 201.0, 202.0],
+            "scan": [1, 1, 1],
+            "intensity": list(_INTENSITIES),
+        }
+    )
+    scan_df = pd.DataFrame({"scan": [1], "scan_time": [1.0]})
+    return data, scan_df
+
+
+def _target_window():
+    return {
+        "target_mz_list": [201.0],
+        "target_rt_list": [1.0],
+        "mz_tolerance_ppm": 20_000,
+        "rt_tolerance": 0.5,
+    }
+
+
+class _ProfileHost(PHCalculations):
+    def __init__(self, abs_floor):
+        data, scan_df = _frame()
+        self.parameters = LCMSParameters()
+        self.parameters.lc_ms.ph_inten_min_rel = 0.1
+        self.parameters.lc_ms.ph_inten_min_abs = abs_floor
+        self.parameters.lc_ms.verbose_processing = False
+        self._ms_unprocessed = {1: data}
+        self.scan_df = scan_df
+        self.mass_features = {}
+        self.kept = None
+
+    def check_if_grid(self, data):
+        return True
+
+    def _record(self, data_thres):
+        self.kept = sorted(float(value) for value in data_thres["intensity"])
+
+    def _find_mass_features_ph_single(self, data_thres, *args, **kwargs):
+        self._record(data_thres)
+
+    def _find_mass_features_ph_partition(self, data_thres, *args, **kwargs):
+        self._record(data_thres)
+
+
+class _CentroidHost(PHCalculations):
+    def __init__(self, abs_floor):
+        data, scan_df = _frame()
+        self.parameters = LCMSParameters()
+        self.parameters.lc_ms.ph_inten_min_rel = 0.1
+        self.parameters.lc_ms.ph_inten_min_abs = abs_floor
+        self.parameters.lc_ms.verbose_processing = False
+        self._ms_unprocessed = {1: data}
+        self.scan_df = scan_df
+        self.mass_features = {}
+        self.kept = None
+
+    def roll_up_dataframe(self, df, **kwargs):
+        if len(df):
+            self.kept = sorted(float(value) for value in df["intensity"])
+        return df.iloc[0:0].copy()
+
+
+def test_ph_inten_min_abs_defaults_to_zero_and_is_exported():
+    settings = LiquidChromatographSetting()
+    assert settings.ph_inten_min_abs == 0.0
+    assert settings_to_export_dict(settings)["ph_inten_min_abs"] == 0.0
+
+
+def test_zero_absolute_floor_keeps_the_relative_mask():
+    profile = _ProfileHost(0.0)
+    profile.find_mass_features_ph()
+    centroid = _CentroidHost(0.0)
+    centroid.find_mass_features_ph_centroid()
+    assert profile.kept == [200.0, 1000.0]
+    assert centroid.kept == [200.0, 1000.0]
+
+
+def test_absolute_floor_drops_points_the_relative_floor_would_keep():
+    profile = _ProfileHost(500.0)
+    profile.find_mass_features_ph()
+    centroid = _CentroidHost(500.0)
+    centroid.find_mass_features_ph_centroid()
+    assert profile.kept == [1000.0]
+    assert centroid.kept == [1000.0]
+
+
+def test_relative_floor_still_drops_points_above_the_absolute_floor():
+    profile = _ProfileHost(10.0)
+    profile.find_mass_features_ph()
+    centroid = _CentroidHost(10.0)
+    centroid.find_mass_features_ph_centroid()
+    assert profile.kept == [200.0, 1000.0]
+    assert centroid.kept == [200.0, 1000.0]
+
+
+def test_targeted_search_ignores_the_absolute_floor():
+    profile = _ProfileHost(500.0)
+    profile.find_mass_features_ph(
+        targeted_search=True, target_search_dict=_target_window()
+    )
+    centroid = _CentroidHost(500.0)
+    centroid.find_mass_features_ph_centroid(
+        targeted_search=True, target_search_dict=_target_window()
+    )
+    assert profile.kept == [50.0, 200.0, 1000.0]
+    assert centroid.kept == [50.0, 200.0, 1000.0]