corems.mass_spectrum.calc.Calibration

  1# -*- coding: utf-8 -*-
  2"""
  3Created on Wed May 13 02:16:09 2020
  4
  5@author: Will Kew
  6"""
  7
  8# import modules
  9import csv
 10import warnings
 11from io import BytesIO, StringIO
 12from pathlib import Path
 13
 14import numpy as np
 15import pandas as pd
 16from s3path import S3Path
 17
 18# import scipy modules for calibration
 19from scipy.optimize import minimize
 20
 21
 22class MzDomainCalibration:
 23    """MzDomainCalibration class for recalibrating mass spectra
 24
 25    Parameters
 26    ----------
 27    mass_spectrum : CoreMS MassSpectrum Object
 28        The mass spectrum to be calibrated.
 29    ref_masslist : str
 30        The path to a reference mass list.
 31    mzsegment : tuple of floats, optional
 32        The mz range to recalibrate, or None. Used for calibration of specific parts of the mz domain at a time.
 33        Future work - allow multiple mzsegments to be passed.
 34
 35    Attributes
 36    ----------
 37    mass_spectrum : CoreMS MassSpectrum Object
 38        The mass spectrum to be calibrated.
 39    mzsegment : tuple of floats or None
 40        The mz range to recalibrate, or None.
 41    ref_mass_list_path : str or Path
 42        The path to the reference mass list.
 43
 44    Methods
 45    -------
 46    * run().
 47        Main function to run this class.
 48    * load_ref_mass_list().
 49        Load reference mass list (Bruker format).
 50    * gen_ref_mass_list_from_assigned(min_conf=0.7).
 51        Generate reference mass list from assigned masses.
 52    * find_calibration_points(df_ref, calib_ppm_error_threshold=(-1, 1), calib_snr_threshold=5).
 53        Find calibration points in the mass spectrum based on the reference mass list.
 54    * robust_calib(param, cal_peaks_mz, cal_refs_mz, order=1).
 55        Recalibration function.
 56    * recalibrate_mass_spectrum(cal_peaks_mz, cal_refs_mz, order=1, diagnostic=False).
 57        Main recalibration function which uses a robust linear regression.
 58
 59
 60    """
 61
 62    def __init__(self, mass_spectrum, ref_masslist, mzsegment=None):
 63        self.mass_spectrum = mass_spectrum
 64        self.mzsegment = mzsegment
 65
 66        # define reference mass list - bruker .ref format
 67        self.ref_mass_list_path = ref_masslist
 68        if self.mass_spectrum.percentile_assigned(mute_output=True)[0] != 0:
 69            warnings.warn(
 70                "Warning: calibrating spectra which have already been assigned may yield erroneous results"
 71            )
 72        self.mass_spectrum.mz_cal = self.mass_spectrum.mz_exp
 73        self.mass_spectrum.mz_cal_profile = self.mass_spectrum._mz_exp
 74
 75        if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
 76            print(
 77                "MS Obj loaded - " + str(len(mass_spectrum.mspeaks)) + " peaks found."
 78            )
 79
 80    def load_ref_mass_list(self):
 81        """
 82        Load reference mass list (Bruker format).
 83
 84        Loads in a reference mass list from a .ref file. Some versions of
 85        Bruker's software produce .ref files with a different format where
 86        the header lines (starting with '#' or '##') and delimiters may vary.
 87        The file may be located locally or on S3 and will be handled accordingly.
 88
 89        Returns
 90        -------
 91        df_ref : Pandas DataFrame
 92            Reference mass list object.
 93        """
 94        # Get a Path-like object from the input path string or S3Path
 95        refmasslist = (
 96            Path(self.ref_mass_list_path)
 97            if isinstance(self.ref_mass_list_path, str)
 98            else self.ref_mass_list_path
 99        )
100
101        # Make sure the file exists
102        if not refmasslist.exists():
103            raise FileExistsError("File does not exist: %s" % refmasslist)
104
105        # Read all lines from the file (handling S3 vs local differently)
106        if isinstance(refmasslist, S3Path):
107            # For S3, read the file in binary, then decode to string and split into lines.
108            content = refmasslist.open("rb").read()
109            all_lines = content.decode("utf-8").splitlines(keepends=True)
110        else:
111            # For a local file, open in text mode and read lines.
112            with refmasslist.open("r") as f:
113                all_lines = f.readlines()
114
115        # Identify the index of the first line of the actual data.
116        # We assume header lines start with '#' (or '##') and ignore blank lines.
117        data_start_idx = 0
118        for idx, line in enumerate(all_lines):
119            if line.strip() and not line.lstrip().startswith("#"):
120                data_start_idx = idx
121                break
122
123        # If there are not enough lines to guess the dialect, throw an error
124        if data_start_idx >= len(all_lines):
125            raise ValueError("The file does not appear to contain any data lines.")
126
127        # Use a couple of the data lines to let csv.Sniffer detect the delimiter
128        sample_lines = "".join(all_lines[data_start_idx : data_start_idx + 2])
129        try:
130            dialect = csv.Sniffer().sniff(sample_lines)
131            delimiter = dialect.delimiter
132        except csv.Error:
133            # If csv.Sniffer fails, default to a common delimiter (e.g., comma)
134            delimiter = ","
135
136        # Join the lines from the beginning of data (this might include a blank line)
137        joined_data = "".join(all_lines[data_start_idx:])
138
139        # Depending on whether the file is S3 or local, wrap the data as needed for pandas
140        if isinstance(refmasslist, S3Path):
141            data_buffer = BytesIO(joined_data.encode("utf-8"))
142        else:
143            data_buffer = StringIO(joined_data)
144
145        # Read data into a DataFrame.
146        # Adjust columns and names as needed – here we assume at least 2 columns:
147        df_ref = pd.read_csv(
148            data_buffer,
149            sep=delimiter,
150            header=None,
151            usecols=[0, 1],  # Modify if more columns are required.
152            names=["Formula", "m/z"],
153        )
154
155        df_ref.sort_values(by="m/z", ascending=True, inplace=True)
156
157        if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
158            print(
159                "Reference mass list loaded - {} calibration masses loaded.".format(
160                    len(df_ref)
161                )
162            )
163
164        return df_ref
165
166    def gen_ref_mass_list_from_assigned(self, min_conf: float = 0.7):
167        """Generate reference mass list from assigned masses
168
169        This function will generate a ref mass dataframe object from an assigned corems mass spec obj
170        using assigned masses above a certain minimum confidence threshold.
171
172        This function needs to be retested and check it is covered in the unit tests.
173
174        Parameters
175        ----------
176        min_conf : float, optional
177            minimum confidence score. The default is 0.7.
178
179        Returns
180        -------
181        df_ref : Pandas DataFrame
182            reference mass list - based on calculated masses.
183
184        """
185        # TODO this function needs to be retested and check it is covered in the unit tests
186        df = self.mass_spectrum.to_dataframe()
187        df = df[df["Confidence Score"] > min_conf]
188        df_ref = pd.DataFrame(columns=["m/z"])
189        df_ref["m/z"] = df["Calculated m/z"]
190        if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
191            print(
192                "Reference mass list generated - "
193                + str(len(df_ref))
194                + " calibration masses."
195            )
196        return df_ref
197
198    def find_calibration_points(
199        self,
200        df_ref,
201        calib_ppm_error_threshold: tuple[float, float] = (-1, 1),
202        calib_snr_threshold: float = 5,
203        calibration_ref_match_method: str = "merged",
204        calibration_ref_match_tolerance: float = 0.003,
205        calibration_ref_match_std_raw_error_limit: float = 1.5,
206    ):
207        """Function to find calibration points in the mass spectrum
208
209        Based on the reference mass list.
210
211        Parameters
212        ----------
213        df_ref : Pandas DataFrame
214            reference mass list for recalibration.
215        calib_ppm_error_threshold : tuple of floats, optional
216            ppm error for finding calibration masses in the spectrum. The default is -1,1.
217            Note: This is based on the calculation of ppm = ((mz_measure - mz_theoretical)/mz_theoretical)*1e6.
218                Some software does this the other way around and value signs must be inverted for that to work.
219        calib_snr_threshold : float, optional
220            snr threshold for finding calibration masses in the spectrum. The default is 5.
221            If SNR data is unavailable, peaks are filtered by intensity percentile using the formula:
222            percentile = max(5, 100 - calib_snr_threshold)
223        calibration_ref_match_method : str, optional
224            method for matching calibration references. The default is "merged".
225        calibration_ref_match_tolerance : float, optional
226            tolerance for matching calibration references. The default is 0.003.
227        calibration_ref_match_std_raw_error_limit : float, optional
228            standard deviation raw error limit for calibration references. The default is 1.5.
229
230        Returns
231        -------
232        cal_peaks_mz : list of floats
233            masses of measured ions to use in calibration routine
234        cal_refs_mz : list of floats
235            reference mz values of found calibration points.
236
237        """
238
239        # Check if SNR data is available by testing the first peak
240        use_snr = False
241        if len(self.mass_spectrum.mspeaks) > 0:
242            first_peak = self.mass_spectrum.mspeaks[0]
243            if (
244                hasattr(first_peak, "signal_to_noise")
245                and first_peak.signal_to_noise is not None
246                and not np.isnan(first_peak.signal_to_noise)
247                and first_peak.signal_to_noise > 0
248            ):
249                use_snr = True
250
251        # This approach is much more efficient and expedient than the original implementation.
252        peaks_mz = []
253        peaks_intensity = []
254
255        if use_snr:
256            # Use SNR filtering
257            for x in self.mass_spectrum.mspeaks:
258                if x.signal_to_noise > calib_snr_threshold:
259                    if self.mzsegment:
260                        if min(self.mzsegment) <= x.mz_exp <= max(self.mzsegment):
261                            peaks_mz.append(x.mz_exp)
262                    else:
263                        peaks_mz.append(x.mz_exp)
264        else:
265            # Fallback to intensity percentile filtering
266            intensity_percentile = max(5, 100 - calib_snr_threshold)
267            warnings.warn(
268                f"SNR data unavailable for calibration. Using intensity-based filtering instead. "
269                f"SNR threshold of {calib_snr_threshold} corresponds to intensity percentile >= {intensity_percentile}%."
270            )
271
272            # Collect all peaks and their intensities
273            all_peaks_data = []
274            for x in self.mass_spectrum.mspeaks:
275                if self.mzsegment:
276                    if min(self.mzsegment) <= x.mz_exp <= max(self.mzsegment):
277                        all_peaks_data.append((x.mz_exp, x.abundance))
278                else:
279                    all_peaks_data.append((x.mz_exp, x.abundance))
280
281            if all_peaks_data:
282                peaks_mz_list, intensities = zip(*all_peaks_data)
283                intensity_threshold = np.percentile(intensities, intensity_percentile)
284
285                for mz, intensity in all_peaks_data:
286                    if intensity >= intensity_threshold:
287                        peaks_mz.append(mz)
288
289        peaks_mz = np.asarray(peaks_mz)
290
291        if calibration_ref_match_method == "legacy":
292            # This legacy approach iterates through each reference match and finds the entries within 1 mz and within the user defined PPM error threshold
293            # Then it removes ambiguities - which means the calibration threshold hasto be very tight.
294            cal_peaks_mz = []
295            cal_refs_mz = []
296            for mzref in df_ref["m/z"]:
297                tmp_peaks_mz = peaks_mz[abs(peaks_mz - mzref) < 1]
298                for mzmeas in tmp_peaks_mz:
299                    delta_mass = ((mzmeas - mzref) / mzref) * 1e6
300                    if delta_mass < max(calib_ppm_error_threshold):
301                        if delta_mass > min(calib_ppm_error_threshold):
302                            cal_peaks_mz.append(mzmeas)
303                            cal_refs_mz.append(mzref)
304
305            # To remove entries with duplicated indices (reference masses matching multiple peaks)
306            tmpdf = pd.Series(index=cal_refs_mz, data=cal_peaks_mz, dtype=float)
307            tmpdf = tmpdf[~tmpdf.index.duplicated(keep=False)]
308
309            cal_peaks_mz = list(tmpdf.to_numpy())
310            cal_refs_mz = list(tmpdf.index)
311        elif calibration_ref_match_method == "merged":
312            # This is a new approach (August 2024) which uses Pandas 'merged_asof' to find the peaks closest in m/z between
313            # reference and measured masses. This is a quicker way to match, and seems to get more matches.
314            # It may not work as well when the data are far from correc initial mass
315            # e.g. if the correct peak is further from the reference than an incorrect peak.
316            meas_df = pd.DataFrame(columns=["meas_m/z"], data=peaks_mz)
317            tolerance = calibration_ref_match_tolerance
318            merged_df = pd.merge_asof(
319                df_ref,
320                meas_df,
321                left_on="m/z",
322                right_on="meas_m/z",
323                tolerance=tolerance,
324                direction="nearest",
325            )
326            merged_df.dropna(how="any", inplace=True)
327            merged_df["Error_ppm"] = (
328                (merged_df["meas_m/z"] - merged_df["m/z"]) / merged_df["m/z"]
329            ) * 1e6
330            median_raw_error = merged_df["Error_ppm"].median()
331            std_raw_error = merged_df["Error_ppm"].std()
332            if std_raw_error > calibration_ref_match_std_raw_error_limit:
333                std_raw_error = calibration_ref_match_std_raw_error_limit
334            self.mass_spectrum.calibration_raw_error_median = median_raw_error
335            self.mass_spectrum.calibration_raw_error_stdev = std_raw_error
336            merged_df = merged_df[
337                (merged_df["Error_ppm"] > (median_raw_error - 1.5 * std_raw_error))
338                & (merged_df["Error_ppm"] < (median_raw_error + 1.5 * std_raw_error))
339            ]
340            # merged_df= merged_df[(merged_df['Error_ppm']>min(calib_ppm_error_threshold))&(merged_df['Error_ppm']<max(calib_ppm_error_threshold))]
341            cal_peaks_mz = list(merged_df["meas_m/z"])
342            cal_refs_mz = list(merged_df["m/z"])
343        else:
344            raise ValueError(f"{calibration_ref_match_method} not allowed.")
345
346        # it is crucial the mass lists are in same order
347        # corems likes to do masses from high to low.
348        cal_refs_mz.sort(reverse=False)
349        cal_peaks_mz.sort(reverse=False)
350        if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
351            print(
352                str(len(cal_peaks_mz))
353                + " calibration points matched within thresholds."
354            )
355        return cal_peaks_mz, cal_refs_mz
356
357    def robust_calib(
358        self,
359        param: list[float],
360        cal_peaks_mz: list[float],
361        cal_refs_mz: list[float],
362        order: int = 1,
363    ):
364        """Recalibration function
365
366        Computes the rms of m/z errors to minimize when calibrating.
367        This is adapted from from spike.
368
369        Parameters
370        ----------
371        param : list of floats
372            generated by minimize function from scipy optimize.
373        cal_peaks_mz : list of floats
374            masses of measured peaks to use in mass calibration.
375        cal_peaks_mz : list of floats
376            reference mz values of found calibration points.
377        order : int, optional
378            order of the recalibration function. 1 = linear, 2 = quadratic. The default is 1.
379
380        Returns
381        -------
382        rmserror : float
383            root mean square mass error for calibration points.
384
385        """
386        Aterm = param[0]
387        Bterm = param[1]
388        try:
389            Cterm = param[2]
390        except IndexError:
391            pass
392
393        # get the mspeaks from the mass spectrum object which were calibration points
394        # mspeaks = [self.mass_spectrum.mspeaks[x] for x in imzmeas]
395        # get their calibrated mass values
396        # mspeakmzs = [x.mz_cal for x in mspeaks]
397        cal_peaks_mz = np.asarray(cal_peaks_mz)
398
399        # linearz
400        if order == 1:
401            ref_recal_points = (Aterm * cal_peaks_mz) + Bterm
402        # quadratic
403        elif order == 2:
404            ref_recal_points = (Aterm * (cal_peaks_mz)) + (
405                Bterm * np.power((cal_peaks_mz), 2) + Cterm
406            )
407
408        # sort both the calibration points (measured, recalibrated)
409        ref_recal_points.sort()
410        # and sort the calibration points (theoretical, predefined)
411        cal_refs_mz.sort()
412
413        # calculate the ppm error for each calibration point
414        error = ((ref_recal_points - cal_refs_mz) / cal_refs_mz) * 1e6
415        # calculate the root mean square error - this is our target to minimize
416        rmserror = np.sqrt(np.mean(error**2))
417        return rmserror
418
419    def recalibrate_mass_spectrum(
420        self,
421        cal_peaks_mz: list[float],
422        cal_refs_mz: list[float],
423        order: int = 1,
424        diagnostic: bool = False,
425    ):
426        """Main recalibration function which uses a robust linear regression
427
428        This function performs the recalibration of the mass spectrum object.
429        It iteratively applies
430
431        Parameters
432        ----------
433        cal_peaks_mz : list of float
434            masses of measured peaks to use in mass calibration.
435        cal_refs_mz : list of float
436            reference mz values of found calibration points.
437        order : int, optional
438            order of the recalibration function. 1 = linear, 2 = quadratic. The default is 1.
439
440        Returns
441        -------
442        mass_spectrum : CoreMS mass spectrum object
443            Calibrated mass spectrum object
444
445
446        Notes
447        -----
448        This function is adapted, in part, from the SPIKE project [1,2] and is based on the robust linear regression method.
449
450        References
451        ----------
452        1.  Chiron L., Coutouly M-A., Starck J-P., Rolando C., Delsuc M-A.
453            SPIKE a Processing Software dedicated to Fourier Spectroscopies
454            https://arxiv.org/abs/1608.06777 (2016)
455        2.  SPIKE - https://github.com/spike-project/spike
456
457        """
458        # initialise parameters for recalibration
459        # these are the 'Aterm, Bterm, Cterm'
460        # as spectra are already freq->mz calibrated, these terms are very small
461        # may be beneficial to formally separate them from the freq->mz terms
462        if order == 1:
463            Po = [1, 0]
464        elif order == 2:
465            Po = [1, 0, 0]
466
467        if len(cal_peaks_mz) >= 2:
468            if self.mzsegment:  # If only part of the spectrum is to be recalibrated
469                mz_exp_peaks = np.array(
470                    [mspeak.mz_exp for mspeak in self.mass_spectrum]
471                )
472                # Split the array into two parts - one to recailbrate, one to keep unchanged.
473                mz_exp_peaks_tocal = mz_exp_peaks[
474                    (mz_exp_peaks >= min(self.mzsegment))
475                    & (mz_exp_peaks <= max(self.mzsegment))
476                ]
477                mz_exp_peaks_unchanged = mz_exp_peaks[
478                    ~(mz_exp_peaks >= min(self.mzsegment))
479                    | ~(mz_exp_peaks <= max(self.mzsegment))
480                ]
481                # TODO: - segmented calibration needs a way to better track the calibration args/values...
482                if not self.mass_spectrum.is_centroid:
483                    mz_exp_profile = np.array(self.mass_spectrum.mz_exp_profile)
484                    # Split the array into two parts - one to recailbrate, one to keep unchanged.
485                    mz_exp_profile_tocal = mz_exp_profile[
486                        (mz_exp_profile >= min(self.mzsegment))
487                        & (mz_exp_profile <= max(self.mzsegment))
488                    ]
489                    mz_exp_profile_unchanged = mz_exp_profile[
490                        ~(mz_exp_profile >= min(self.mzsegment))
491                        | ~(mz_exp_profile <= max(self.mzsegment))
492                    ]
493            else:  # if just recalibrating the whole spectrum
494                mz_exp_peaks_tocal = np.array(
495                    [mspeak.mz_exp for mspeak in self.mass_spectrum]
496                )
497                if not self.mass_spectrum.is_centroid:
498                    mz_exp_profile_tocal = np.array(self.mass_spectrum.mz_exp_profile)
499
500            minimize_method = self.mass_spectrum.settings.calib_minimize_method
501            res = minimize(
502                self.robust_calib,
503                Po,
504                args=(cal_peaks_mz, cal_refs_mz, order),
505                method=minimize_method,
506            )
507            if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
508                print(
509                    "minimize function completed with RMS error of: {:0.3f} ppm".format(
510                        res["fun"]
511                    )
512                )
513                print(
514                    "minimize function performed {:1d} fn evals and {:1d} iterations".format(
515                        res["nfev"], res["nit"]
516                    )
517                )
518            Pn = res.x
519
520            # mz_exp_ms = np.array([mspeak.mz_exp for mspeak in self.mass_spectrum])
521
522            if order == 1:
523                mz_domain = (Pn[0] * mz_exp_peaks_tocal) + Pn[1]
524                if not self.mass_spectrum.is_centroid:
525                    mz_profile_calc = (Pn[0] * mz_exp_profile_tocal) + Pn[1]
526
527            elif order == 2:
528                mz_domain = (Pn[0] * (mz_exp_peaks_tocal)) + (
529                    Pn[1] * np.power((mz_exp_peaks_tocal), 2) + Pn[2]
530                )
531
532                if not self.mass_spectrum.is_centroid:
533                    mz_profile_calc = (Pn[0] * (mz_exp_profile_tocal)) + (
534                        Pn[1] * np.power((mz_exp_profile_tocal), 2) + Pn[2]
535                    )
536
537            if self.mzsegment:
538                # Recombine the mass domains
539                mz_domain = np.concatenate([mz_domain, mz_exp_peaks_unchanged])
540                mz_domain.sort()
541                if not self.mass_spectrum.is_centroid:
542                    mz_profile_calc = np.concatenate(
543                        [mz_profile_calc, mz_exp_profile_unchanged]
544                    )
545                    mz_profile_calc.sort()
546                # Sort them
547                if (
548                    mz_exp_peaks[0] > mz_exp_peaks[1]
549                ):  # If originally descending mass order
550                    mz_domain = mz_domain[::-1]
551                    if not self.mass_spectrum.is_centroid:
552                        mz_profile_calc = mz_profile_calc[::-1]
553
554            self.mass_spectrum.mz_cal = mz_domain
555            if not self.mass_spectrum.is_centroid:
556                self.mass_spectrum.mz_cal_profile = mz_profile_calc
557
558            self.mass_spectrum.calibration_order = order
559            self.mass_spectrum.calibration_RMS = float(res["fun"])
560            self.mass_spectrum.calibration_points = int(len(cal_refs_mz))
561            self.mass_spectrum.calibration_ref_mzs = cal_refs_mz
562            self.mass_spectrum.calibration_meas_mzs = cal_peaks_mz
563
564            self.mass_spectrum.calibration_segment = self.mzsegment
565
566            if diagnostic:
567                return self.mass_spectrum, res
568            return self.mass_spectrum
569        else:
570            warnings.warn("Too few calibration points - aborting.")
571            return self.mass_spectrum
572
573    def run(self):
574        """Run the calibration routine
575
576        This function runs the calibration routine.
577
578        """
579        calib_snr_threshold = self.mass_spectrum.settings.calib_sn_threshold
580        max_calib_ppm_error = self.mass_spectrum.settings.max_calib_ppm_error
581        min_calib_ppm_error = self.mass_spectrum.settings.min_calib_ppm_error
582        calib_pol_order = self.mass_spectrum.settings.calib_pol_order
583        calibration_ref_match_method = (
584            self.mass_spectrum.settings.calibration_ref_match_method
585        )
586        calibration_ref_match_tolerance = (
587            self.mass_spectrum.settings.calibration_ref_match_tolerance
588        )
589        calibration_ref_match_std_raw_error_limit = (
590            self.mass_spectrum.settings.calibration_ref_match_std_raw_error_limit
591        )
592
593        # load reference mass list
594        df_ref = self.load_ref_mass_list()
595
596        # find calibration points
597        cal_peaks_mz, cal_refs_mz = self.find_calibration_points(
598            df_ref,
599            calib_ppm_error_threshold=(min_calib_ppm_error, max_calib_ppm_error),
600            calib_snr_threshold=calib_snr_threshold,
601            calibration_ref_match_method=calibration_ref_match_method,
602            calibration_ref_match_tolerance=calibration_ref_match_tolerance,
603            calibration_ref_match_std_raw_error_limit=calibration_ref_match_std_raw_error_limit,
604        )
605        if len(cal_peaks_mz) == 2:
606            self.mass_spectrum.settings.calib_pol_order = 1
607            calib_pol_order = 1
608            if self.mass_spectrum.parameters.mass_spectrum.verbose_processing:
609                print("Only 2 calibration points found, forcing a linear recalibration")
610        elif len(cal_peaks_mz) < 2:
611            warnings.warn("Too few calibration points found, function will fail")
612        self.recalibrate_mass_spectrum(cal_peaks_mz, cal_refs_mz, order=calib_pol_order)