Skip to contents

Detect and correct isolated spikes in spectra

Source: R/correct_spike.R
correct_spike.Rd

correct_spike() detects isolated positive or negative intensity artifacts and replaces only accepted spike intervals with local interpolation. The default MAD-prominence-width method detects narrow positive and negative peaks, estimates noise from the raw median absolute deviation of first differences, and replaces detected intervals from nearby clean values. The legacy residual method compares each point with a wavenumber-aware local interpolation and scales the residual by a local median absolute deviation. Two additional paper-backed methods use peak prominence and width measured in sample (CCD-pixel) units.

Usage

correct_spike(x, ...)

# Default S3 method
correct_spike(x, ...)

# S3 method for class 'OpenSpecy'
correct_spike(
  x,
  method = c("mad_prominence_width", "residual", "prominence_fwhm",
    "prominence_fwhm_ratio"),
  direction = c("both", "positive", "negative"),
  residual_window = 5L,
  residual_threshold = 8,
  residual_max_width = 1L,
  prominence_threshold = NULL,
  width_threshold = NULL,
  noise_multiplier = 10,
  rel_height = 0.8,
  interpolation_points = 5L,
  interpolation = c("linear", "quadratic"),
  z_threshold = 3.5,
  min_peaks = 20L,
  ...
)

Arguments

x

an OpenSpecy object.

method

character; detection method. One of "mad_prominence_width" (default), "residual", "prominence_fwhm", or "prominence_fwhm_ratio".

direction

character; detect "both" positive and negative spikes, only "positive" spikes, or only "negative" spikes.

residual_window

positive integer; points on each side used by the local residual predictor.

residual_threshold

positive numeric; absolute robust residual score required by the residual method.

residual_max_width

positive integer; widest consecutive candidate interval accepted by the residual method. The one-point default is deliberately conservative.

prominence_threshold

positive numeric or NULL; minimum peak prominence. NULL uses noise_multiplier times the raw MAD of first differences for the default method; it remains required for the manual prominence/FWHM method.

width_threshold

positive numeric or NULL; maximum peak FWHM in sample (CCD-pixel) units for the manual prominence/FWHM method. NULL uses 2 points for the default method.

noise_multiplier

positive numeric multiplier for the automatic raw MAD prominence threshold used by method = "mad_prominence_width".

rel_height

numeric in (0, 1]; prominence fraction at which the interval replaced by paper methods is measured. Coca-Lopez used 0.8 for most examples.

interpolation_points

positive integer; finite, unflagged neighboring points required on each side of an accepted interval. This is the paper's m parameter.

interpolation

character; "linear" (default) or "quadratic" local interpolation.

z_threshold

positive numeric; upper Z-score threshold for automated prominence/FWHM-ratio detection. The paper uses values greater than 3.5.

min_peaks

integer of at least two; minimum number of measurable peaks used to estimate automated ratio outliers.

...

must be empty. Unexpected arguments are rejected so detector tuning misspellings cannot be silently ignored.

Value

An OpenSpecy object with accepted spike intervals corrected. The wavenumber axis, spectra dimensions and names, metadata alignment, and existing attributes are preserved. A successful or rejected attempted correction stores an automatic_spike attribute containing the method, parameters, corrected and rejected regions, affected spectra, detector counts, pass count, and transaction reason. If nothing is detected, x is returned unchanged.

Details

method = "mad_prominence_width" reproduces the automated workflow supplied by Nicolas Coca Lopez without requiring pracma. With the defaults, peaks must span no more than 2 points and exceed 10 times the raw MAD of first differences. prominence_threshold may override that automatic threshold. Detected intervals receive a one-point guard on either side and use a 5-point local interpolation window; at an edge the nearest clean value is used rather than extrapolating a slope.

method = "prominence_fwhm" requires user-supplied prominence_threshold and width_threshold values. These thresholds depend on the material, instrument, spectral resolution, and acquisition settings; the graphene values reported by Coca-Lopez are deliberately not universal defaults. method = "prominence_fwhm_ratio" instead treats prominence/FWHM values above z_threshold standard deviations as spikes and requires at least min_peaks measurable peaks.

Peak widths and flagged intervals follow the prominence contour definition used by scipy.signal.peak_widths(): FWHM is measured at rel_height = 0.5, while the interval replaced is measured at the requested rel_height. The paper-backed modes require interpolation_points finite, unflagged samples on both sides; boundary values are never wrapped. The default method searches that many points on either side and uses the nearest clean value when only one side exists. Close spike intervals are merged before interpolation so one spike cannot be used to repair another. Linear interpolation that materially disagrees with a local quadratic reconstruction over a multi-point interval is rejected to avoid silently truncating an underlying broad band.

Correction proceeds through bounded transactional passes while the detector's correctable count strictly decreases. This lets a newly revealed spike be corrected without rolling back safe earlier replacements. Processing stops on no progress; boundary, interpolation, and band-protection safeguards stay in force, and any remaining safeguarded candidates are recorded rather than forced.

No single-spectrum method can always distinguish a cosmic-ray spike from a genuine band with the same shape. Calibrate paper thresholds on representative standards, especially for narrow-band materials such as calcite and polystyrene, and inspect the automatic_spike diagnostic attribute.

References

Coca-Lopez N (2024). "An intuitive approach for spike removal in Raman spectra based on peaks' prominence and width." Analytica Chimica Acta, 1295, 342312. doi:10.1016/j.aca.2024.342312 .

Author

Win Cowger, Nicolas Coca Lopez

Examples

wave <- seq(400, 1800, length.out = 101)
values <- sin(wave / 200)
values[51] <- values[51] + 20
spectrum <- as_OpenSpecy(wave, data.frame(sample = values))
corrected <- correct_spike(spectrum)