corems.molecular_id.search.lcms_spectral_search

View Source

  1import re
  2
  3import numpy as np
  4
  5from corems.molecular_id.factory.spectrum_search_results import SpectrumSearchResults
  6
  7
  8class LCMSSpectralSearch:
  9    """
 10    Methods for searching LCMS spectra.
 11
 12    This class is designed to be a mixin class for the :obj:`~corems.mass_spectra.factory.lc_class.LCMSBase` class.
 13
 14    """
 15
 16    @staticmethod
 17    def get_more_match_quals(
 18        query_mz_arr, lib_entry, mz_tol_da=0.1, include_fragment_types=False
 19    ):
 20        """
 21        Return additional match qualities between query and library entry.
 22
 23        Parameters
 24        ----------
 25        query_mz_arr : np.array
 26            Array of query spectrum. Shape (N, 2), with m/z in the first column
 27            and abundance in the second.
 28        lib_entry : dict
 29            Library spectrum entry, with 'mz' key containing the spectrum in
 30            the format (mz, abundance),(mz, abundance), i.e. from MetabRef.
 31        mz_tol_da : float, optional
 32            Tolerance in Da for matching peaks (in MS2). Default is 0.1.
 33        include_fragment_types : bool, optional
 34            If True, include fragment type comparisons in output.
 35            Defaults to False.
 36
 37        Returns
 38        -------
 39        tuple
 40            Tuple of (query_in_lib, query_in_lib_fract, lib_in_query, lib_in_query_fract, query_frags, lib_frags, lib_precursor_mz).
 41
 42        Notes
 43        -----
 44        query_in_lib : int
 45            Number of peaks in query that are present in the library entry (within mz_tol_da).
 46        query_in_lib_fract : float
 47            Fraction of peaks in query that are present in the library entry (within mz_tol_da).
 48        lib_in_query : int
 49            Number of peaks in the library entry that are present in the query (within mz_tol_da).
 50        lib_in_query_fract : float
 51            Fraction of peaks in the library entry that are present in the query (within mz_tol_da).
 52        query_frags : list
 53            List of unique fragment types present in the query, generally 'MLF' or 'LSF' or both.
 54        lib_frags : list
 55            List of unique fragment types present in the library entry, generally 'MLF' or 'LSF' or both.
 56
 57        Raises
 58        ------
 59        ValueError
 60            If library entry does not have 'fragment_types' key and include_fragment_types is True.
 61
 62        """
 63
 64        if "mz" in lib_entry.keys():
 65            # Get the original mz values from the library entry
 66            lib_mzs = np.array(
 67                re.findall(r"\(([^,]+),([^)]+)\)", lib_entry["mz"]), dtype=float
 68            ).reshape(-1, 2)[:, 0]
 69        elif "peaks" in lib_entry.keys() and lib_entry["peaks"] is not None:
 70            lib_mzs = lib_entry["peaks"][:, 0]
 71
 72        # Get count and fraction of peaks in query that are in lib entry
 73        query_in_lib = 0
 74        for peak in query_mz_arr:
 75            if np.any(np.isclose(lib_mzs, peak, atol=mz_tol_da)):
 76                query_in_lib += 1
 77        query_in_lib_fract = query_in_lib / len(query_mz_arr)
 78
 79        # Get count and fraction of peaks in lib that are in query
 80        lib_in_query = 0
 81        for peak in lib_mzs:
 82            if np.any(np.isclose(query_mz_arr, peak, atol=mz_tol_da)):
 83                lib_in_query += 1
 84        lib_in_query_fract = lib_in_query / len(lib_mzs)
 85
 86        if include_fragment_types:
 87            # Check that fragment types are present in the library entry
 88            if "fragment_types" not in lib_entry.keys():
 89                raise ValueError(
 90                    "Flash entropy library entry must have 'fragment_types' key to include fragment types in output."
 91                )
 92
 93            # Get types of fragments in the lib entry and convert it to a list on ", "
 94            lib_frags = [x.strip() for x in lib_entry["fragment_types"].split(",")]
 95            # make list of the fragment types that are present in the query spectrum
 96            lib_in_query_ids = list(
 97                set(
 98                    [
 99                        ind
100                        for ind, x in enumerate(lib_mzs)
101                        if len(np.where(np.isclose(query_mz_arr, x, atol=mz_tol_da))[0])
102                        > 0
103                    ]
104                )
105            )
106            query_frags = list(set([lib_frags[x] for x in lib_in_query_ids]))
107            lib_frags = list(set(lib_frags))
108
109        else:
110            query_frags = None
111            lib_frags = None
112
113        return (
114            query_in_lib,
115            query_in_lib_fract,
116            lib_in_query,
117            lib_in_query_fract,
118            query_frags,
119            lib_frags,
120        )
121
122    def fe_search(
123        self,
124        scan_list,
125        fe_lib,
126        precursor_mz_list=[],
127        use_mass_features=True,
128        peak_sep_da=0.01,
129        get_additional_metrics=True,
130        accumulate_results=False,
131    ):
132        """
133        Search LCMS spectra using a FlashEntropy approach.
134
135        Parameters
136        ----------
137        scan_list : list
138            List of scan numbers to search.
139        fe_lib : :obj:`~ms_entropy.FlashEntropySearch`
140            FlashEntropy Search instance.
141        precursor_mz_list : list, optional
142            List of precursor m/z values to search, by default [], which implies
143            matched with mass features; to enable this use_mass_features must be True.
144        use_mass_features : bool, optional
145            If True, use mass features to get precursor m/z values, by default True.
146            If True, will add search results to mass features' ms2_similarity_results attribute.
147        peak_sep_da : float, optional
148            Minimum separation between m/z peaks spectra in Da. This needs match the
149            approximate resolution of the search spectra and the FlashEntropySearch
150            instance, by default 0.01.
151        get_additional_metrics : bool, optional
152            If True, get additional metrics from FlashEntropy search, by default True.
153        accumulate_results : bool, optional
154            If True, accumulate results with existing spectral_search_results instead of
155            replacing them. This allows searching the same scans with multiple libraries
156            without overwriting previous results, by default False.
157
158        Returns
159        -------
160        None, but adds results to self.spectral_search_results and associates these
161        spectral_search_results with mass_features within the self.mass_features dictionary.
162
163        """
164        # Retrieve parameters from self
165        # include_fragment_types should used for lipids queries only, not general metabolomics
166        include_fragment_types = self.parameters.lc_ms.include_fragment_types
167        min_match_score = self.parameters.lc_ms.ms2_min_fe_score
168
169        # If precursor_mz_list is empty and use_mass_features is True, get precursor m/z values from mass features for each scan in scan_list
170        if use_mass_features and len(precursor_mz_list) == 0:
171            precursor_mz_list = []
172            for scan in scan_list:
173                mf_ids = [
174                    key
175                    for key, value in self.mass_features.items()
176                    if scan in value.ms2_mass_spectra
177                ]
178                precursor_mz = [
179                    value.mz
180                    for key, value in self.mass_features.items()
181                    if key in mf_ids
182                ]
183                precursor_mz_list.append(precursor_mz)
184
185        # Check that precursor_mz_list same length as scan_list, if not, raise error
186        if len(precursor_mz_list) != len(scan_list):
187            raise ValueError("Length of precursor_mz_list is not equal to scan_list.")
188
189        # Loop through each query spectrum / precursor match and save ids of db spectrum that are decent matches
190        overall_results_dict = {}
191        for i in np.arange(len(scan_list)):
192            scan_oi = scan_list[i]
193            if len(self._ms[scan_oi].mspeaks) > 0:
194                precursor_mzs = precursor_mz_list[i]
195                overall_results_dict[scan_oi] = {}
196                for precursor_mz in precursor_mzs:
197                    query_spectrum = fe_lib.clean_spectrum_for_search(
198                        precursor_mz=precursor_mz,
199                        peaks=np.vstack(
200                            (self._ms[scan_oi].mz_exp, self._ms[scan_oi].abundance)
201                        ).T,
202                        precursor_ions_removal_da=None,
203                        noise_threshold=self._ms[
204                            scan_oi
205                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
206                        / 100,
207                        min_ms2_difference_in_da=peak_sep_da,
208                    )
209                    search_results = fe_lib.search(
210                        precursor_mz=precursor_mz,
211                        peaks=query_spectrum,
212                        ms1_tolerance_in_da=self.parameters.mass_spectrum[
213                            "ms1"
214                        ].molecular_search.max_ppm_error
215                        * 10**-6
216                        * precursor_mz,
217                        ms2_tolerance_in_da=peak_sep_da * 0.5,
218                        method={"identity"},
219                        precursor_ions_removal_da=None,
220                        noise_threshold=self._ms[
221                            scan_oi
222                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
223                        / 100,
224                        target="cpu",
225                    )["identity_search"]
226                    match_inds = np.where(search_results > min_match_score)[0]
227
228                    # If any decent matches are found, add them to the results dictionary
229                    if len(match_inds) > 0:
230                        match_scores = search_results[match_inds]
231                        ref_ms_ids = [fe_lib[x]["id"] for x in match_inds]
232                        ref_mol_ids = [
233                            fe_lib[x]["molecular_data_id"] for x in match_inds
234                        ]
235                        ref_precursor_mzs = [
236                            fe_lib[x]["precursor_mz"] for x in match_inds
237                        ]
238                        ion_types = [fe_lib[x]["ion_type"] for x in match_inds]
239                        overall_results_dict[scan_oi][precursor_mz] = {
240                            "ref_mol_id": ref_mol_ids,
241                            "ref_ms_id": ref_ms_ids,
242                            "ref_precursor_mz": ref_precursor_mzs,
243                            "precursor_mz_error_ppm": [
244                                (precursor_mz - x) / precursor_mz * 10**6
245                                for x in ref_precursor_mzs
246                            ],
247                            "entropy_similarity": match_scores,
248                            "ref_ion_type": ion_types,
249                        }
250                        # Add database name, if present
251                        db_name = [
252                            fe_lib[x].get("database_name") for x in match_inds
253                        ]
254                        if db_name is not None:
255                            overall_results_dict[scan_oi][precursor_mz].update(
256                                {"database_name": db_name}
257                            )
258                        if get_additional_metrics:
259                            more_match_quals = [
260                                self.get_more_match_quals(
261                                    self._ms[scan_oi].mz_exp,
262                                    fe_lib[x],
263                                    mz_tol_da=peak_sep_da,
264                                    include_fragment_types=include_fragment_types,
265                                )
266                                for x in match_inds
267                            ]
268                            overall_results_dict[scan_oi][precursor_mz].update(
269                                {
270                                    "query_mz_in_ref_n": [
271                                        x[0] for x in more_match_quals
272                                    ],
273                                    "query_mz_in_ref_fract": [
274                                        x[1] for x in more_match_quals
275                                    ],
276                                    "ref_mz_in_query_n": [
277                                        x[2] for x in more_match_quals
278                                    ],
279                                    "ref_mz_in_query_fract": [
280                                        x[3] for x in more_match_quals
281                                    ],
282                                }
283                            )
284                            if include_fragment_types:
285                                overall_results_dict[scan_oi][precursor_mz].update(
286                                    {
287                                        "query_frag_types": [
288                                            x[4] for x in more_match_quals
289                                        ],
290                                        "ref_frag_types": [
291                                            x[5] for x in more_match_quals
292                                        ],
293                                    }
294                                )
295
296        # Drop scans with no results from dictionary
297        overall_results_dict = {k: v for k, v in overall_results_dict.items() if v}
298
299        # Cast each entry as a MS2SearchResults object
300        for scan_id in overall_results_dict.keys():
301            for precursor_mz in overall_results_dict[scan_id].keys():
302                ms2_spectrum = self._ms[scan_id]
303                ms2_search_results = overall_results_dict[scan_id][precursor_mz]
304                overall_results_dict[scan_id][precursor_mz] = SpectrumSearchResults(
305                    ms2_spectrum, precursor_mz, ms2_search_results
306                )
307
308        # Add MS2SearchResults to the existing spectral search results dictionary
309        if accumulate_results:
310            # Merge results with existing spectral_search_results
311            for scan_id, precursor_dict in overall_results_dict.items():
312                if scan_id in self.spectral_search_results:
313                    # Scan already has results, merge precursor_mz dictionaries
314                    self.spectral_search_results[scan_id].update(precursor_dict)
315                else:
316                    # New scan, add entire dictionary
317                    self.spectral_search_results[scan_id] = precursor_dict
318        else:
319            # Replace existing results (original behavior)
320            self.spectral_search_results.update(overall_results_dict)
321
322        # If there are mass features, associate the results with each mass feature
323        if len(self.mass_features) > 0:
324            # Determine which results to associate with mass features
325            if accumulate_results:
326                # When accumulating, only associate new results from this search
327                # to avoid duplicating previously associated results
328                results_to_associate = overall_results_dict
329            else:
330                # When not accumulating, clear existing associations and re-associate all results
331                for mass_feature_id in self.mass_features.keys():
332                    self.mass_features[mass_feature_id].ms2_similarity_results = []
333                results_to_associate = self.spectral_search_results
334            
335            for mass_feature_id, mass_feature in self.mass_features.items():
336                scan_ids = mass_feature.ms2_scan_numbers
337                for ms2_scan_id in scan_ids:
338                    precursor_mz = mass_feature.mz
339                    try:
340                        results_to_associate[ms2_scan_id][precursor_mz]
341                    except KeyError:
342                        pass
343                    else:
344                        self.mass_features[
345                            mass_feature_id
346                        ].ms2_similarity_results.append(
347                            results_to_associate[ms2_scan_id][precursor_mz]
348                        )

class LCMSSpectralSearch: View Source

  9class LCMSSpectralSearch:
 10    """
 11    Methods for searching LCMS spectra.
 12
 13    This class is designed to be a mixin class for the :obj:`~corems.mass_spectra.factory.lc_class.LCMSBase` class.
 14
 15    """
 16
 17    @staticmethod
 18    def get_more_match_quals(
 19        query_mz_arr, lib_entry, mz_tol_da=0.1, include_fragment_types=False
 20    ):
 21        """
 22        Return additional match qualities between query and library entry.
 23
 24        Parameters
 25        ----------
 26        query_mz_arr : np.array
 27            Array of query spectrum. Shape (N, 2), with m/z in the first column
 28            and abundance in the second.
 29        lib_entry : dict
 30            Library spectrum entry, with 'mz' key containing the spectrum in
 31            the format (mz, abundance),(mz, abundance), i.e. from MetabRef.
 32        mz_tol_da : float, optional
 33            Tolerance in Da for matching peaks (in MS2). Default is 0.1.
 34        include_fragment_types : bool, optional
 35            If True, include fragment type comparisons in output.
 36            Defaults to False.
 37
 38        Returns
 39        -------
 40        tuple
 41            Tuple of (query_in_lib, query_in_lib_fract, lib_in_query, lib_in_query_fract, query_frags, lib_frags, lib_precursor_mz).
 42
 43        Notes
 44        -----
 45        query_in_lib : int
 46            Number of peaks in query that are present in the library entry (within mz_tol_da).
 47        query_in_lib_fract : float
 48            Fraction of peaks in query that are present in the library entry (within mz_tol_da).
 49        lib_in_query : int
 50            Number of peaks in the library entry that are present in the query (within mz_tol_da).
 51        lib_in_query_fract : float
 52            Fraction of peaks in the library entry that are present in the query (within mz_tol_da).
 53        query_frags : list
 54            List of unique fragment types present in the query, generally 'MLF' or 'LSF' or both.
 55        lib_frags : list
 56            List of unique fragment types present in the library entry, generally 'MLF' or 'LSF' or both.
 57
 58        Raises
 59        ------
 60        ValueError
 61            If library entry does not have 'fragment_types' key and include_fragment_types is True.
 62
 63        """
 64
 65        if "mz" in lib_entry.keys():
 66            # Get the original mz values from the library entry
 67            lib_mzs = np.array(
 68                re.findall(r"\(([^,]+),([^)]+)\)", lib_entry["mz"]), dtype=float
 69            ).reshape(-1, 2)[:, 0]
 70        elif "peaks" in lib_entry.keys() and lib_entry["peaks"] is not None:
 71            lib_mzs = lib_entry["peaks"][:, 0]
 72
 73        # Get count and fraction of peaks in query that are in lib entry
 74        query_in_lib = 0
 75        for peak in query_mz_arr:
 76            if np.any(np.isclose(lib_mzs, peak, atol=mz_tol_da)):
 77                query_in_lib += 1
 78        query_in_lib_fract = query_in_lib / len(query_mz_arr)
 79
 80        # Get count and fraction of peaks in lib that are in query
 81        lib_in_query = 0
 82        for peak in lib_mzs:
 83            if np.any(np.isclose(query_mz_arr, peak, atol=mz_tol_da)):
 84                lib_in_query += 1
 85        lib_in_query_fract = lib_in_query / len(lib_mzs)
 86
 87        if include_fragment_types:
 88            # Check that fragment types are present in the library entry
 89            if "fragment_types" not in lib_entry.keys():
 90                raise ValueError(
 91                    "Flash entropy library entry must have 'fragment_types' key to include fragment types in output."
 92                )
 93
 94            # Get types of fragments in the lib entry and convert it to a list on ", "
 95            lib_frags = [x.strip() for x in lib_entry["fragment_types"].split(",")]
 96            # make list of the fragment types that are present in the query spectrum
 97            lib_in_query_ids = list(
 98                set(
 99                    [
100                        ind
101                        for ind, x in enumerate(lib_mzs)
102                        if len(np.where(np.isclose(query_mz_arr, x, atol=mz_tol_da))[0])
103                        > 0
104                    ]
105                )
106            )
107            query_frags = list(set([lib_frags[x] for x in lib_in_query_ids]))
108            lib_frags = list(set(lib_frags))
109
110        else:
111            query_frags = None
112            lib_frags = None
113
114        return (
115            query_in_lib,
116            query_in_lib_fract,
117            lib_in_query,
118            lib_in_query_fract,
119            query_frags,
120            lib_frags,
121        )
122
123    def fe_search(
124        self,
125        scan_list,
126        fe_lib,
127        precursor_mz_list=[],
128        use_mass_features=True,
129        peak_sep_da=0.01,
130        get_additional_metrics=True,
131        accumulate_results=False,
132    ):
133        """
134        Search LCMS spectra using a FlashEntropy approach.
135
136        Parameters
137        ----------
138        scan_list : list
139            List of scan numbers to search.
140        fe_lib : :obj:`~ms_entropy.FlashEntropySearch`
141            FlashEntropy Search instance.
142        precursor_mz_list : list, optional
143            List of precursor m/z values to search, by default [], which implies
144            matched with mass features; to enable this use_mass_features must be True.
145        use_mass_features : bool, optional
146            If True, use mass features to get precursor m/z values, by default True.
147            If True, will add search results to mass features' ms2_similarity_results attribute.
148        peak_sep_da : float, optional
149            Minimum separation between m/z peaks spectra in Da. This needs match the
150            approximate resolution of the search spectra and the FlashEntropySearch
151            instance, by default 0.01.
152        get_additional_metrics : bool, optional
153            If True, get additional metrics from FlashEntropy search, by default True.
154        accumulate_results : bool, optional
155            If True, accumulate results with existing spectral_search_results instead of
156            replacing them. This allows searching the same scans with multiple libraries
157            without overwriting previous results, by default False.
158
159        Returns
160        -------
161        None, but adds results to self.spectral_search_results and associates these
162        spectral_search_results with mass_features within the self.mass_features dictionary.
163
164        """
165        # Retrieve parameters from self
166        # include_fragment_types should used for lipids queries only, not general metabolomics
167        include_fragment_types = self.parameters.lc_ms.include_fragment_types
168        min_match_score = self.parameters.lc_ms.ms2_min_fe_score
169
170        # If precursor_mz_list is empty and use_mass_features is True, get precursor m/z values from mass features for each scan in scan_list
171        if use_mass_features and len(precursor_mz_list) == 0:
172            precursor_mz_list = []
173            for scan in scan_list:
174                mf_ids = [
175                    key
176                    for key, value in self.mass_features.items()
177                    if scan in value.ms2_mass_spectra
178                ]
179                precursor_mz = [
180                    value.mz
181                    for key, value in self.mass_features.items()
182                    if key in mf_ids
183                ]
184                precursor_mz_list.append(precursor_mz)
185
186        # Check that precursor_mz_list same length as scan_list, if not, raise error
187        if len(precursor_mz_list) != len(scan_list):
188            raise ValueError("Length of precursor_mz_list is not equal to scan_list.")
189
190        # Loop through each query spectrum / precursor match and save ids of db spectrum that are decent matches
191        overall_results_dict = {}
192        for i in np.arange(len(scan_list)):
193            scan_oi = scan_list[i]
194            if len(self._ms[scan_oi].mspeaks) > 0:
195                precursor_mzs = precursor_mz_list[i]
196                overall_results_dict[scan_oi] = {}
197                for precursor_mz in precursor_mzs:
198                    query_spectrum = fe_lib.clean_spectrum_for_search(
199                        precursor_mz=precursor_mz,
200                        peaks=np.vstack(
201                            (self._ms[scan_oi].mz_exp, self._ms[scan_oi].abundance)
202                        ).T,
203                        precursor_ions_removal_da=None,
204                        noise_threshold=self._ms[
205                            scan_oi
206                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
207                        / 100,
208                        min_ms2_difference_in_da=peak_sep_da,
209                    )
210                    search_results = fe_lib.search(
211                        precursor_mz=precursor_mz,
212                        peaks=query_spectrum,
213                        ms1_tolerance_in_da=self.parameters.mass_spectrum[
214                            "ms1"
215                        ].molecular_search.max_ppm_error
216                        * 10**-6
217                        * precursor_mz,
218                        ms2_tolerance_in_da=peak_sep_da * 0.5,
219                        method={"identity"},
220                        precursor_ions_removal_da=None,
221                        noise_threshold=self._ms[
222                            scan_oi
223                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
224                        / 100,
225                        target="cpu",
226                    )["identity_search"]
227                    match_inds = np.where(search_results > min_match_score)[0]
228
229                    # If any decent matches are found, add them to the results dictionary
230                    if len(match_inds) > 0:
231                        match_scores = search_results[match_inds]
232                        ref_ms_ids = [fe_lib[x]["id"] for x in match_inds]
233                        ref_mol_ids = [
234                            fe_lib[x]["molecular_data_id"] for x in match_inds
235                        ]
236                        ref_precursor_mzs = [
237                            fe_lib[x]["precursor_mz"] for x in match_inds
238                        ]
239                        ion_types = [fe_lib[x]["ion_type"] for x in match_inds]
240                        overall_results_dict[scan_oi][precursor_mz] = {
241                            "ref_mol_id": ref_mol_ids,
242                            "ref_ms_id": ref_ms_ids,
243                            "ref_precursor_mz": ref_precursor_mzs,
244                            "precursor_mz_error_ppm": [
245                                (precursor_mz - x) / precursor_mz * 10**6
246                                for x in ref_precursor_mzs
247                            ],
248                            "entropy_similarity": match_scores,
249                            "ref_ion_type": ion_types,
250                        }
251                        # Add database name, if present
252                        db_name = [
253                            fe_lib[x].get("database_name") for x in match_inds
254                        ]
255                        if db_name is not None:
256                            overall_results_dict[scan_oi][precursor_mz].update(
257                                {"database_name": db_name}
258                            )
259                        if get_additional_metrics:
260                            more_match_quals = [
261                                self.get_more_match_quals(
262                                    self._ms[scan_oi].mz_exp,
263                                    fe_lib[x],
264                                    mz_tol_da=peak_sep_da,
265                                    include_fragment_types=include_fragment_types,
266                                )
267                                for x in match_inds
268                            ]
269                            overall_results_dict[scan_oi][precursor_mz].update(
270                                {
271                                    "query_mz_in_ref_n": [
272                                        x[0] for x in more_match_quals
273                                    ],
274                                    "query_mz_in_ref_fract": [
275                                        x[1] for x in more_match_quals
276                                    ],
277                                    "ref_mz_in_query_n": [
278                                        x[2] for x in more_match_quals
279                                    ],
280                                    "ref_mz_in_query_fract": [
281                                        x[3] for x in more_match_quals
282                                    ],
283                                }
284                            )
285                            if include_fragment_types:
286                                overall_results_dict[scan_oi][precursor_mz].update(
287                                    {
288                                        "query_frag_types": [
289                                            x[4] for x in more_match_quals
290                                        ],
291                                        "ref_frag_types": [
292                                            x[5] for x in more_match_quals
293                                        ],
294                                    }
295                                )
296
297        # Drop scans with no results from dictionary
298        overall_results_dict = {k: v for k, v in overall_results_dict.items() if v}
299
300        # Cast each entry as a MS2SearchResults object
301        for scan_id in overall_results_dict.keys():
302            for precursor_mz in overall_results_dict[scan_id].keys():
303                ms2_spectrum = self._ms[scan_id]
304                ms2_search_results = overall_results_dict[scan_id][precursor_mz]
305                overall_results_dict[scan_id][precursor_mz] = SpectrumSearchResults(
306                    ms2_spectrum, precursor_mz, ms2_search_results
307                )
308
309        # Add MS2SearchResults to the existing spectral search results dictionary
310        if accumulate_results:
311            # Merge results with existing spectral_search_results
312            for scan_id, precursor_dict in overall_results_dict.items():
313                if scan_id in self.spectral_search_results:
314                    # Scan already has results, merge precursor_mz dictionaries
315                    self.spectral_search_results[scan_id].update(precursor_dict)
316                else:
317                    # New scan, add entire dictionary
318                    self.spectral_search_results[scan_id] = precursor_dict
319        else:
320            # Replace existing results (original behavior)
321            self.spectral_search_results.update(overall_results_dict)
322
323        # If there are mass features, associate the results with each mass feature
324        if len(self.mass_features) > 0:
325            # Determine which results to associate with mass features
326            if accumulate_results:
327                # When accumulating, only associate new results from this search
328                # to avoid duplicating previously associated results
329                results_to_associate = overall_results_dict
330            else:
331                # When not accumulating, clear existing associations and re-associate all results
332                for mass_feature_id in self.mass_features.keys():
333                    self.mass_features[mass_feature_id].ms2_similarity_results = []
334                results_to_associate = self.spectral_search_results
335            
336            for mass_feature_id, mass_feature in self.mass_features.items():
337                scan_ids = mass_feature.ms2_scan_numbers
338                for ms2_scan_id in scan_ids:
339                    precursor_mz = mass_feature.mz
340                    try:
341                        results_to_associate[ms2_scan_id][precursor_mz]
342                    except KeyError:
343                        pass
344                    else:
345                        self.mass_features[
346                            mass_feature_id
347                        ].ms2_similarity_results.append(
348                            results_to_associate[ms2_scan_id][precursor_mz]
349                        )

Methods for searching LCMS spectra.

This class is designed to be a mixin class for the ~corems.mass_spectra.factory.lc_class.LCMSBase class.

@staticmethod

def get_more_match_quals(query_mz_arr, lib_entry, mz_tol_da=0.1, include_fragment_types=False): View Source

 17    @staticmethod
 18    def get_more_match_quals(
 19        query_mz_arr, lib_entry, mz_tol_da=0.1, include_fragment_types=False
 20    ):
 21        """
 22        Return additional match qualities between query and library entry.
 23
 24        Parameters
 25        ----------
 26        query_mz_arr : np.array
 27            Array of query spectrum. Shape (N, 2), with m/z in the first column
 28            and abundance in the second.
 29        lib_entry : dict
 30            Library spectrum entry, with 'mz' key containing the spectrum in
 31            the format (mz, abundance),(mz, abundance), i.e. from MetabRef.
 32        mz_tol_da : float, optional
 33            Tolerance in Da for matching peaks (in MS2). Default is 0.1.
 34        include_fragment_types : bool, optional
 35            If True, include fragment type comparisons in output.
 36            Defaults to False.
 37
 38        Returns
 39        -------
 40        tuple
 41            Tuple of (query_in_lib, query_in_lib_fract, lib_in_query, lib_in_query_fract, query_frags, lib_frags, lib_precursor_mz).
 42
 43        Notes
 44        -----
 45        query_in_lib : int
 46            Number of peaks in query that are present in the library entry (within mz_tol_da).
 47        query_in_lib_fract : float
 48            Fraction of peaks in query that are present in the library entry (within mz_tol_da).
 49        lib_in_query : int
 50            Number of peaks in the library entry that are present in the query (within mz_tol_da).
 51        lib_in_query_fract : float
 52            Fraction of peaks in the library entry that are present in the query (within mz_tol_da).
 53        query_frags : list
 54            List of unique fragment types present in the query, generally 'MLF' or 'LSF' or both.
 55        lib_frags : list
 56            List of unique fragment types present in the library entry, generally 'MLF' or 'LSF' or both.
 57
 58        Raises
 59        ------
 60        ValueError
 61            If library entry does not have 'fragment_types' key and include_fragment_types is True.
 62
 63        """
 64
 65        if "mz" in lib_entry.keys():
 66            # Get the original mz values from the library entry
 67            lib_mzs = np.array(
 68                re.findall(r"\(([^,]+),([^)]+)\)", lib_entry["mz"]), dtype=float
 69            ).reshape(-1, 2)[:, 0]
 70        elif "peaks" in lib_entry.keys() and lib_entry["peaks"] is not None:
 71            lib_mzs = lib_entry["peaks"][:, 0]
 72
 73        # Get count and fraction of peaks in query that are in lib entry
 74        query_in_lib = 0
 75        for peak in query_mz_arr:
 76            if np.any(np.isclose(lib_mzs, peak, atol=mz_tol_da)):
 77                query_in_lib += 1
 78        query_in_lib_fract = query_in_lib / len(query_mz_arr)
 79
 80        # Get count and fraction of peaks in lib that are in query
 81        lib_in_query = 0
 82        for peak in lib_mzs:
 83            if np.any(np.isclose(query_mz_arr, peak, atol=mz_tol_da)):
 84                lib_in_query += 1
 85        lib_in_query_fract = lib_in_query / len(lib_mzs)
 86
 87        if include_fragment_types:
 88            # Check that fragment types are present in the library entry
 89            if "fragment_types" not in lib_entry.keys():
 90                raise ValueError(
 91                    "Flash entropy library entry must have 'fragment_types' key to include fragment types in output."
 92                )
 93
 94            # Get types of fragments in the lib entry and convert it to a list on ", "
 95            lib_frags = [x.strip() for x in lib_entry["fragment_types"].split(",")]
 96            # make list of the fragment types that are present in the query spectrum
 97            lib_in_query_ids = list(
 98                set(
 99                    [
100                        ind
101                        for ind, x in enumerate(lib_mzs)
102                        if len(np.where(np.isclose(query_mz_arr, x, atol=mz_tol_da))[0])
103                        > 0
104                    ]
105                )
106            )
107            query_frags = list(set([lib_frags[x] for x in lib_in_query_ids]))
108            lib_frags = list(set(lib_frags))
109
110        else:
111            query_frags = None
112            lib_frags = None
113
114        return (
115            query_in_lib,
116            query_in_lib_fract,
117            lib_in_query,
118            lib_in_query_fract,
119            query_frags,
120            lib_frags,
121        )

Return additional match qualities between query and library entry.

Parameters

query_mz_arr (np.array): Array of query spectrum. Shape (N, 2), with m/z in the first column and abundance in the second.
lib_entry (dict): Library spectrum entry, with 'mz' key containing the spectrum in the format (mz, abundance),(mz, abundance), i.e. from MetabRef.
mz_tol_da (float, optional): Tolerance in Da for matching peaks (in MS2). Default is 0.1.
include_fragment_types (bool, optional): If True, include fragment type comparisons in output. Defaults to False.

Returns

tuple: Tuple of (query_in_lib, query_in_lib_fract, lib_in_query, lib_in_query_fract, query_frags, lib_frags, lib_precursor_mz).

Notes

query_in_lib : int Number of peaks in query that are present in the library entry (within mz_tol_da). query_in_lib_fract : float Fraction of peaks in query that are present in the library entry (within mz_tol_da). lib_in_query : int Number of peaks in the library entry that are present in the query (within mz_tol_da). lib_in_query_fract : float Fraction of peaks in the library entry that are present in the query (within mz_tol_da). query_frags : list List of unique fragment types present in the query, generally 'MLF' or 'LSF' or both. lib_frags : list List of unique fragment types present in the library entry, generally 'MLF' or 'LSF' or both.

Raises

ValueError: If library entry does not have 'fragment_types' key and include_fragment_types is True.

def fe_search( self, scan_list, fe_lib, precursor_mz_list=[], use_mass_features=True, peak_sep_da=0.01, get_additional_metrics=True, accumulate_results=False): View Source

123    def fe_search(
124        self,
125        scan_list,
126        fe_lib,
127        precursor_mz_list=[],
128        use_mass_features=True,
129        peak_sep_da=0.01,
130        get_additional_metrics=True,
131        accumulate_results=False,
132    ):
133        """
134        Search LCMS spectra using a FlashEntropy approach.
135
136        Parameters
137        ----------
138        scan_list : list
139            List of scan numbers to search.
140        fe_lib : :obj:`~ms_entropy.FlashEntropySearch`
141            FlashEntropy Search instance.
142        precursor_mz_list : list, optional
143            List of precursor m/z values to search, by default [], which implies
144            matched with mass features; to enable this use_mass_features must be True.
145        use_mass_features : bool, optional
146            If True, use mass features to get precursor m/z values, by default True.
147            If True, will add search results to mass features' ms2_similarity_results attribute.
148        peak_sep_da : float, optional
149            Minimum separation between m/z peaks spectra in Da. This needs match the
150            approximate resolution of the search spectra and the FlashEntropySearch
151            instance, by default 0.01.
152        get_additional_metrics : bool, optional
153            If True, get additional metrics from FlashEntropy search, by default True.
154        accumulate_results : bool, optional
155            If True, accumulate results with existing spectral_search_results instead of
156            replacing them. This allows searching the same scans with multiple libraries
157            without overwriting previous results, by default False.
158
159        Returns
160        -------
161        None, but adds results to self.spectral_search_results and associates these
162        spectral_search_results with mass_features within the self.mass_features dictionary.
163
164        """
165        # Retrieve parameters from self
166        # include_fragment_types should used for lipids queries only, not general metabolomics
167        include_fragment_types = self.parameters.lc_ms.include_fragment_types
168        min_match_score = self.parameters.lc_ms.ms2_min_fe_score
169
170        # If precursor_mz_list is empty and use_mass_features is True, get precursor m/z values from mass features for each scan in scan_list
171        if use_mass_features and len(precursor_mz_list) == 0:
172            precursor_mz_list = []
173            for scan in scan_list:
174                mf_ids = [
175                    key
176                    for key, value in self.mass_features.items()
177                    if scan in value.ms2_mass_spectra
178                ]
179                precursor_mz = [
180                    value.mz
181                    for key, value in self.mass_features.items()
182                    if key in mf_ids
183                ]
184                precursor_mz_list.append(precursor_mz)
185
186        # Check that precursor_mz_list same length as scan_list, if not, raise error
187        if len(precursor_mz_list) != len(scan_list):
188            raise ValueError("Length of precursor_mz_list is not equal to scan_list.")
189
190        # Loop through each query spectrum / precursor match and save ids of db spectrum that are decent matches
191        overall_results_dict = {}
192        for i in np.arange(len(scan_list)):
193            scan_oi = scan_list[i]
194            if len(self._ms[scan_oi].mspeaks) > 0:
195                precursor_mzs = precursor_mz_list[i]
196                overall_results_dict[scan_oi] = {}
197                for precursor_mz in precursor_mzs:
198                    query_spectrum = fe_lib.clean_spectrum_for_search(
199                        precursor_mz=precursor_mz,
200                        peaks=np.vstack(
201                            (self._ms[scan_oi].mz_exp, self._ms[scan_oi].abundance)
202                        ).T,
203                        precursor_ions_removal_da=None,
204                        noise_threshold=self._ms[
205                            scan_oi
206                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
207                        / 100,
208                        min_ms2_difference_in_da=peak_sep_da,
209                    )
210                    search_results = fe_lib.search(
211                        precursor_mz=precursor_mz,
212                        peaks=query_spectrum,
213                        ms1_tolerance_in_da=self.parameters.mass_spectrum[
214                            "ms1"
215                        ].molecular_search.max_ppm_error
216                        * 10**-6
217                        * precursor_mz,
218                        ms2_tolerance_in_da=peak_sep_da * 0.5,
219                        method={"identity"},
220                        precursor_ions_removal_da=None,
221                        noise_threshold=self._ms[
222                            scan_oi
223                        ].parameters.mass_spectrum.noise_threshold_min_relative_abundance
224                        / 100,
225                        target="cpu",
226                    )["identity_search"]
227                    match_inds = np.where(search_results > min_match_score)[0]
228
229                    # If any decent matches are found, add them to the results dictionary
230                    if len(match_inds) > 0:
231                        match_scores = search_results[match_inds]
232                        ref_ms_ids = [fe_lib[x]["id"] for x in match_inds]
233                        ref_mol_ids = [
234                            fe_lib[x]["molecular_data_id"] for x in match_inds
235                        ]
236                        ref_precursor_mzs = [
237                            fe_lib[x]["precursor_mz"] for x in match_inds
238                        ]
239                        ion_types = [fe_lib[x]["ion_type"] for x in match_inds]
240                        overall_results_dict[scan_oi][precursor_mz] = {
241                            "ref_mol_id": ref_mol_ids,
242                            "ref_ms_id": ref_ms_ids,
243                            "ref_precursor_mz": ref_precursor_mzs,
244                            "precursor_mz_error_ppm": [
245                                (precursor_mz - x) / precursor_mz * 10**6
246                                for x in ref_precursor_mzs
247                            ],
248                            "entropy_similarity": match_scores,
249                            "ref_ion_type": ion_types,
250                        }
251                        # Add database name, if present
252                        db_name = [
253                            fe_lib[x].get("database_name") for x in match_inds
254                        ]
255                        if db_name is not None:
256                            overall_results_dict[scan_oi][precursor_mz].update(
257                                {"database_name": db_name}
258                            )
259                        if get_additional_metrics:
260                            more_match_quals = [
261                                self.get_more_match_quals(
262                                    self._ms[scan_oi].mz_exp,
263                                    fe_lib[x],
264                                    mz_tol_da=peak_sep_da,
265                                    include_fragment_types=include_fragment_types,
266                                )
267                                for x in match_inds
268                            ]
269                            overall_results_dict[scan_oi][precursor_mz].update(
270                                {
271                                    "query_mz_in_ref_n": [
272                                        x[0] for x in more_match_quals
273                                    ],
274                                    "query_mz_in_ref_fract": [
275                                        x[1] for x in more_match_quals
276                                    ],
277                                    "ref_mz_in_query_n": [
278                                        x[2] for x in more_match_quals
279                                    ],
280                                    "ref_mz_in_query_fract": [
281                                        x[3] for x in more_match_quals
282                                    ],
283                                }
284                            )
285                            if include_fragment_types:
286                                overall_results_dict[scan_oi][precursor_mz].update(
287                                    {
288                                        "query_frag_types": [
289                                            x[4] for x in more_match_quals
290                                        ],
291                                        "ref_frag_types": [
292                                            x[5] for x in more_match_quals
293                                        ],
294                                    }
295                                )
296
297        # Drop scans with no results from dictionary
298        overall_results_dict = {k: v for k, v in overall_results_dict.items() if v}
299
300        # Cast each entry as a MS2SearchResults object
301        for scan_id in overall_results_dict.keys():
302            for precursor_mz in overall_results_dict[scan_id].keys():
303                ms2_spectrum = self._ms[scan_id]
304                ms2_search_results = overall_results_dict[scan_id][precursor_mz]
305                overall_results_dict[scan_id][precursor_mz] = SpectrumSearchResults(
306                    ms2_spectrum, precursor_mz, ms2_search_results
307                )
308
309        # Add MS2SearchResults to the existing spectral search results dictionary
310        if accumulate_results:
311            # Merge results with existing spectral_search_results
312            for scan_id, precursor_dict in overall_results_dict.items():
313                if scan_id in self.spectral_search_results:
314                    # Scan already has results, merge precursor_mz dictionaries
315                    self.spectral_search_results[scan_id].update(precursor_dict)
316                else:
317                    # New scan, add entire dictionary
318                    self.spectral_search_results[scan_id] = precursor_dict
319        else:
320            # Replace existing results (original behavior)
321            self.spectral_search_results.update(overall_results_dict)
322
323        # If there are mass features, associate the results with each mass feature
324        if len(self.mass_features) > 0:
325            # Determine which results to associate with mass features
326            if accumulate_results:
327                # When accumulating, only associate new results from this search
328                # to avoid duplicating previously associated results
329                results_to_associate = overall_results_dict
330            else:
331                # When not accumulating, clear existing associations and re-associate all results
332                for mass_feature_id in self.mass_features.keys():
333                    self.mass_features[mass_feature_id].ms2_similarity_results = []
334                results_to_associate = self.spectral_search_results
335            
336            for mass_feature_id, mass_feature in self.mass_features.items():
337                scan_ids = mass_feature.ms2_scan_numbers
338                for ms2_scan_id in scan_ids:
339                    precursor_mz = mass_feature.mz
340                    try:
341                        results_to_associate[ms2_scan_id][precursor_mz]
342                    except KeyError:
343                        pass
344                    else:
345                        self.mass_features[
346                            mass_feature_id
347                        ].ms2_similarity_results.append(
348                            results_to_associate[ms2_scan_id][precursor_mz]
349                        )

Search LCMS spectra using a FlashEntropy approach.

Parameters

scan_list (list): List of scan numbers to search.
fe_lib (~ms_entropy.FlashEntropySearch): FlashEntropy Search instance.
precursor_mz_list (list, optional): List of precursor m/z values to search, by default [], which implies matched with mass features; to enable this use_mass_features must be True.
use_mass_features (bool, optional): If True, use mass features to get precursor m/z values, by default True. If True, will add search results to mass features' ms2_similarity_results attribute.
peak_sep_da (float, optional): Minimum separation between m/z peaks spectra in Da. This needs match the approximate resolution of the search spectra and the FlashEntropySearch instance, by default 0.01.
get_additional_metrics (bool, optional): If True, get additional metrics from FlashEntropy search, by default True.
accumulate_results (bool, optional): If True, accumulate results with existing spectral_search_results instead of replacing them. This allows searching the same scans with multiple libraries without overwriting previous results, by default False.

Returns

None, but adds results to self.spectral_search_results and associates these
spectral_search_results with mass_features within the self.mass_features dictionary.