Exploring QC: RunScatter and PCA

Tutorial

This tutorial covers two complementary QC-exploration views on the built-in lipidomics dataset. plot_runscatter() visualises a feature’s signal across the analysis sequence — helping to identify trends, detect outliers, and assess analytical performance — while Principal Component Analysis (plot_pca(), plot_pca_loading(), detect_outlier_pca()) provides a multivariate view for spotting injection outliers and residual batch effects. The RunScatter examples use a single feature (TG 48:2 [-18:1]).

Time ~20 min  ·  Level Intermediate  ·  Prerequisites Basic workflow

Setup

library(mrmhub)

mexp <- lipidomics_dataset
mexp <- normalize_by_istd(mexp)
mexp <- calc_qc_metrics(mexp)

Basic plot

A minimal call plots all QC and sample types for one feature.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                rows_page = 1, cols_page = 1)

Selecting QC types

Use qc_types to display only specific sample types.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                qc_types = c("BQC", "SPL"),
                rows_page = 1, cols_page = 1)

Batch display

Set show_batches = TRUE and batch_zebra_stripe = TRUE to highlight batch boundaries with alternating shaded areas.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                show_batches = TRUE,
                batch_zebra_stripe = TRUE,
                rows_page = 1, cols_page = 1)

Reference lines

Show mean ± k × SD reference lines with show_reference_lines = TRUE. Use reference_sd_shade to display the range as a shaded band instead of lines.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                qc_types = c("BQC", "SPL"),
                show_reference_lines = TRUE,
                ref_qc_types = "BQC",
                reference_k_sd = 2,
                reference_sd_shade = TRUE,
                rows_page = 1, cols_page = 1)

Outlier capping

Use cap_outliers = TRUE to cap extreme values based on MAD fences. This is useful when outliers obscure the trends of QC or study samples.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                cap_outliers = TRUE,
                cap_sample_k_mad = 3,
                rows_page = 1, cols_page = 1)

Log scale

Set log_scale = TRUE to apply a log10 transformation to the y-axis. Zero or negative values are replaced with the minimum positive value divided by 5.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                log_scale = TRUE,
                rows_page = 1, cols_page = 1)

Removing gaps in the analysis sequence

Gaps in the x-axis can arise in two ways:

1. Filtered QC types — when only a subset of QC types is selected, the unselected positions leave gaps. Use collapse_excluded = TRUE to close them.
2. Unannotated analyses — when some runs in the analytical sequence have no corresponding entry in @dataset (e.g. solvent blanks or system suitability injections that were not imported), the analysis order contains discontinuities. Use remove_gaps = TRUE to collapse these and mark the former gap positions with vertical indicator lines.

Both options can be combined.

Collapsing excluded QC types

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                qc_types = c("LTR"),
                collapse_excluded = TRUE,
                rows_page = 1, cols_page = 1)

Removing gaps from unannotated analyses

To demonstrate this feature, we first simulate a dataset with unannotated analyses by excluding a block of runs from batch 2. This creates a contiguous gap in the analysis order, as would occur when solvent blanks or system suitability injections are not imported.

mexp_gaps <- exclude_analyses(
  mexp,
  analyses = c(
    paste0("Longit_batch2_", 1:45),
    "Longit_batch2_PQC 12",
    "Longit_batch2_PQC 13",
    "Longit_batch2_PQC 14",
    "Longit_batch2_PQC 15",
    "Longit_batch2_TQC13",
    "Longit_batch2_TQC14",
    "Longit_batch2_TQC15"
  ),
  clear_existing = TRUE
)
mexp_gaps <- normalize_by_istd(mexp_gaps)
mexp_gaps <- calc_qc_metrics(mexp_gaps)

plot_runscatter(mexp_gaps, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                remove_gaps = TRUE,
                rows_page = 1, cols_page = 1)

Combining both

When filtering to a single QC type and the sequence contains unannotated runs, combine collapse_excluded and remove_gaps to produce a compact, fully annotated plot.

plot_runscatter(mexp_gaps, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                qc_types = c("BQC"),
                collapse_excluded = TRUE,
                remove_gaps = TRUE,
                rows_page = 1, cols_page = 1)

Label wrapping

When feature names are long, strip labels can overflow. Set label_wrap = TRUE to wrap labels across multiple lines, controlled by label_wrap_width. Here we show three features to illustrate the effect on strip text.

plot_runscatter(mexp, variable = "intensity",
                include_feature_filter = "TG 48\\:2 \\[\\-18:1\\]",
                rows_page = 1, cols_page = 3,
                specific_page = 1,
                label_wrap = TRUE,
                label_wrap_width = 8)

Trend curves

Trend curves visualise the fitted signal used during drift or batch correction. They require a prior correction step — show_trend = TRUE will error if no drift or batch correction has been applied. To plot the values before the last correction, append _before to the variable name (e.g. variable = "intensity_before"); for the original uncorrected values use _raw (e.g. variable = "intensity_raw").

See the Drift Correction tutorial for a worked example.

Principal component analysis (PCA)

PCA is a routine multivariate check in targeted MS workflows: it helps spot injection outliers, visualise residual batch effects, and confirm that drift or batch corrections reduced unwanted variance. The plots below reuse the mexp prepared in Setup.

Score plot by QC type

The score plot summarises sample variance in two dimensions. Biological QCs (BQC) should cluster tightly near the centre if normalisation and any corrections succeeded, while study samples (SPL) typically show wider biological spread. The ellipse_variable controls the grouping of the confidence ellipses (qc_type, batch_id, or "none"); a markedly dispersed QC cluster suggests insufficient normalisation or remaining instrument drift.

plot_pca(mexp,
         variable = "norm_intensity",
         qc_types = c("BQC", "SPL"),
         ellipse_variable = "qc_type")

Score plot by batch

If samples separate along PC1 or PC2 by batch_id, batch effects persist and a centering-based correction is likely warranted. Running the same plot before and after correct_batch_centering() confirms that the batch ellipses overlap after correction.

plot_pca(mexp,
         variable = "norm_intensity",
         qc_types = c("BQC", "SPL"),
         ellipse_variable = "batch_id")

Loadings

The loading plot identifies the features driving the principal components; features at the extremes contribute most to sample separation. A single feature dominating PC1 should be inspected — it is often a saturated or contaminated transition rather than a biological signal.

plot_pca_loading(mexp,
                 variable = "norm_intensity",
                 qc_types = c("BQC", "SPL"))

Detecting outlier injections

detect_outlier_pca() flags analyses whose score on a chosen principal component lies outside a fence, defined as either mean ± k·SD (outlier_detection = "sd") or median ± k·MAD (outlier_detection = "mad"), with fence_multiplicator setting k. It returns the analysis_id values that exceed the fence on the selected component, or NULL if none are flagged.

outlier_ids <- detect_outlier_pca(
  mexp,
  variable = "intensity",
  filter_data = FALSE,
  pca_component = 1,
  qc_types = c("BQC", "SPL"),
  outlier_detection = "mad",
  fence_multiplicator = 4
)

outlier_ids
#> [1] "Longit_batch6_51"

Choice of fence method

The SD-based fence is sensitive to the very outliers it tries to detect: a single extreme observation inflates the SD estimate. The MAD-based fence is robust to a few large deviations and is generally preferable for QC screening. A typical multiplier is fence_multiplicator = 3 (≈99.7% under normality) or 4 for a more permissive screen. The function evaluates one principal component at a time — re-run with pca_component = 2 to screen PC2.

A sample flagged by PCA is a candidate for investigation, not a verdict. Outlier patterns frequently reflect genuine biology (e.g. a disease group or a sex difference), so an injection should be excluded only where a documented technical cause is identified — a failed injection, contamination, instrument fault, or sample-handling error.

Excluding analyses and features

After visual and documented confirmation, remove the offending injections or features. exclude_analyses() / exclude_features() set the affected rows aside for downstream steps while the original data remain in mexp@dataset_orig; setting clear_existing = TRUE replaces any previous exclusion list rather than appending.

# Review candidates before excluding
mexp@annot_analyses |>
  dplyr::filter(analysis_id %in% outlier_ids) |>
  dplyr::select(analysis_id, qc_type, batch_id, analysis_order)

# Exclude after technical confirmation
mexp <- exclude_analyses(mexp, analyses = outlier_ids, clear_existing = FALSE)

# Exclude a known-problematic feature (e.g. a saturated transition)
mexp <- exclude_features(mexp, features = c("PC 32:0"), clear_existing = FALSE)

Interpretation guide

PCA pattern	Likely cause	Action
BQC samples dispersed across the score plot	Poor precision; ISTD assignment problem	Inspect plot_normalization_qc(); verify ISTD pairing
Clear separation by batch_id on PC1 or PC2	Uncorrected batch effect	Apply correct_batch_centering()
Single injection isolated from all groups	Technical outlier	Investigate cause; exclude only with documented reason
One feature dominates loadings on PC1	Saturation, contamination, or single-transition artefact	Inspect with plot_runscatter(); exclude feature if confirmed
BQC tight, SPL spread	Genuine biological variability	Proceed

Next steps

• Drift and Batch Correction — apply and inspect corrections
• Visualisation Functions — the full list of plotting functions
• Basic MRMhub Workflow — where these checks fit in the pipeline