# Bio Applied Lc Ms Preprocessing

> Preprocess raw LC-MS mzML with XCMS centWave peak picking, obiwarp RT alignment, gap filling, PQN/QC normalization, adduct grouping. Use when building an XCMS pipeline or preprocessing untargeted metabolomics runs.

- Skill: `pavel-kravchenko/bio-applied-lc-ms-preprocessing` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/bio-applied-lc-ms-preprocessing`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/bio-applied-lc-ms-preprocessing/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: DevOps & Infra
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/bio-applied-lc-ms-preprocessing

---


# LC-MS Metabolomics Data Preprocessing

## When to Use

- Converting a batch of raw or centroided mzML files into an aligned features x samples intensity matrix
- Picking chromatographic peaks (centWave) and correcting retention-time drift across injections (obiwarp)
- Normalizing untargeted metabolomics data with pooled QC injections (PQN, LOESS batch correction, CV filtering)
- Grouping adduct/isotope features ([M+H]+, [M+Na]+, [M+NH4]+) that belong to the same metabolite
- Preparing a feature table for downstream PCA, OPLS-DA, or differential abundance testing

## Version Compatibility

- R >= 4.3, xcms >= 4.0 (Bioconductor >= 3.18), MSnbase >= 2.28, CAMERA >= 1.58
- Python >= 3.10, pyopenms >= 3.1, pandas >= 2.0, scikit-learn >= 1.4
- Raw data in mzML (centroided); vendor `.raw`/`.d`/`.wiff` must be converted first (ProteoWizard `msconvert`)

## Prerequisites

- `BiocManager::install(c("xcms", "CAMERA", "MSnbase"))` in R; `pip install pyopenms pandas scikit-learn` in Python
- Concepts: m/z, retention time (RT), ESI polarity, MS1 vs MS/MS, pooled QC-sample design
- Related skills: `bio-applied-metabolite-identification` (spectral matching once you have a feature table), `bio-applied-proteomics` (shared MS instrumentation concepts)

**Goal:** Turn a directory of mzML files into an aligned, gap-filled feature table (rows = features, columns = samples).
**Approach:** XCMS3 pipeline — `readMSData` → `findChromPeaks` (centWave) → `adjustRtime` (obiwarp) → `groupChromPeaks` (peak density correspondence) → `fillChromPeaks`.

```r
library(xcms)
library(MSnbase)

#' Run the standard XCMS preprocessing pipeline on a set of mzML files
#'
#' @param mzml_files character vector of paths to centroided mzML files
#' @param sample_group character vector of group labels, same length/order as mzml_files
#' @return an XcmsExperiment with peaks picked, RT-aligned, and grouped across samples
run_xcms_pipeline <- function(mzml_files, sample_group) {
  pheno <- data.frame(
    sample_name  = sub("\\.mzML$", "", basename(mzml_files)),
    sample_group = sample_group
  )
  raw_data <- readMSData(mzml_files, pdata = new("NAnnotatedDataFrame", pheno), mode = "onDisk")

  # 1. Peak picking: centWave detects chromatographic peaks in profile-mode EIC traces
  cwp <- CentWaveParam(ppm = 15, peakwidth = c(5, 30), snthresh = 10, prefilter = c(3, 1000))
  xdata <- findChromPeaks(raw_data, param = cwp)

  # 2. Retention-time alignment: obiwarp warps each sample's RT axis onto a reference
  xdata <- adjustRtime(xdata, param = ObiwarpParam(binSize = 0.6))

  # 3. Correspondence: group peaks across samples into shared "features" by density in m/z-RT space
  pdp <- PeakDensityParam(sampleGroups = sample_group, bw = 5, minFraction = 0.5)
  xdata <- groupChromPeaks(xdata, param = pdp)

  # 4. Gap filling: integrate raw signal at feature coordinates where no peak was detected
  fillChromPeaks(xdata, param = ChromPeakAreaParam())
}

# xdata <- run_xcms_pipeline(list.files("data/mzML", "\\.mzML$", full.names = TRUE),
#                            rep(c("control", "treatment"), each = 6))
# feature_table <- featureValues(xdata, value = "into", method = "maxint")
# write.csv(cbind(featureDefinitions(xdata)[, c("mzmed", "rtmed")], feature_table), "feature_table.csv")
```

**Goal:** Correct systematic intensity drift and drop unreliable features before statistics.
**Approach:** Probabilistic Quotient Normalization (PQN) removes sample-wise dilution effects; CV filtering on pooled QC injections removes features with unstable measurement.

```python
import pandas as pd


def pqn_normalize(feature_table: pd.DataFrame) -> pd.DataFrame:
    """Probabilistic Quotient Normalization for an LC-MS feature table.

    feature_table: samples (rows) x features (columns) of raw intensities.
    Each sample is first total-intensity-normalized, then scaled by the median
    ratio ("quotient") of its features to a reference spectrum (the per-feature
    median across samples) -- corrects dilution without being skewed by a few
    high-intensity outlier features.
    """
    row_sums = feature_table.sum(axis=1)
    integral_norm = feature_table.div(row_sums, axis=0)
    reference_spectrum = integral_norm.median(axis=0)
    quotients = integral_norm.div(reference_spectrum, axis=1)
    median_quotient = quotients.median(axis=1)
    return feature_table.div(median_quotient, axis=0)


def filter_by_qc_cv(feature_table: pd.DataFrame, qc_samples: list, cv_threshold: float = 0.30) -> pd.DataFrame:
    """Drop features whose coefficient of variation in pooled QC samples exceeds cv_threshold.

    Standard metabolomics QC: CV > 30% in repeated QC injections marks a feature
    as too noisy/unstable to trust for biological comparisons.
    """
    qc_data = feature_table.loc[qc_samples]
    cv = qc_data.std(axis=0) / qc_data.mean(axis=0)
    keep_features = cv[cv <= cv_threshold].index
    return feature_table[keep_features]
```

**Goal:** Collapse redundant adduct/isotope features into one entry per metabolite before annotation.
**Approach:** Features that co-elute (same RT) and differ by a known adduct mass shift almost always come from the same underlying molecule.

```python
import pandas as pd

# Mass shift of each adduct relative to [M+H]+, in Da
ADDUCT_MASS_SHIFT = {"[M+Na]+": 21.9819, "[M+K]+": 37.9559, "[M+NH4]+": 17.0265}


def group_adducts(features: pd.DataFrame, mz_col: str = "mz", rt_col: str = "rt",
                   rt_tolerance: float = 0.05, ppm_tolerance: float = 10) -> pd.DataFrame:
    """Assign an adduct_group id to co-eluting, mass-shift-related features.

    features: DataFrame indexed by feature id with numeric mz_col/rt_col columns.
    rt_tolerance is in the same units as rt_col (minutes for XCMS 'rtmed').
    Returns a copy of features with an added 'adduct_group' integer column.
    """
    features = features.copy()
    features["adduct_group"] = -1
    assigned, group_id = set(), 0
    for i, row in features.iterrows():
        if i in assigned:
            continue
        features.loc[i, "adduct_group"] = group_id
        assigned.add(i)
        for shift in ADDUCT_MASS_SHIFT.values():
            expected_mz = row[mz_col] + shift
            ppm_window = expected_mz * ppm_tolerance / 1e6
            match = features[
                (features[rt_col].sub(row[rt_col]).abs() <= rt_tolerance)
                & (features[mz_col].sub(expected_mz).abs() <= ppm_window)
            ]
            for j in match.index:
                if j not in assigned:
                    features.loc[j, "adduct_group"] = group_id
                    assigned.add(j)
        group_id += 1
    return features
```

## Pitfalls

- **Mass accuracy drift**: calibrate/recalibrate m/z before peak picking (`ppm` in `CentWaveParam` too tight/loose otherwise causes split or merged peaks and false IDs)
- **Adduct confusion**: the same metabolite yields [M+H]+, [M+Na]+, [M+K]+ ions — group adducts (see above) before any statistical test, or you will double-count one metabolite as several
- **Missing value handling**: zeros/NA in a feature table usually mean below-detection-limit, not truly absent — run `fillChromPeaks` first, then impute remaining gaps with kNN or min/2, never zero-fill
- **CentWave parameters are instrument-specific**: `ppm`, `peakwidth`, and `snthresh` tuned for Orbitrap data will misfire on QTOF data — check EICs of known standards before running the full batch
- **Skipping QC injections**: without pooled QC samples interspersed through the run you cannot distinguish batch drift from biology; always include them and use their CV/LOESS trend for correction

## See Also

- `bio-applied-metabolite-identification` — spectral matching, formula assignment, MSEA once you have a feature table
- `bio-applied-proteomics` — shared LC-MS instrumentation and quantification concepts
- `bio-applied-statistics-for-bioinformatics` — t-test/limma/FDR for differential abundance
- `bio-applied-dimensionality-reduction` — PCA/OPLS-DA on the normalized feature table