Peaks.EghPeeler

Peaks.EghPeeler#

Peaks.EghPeeler

Sequential EGH (Exponentially-modified Gaussian Hybrid) peeling for automatic peak recognition and component count estimation.

Algorithm:

  1. Find tallest remaining peak (argmax of smoothed residual)

  2. Fit EGH(H, mu, sigma, tau) via L-BFGS-B

  3. Check shape: min_sigma <= sigma <= max_sigma (see below)

  4. Check significance: area(fitted) / total_area >= min_area_frac

  5. Subtract fitted EGH from residual

  6. Repeat until no significant candidate remains

Physical basis for max_sigma_ratio#

The theoretical plate number N = (tR / sigma)^2 is a column property, approximately constant across all peaks in the same run. This means sigma = tR / sqrt(N), i.e. sigma is proportional to retention time tR.

Our data is in frame space, where frame numbers do not start from injection. However, the unknown offset f_offset >= 0, so the most permissive bound (f_offset = 0) gives:

sigma_k / sigma_1 <= mu_k / mu_1

Since all peaks lie in a narrow frame window, mu_k / mu_1 ~ 1, and any positive f_offset only tightens the bound. Therefore all real peaks must have sigma within a small multiple of the dominant peak’s sigma. The default ratio of 2.0 includes a generous safety margin; empirically, legitimate peaks fall within 1.5x.

Reference:

https://www.shimadzu.com/an/service-support/technical-support/ analysis-basics/basic/theoretical_plate.html

This replaces the legacy recognize_peaks (greedy sequential subtraction with symmetric Gaussian initialization) from molass-legacy.

egh_peel(x, y, num_components=None, min_area_frac=0.02, min_sigma=3.0, max_sigma_ratio=2.0, min_height_frac=0.05, debug=False)#

Sequential EGH peeling to find peaks and estimate component count.

Parameters:

  • x (array-like) – Frame positions (1-D).

  • y (array-like) – Intensity values (1-D, same length as x).

  • num_components (int or None) – If given, peel exactly this many peaks (ignoring area significance). If None, peel until no significant candidate remains.

  • min_area_frac (float) – Minimum area fraction (relative to total curve area) for a peak to be considered significant. Only used when num_components is None. Default 0.02 (2 %).

  • min_sigma (float) – Minimum allowed sigma. Peaks narrower than this are rejected as noise spikes. Default 3.0 frames.

  • max_sigma_ratio (float or None) – Maximum allowed sigma as a multiple of the first (dominant) peak’s sigma. Peaks wider than max_sigma_ratio * sigma_dominant are rejected as physically implausible (plate number constraint). Default 2.0. Set to None to disable.

  • min_height_frac (float) – Minimum allowed height as a fraction of the dominant peak’s height. Candidate peaks shorter than min_height_frac * H_dominant are rejected as ghost components. Only used when num_components is None. Default 0.05 (5 %).

  • debug (bool) – If True, print diagnostic messages.

Returns:

peak_list – Each element is [H, mu, sigma, tau] — the fitted EGH parameters. Sorted by ascending mu (retention time order).

Return type:

list of list