How ProteinBatcher differs from msqrob2

The reviewer specifically asked how ProteinBatcher relates to msqrob2. The two packages overlap in purpose (differential abundance in MS-based proteomics) but occupy clearly different design points, and they are complementary rather than competing.

msqrob2 contributes statistical methodology. It implements a robust linear (mixed) model framework whose parameter estimates can be stabilised by ridge regression, empirical-Bayes variance estimation and robust M-estimation, and it offers a hurdle workflow that handles missing values without requiring imputation. It is built on the QFeatures infrastructure, so it can model the data while consistently linking the PSM, peptide and protein levels, and can start either from raw peptide intensities or from summarised protein values. Its value lies in the modelling engine itself.

ProteinBatcher, by contrast, does not introduce new statistical methods. It orchestrates established limma modelling into a fixed, reproducible downstream pipeline that operates on a plain protein-level SummarizedExperiment. Its distinctive contributions are operational rather than statistical: a single-call pipeline with a fixed annotation schema, a specific organisation of results into main / common / interaction effects, and dedicated interpretation plots (most notably the deregulogram). In short, msqrob2 is primarily a method (and a peptide-to-protein modelling framework), whereas ProteinBatcher is primarily a standardised wrapper and visualisation layer around protein-level limma analysis.

We therefore do not claim that flexible designs, interaction testing or integrated workflows are unavailable elsewhere in Bioconductor: msqrob2, among others, supports them. The motivation for ProteinBatcher is to provide a convention-driven, low-boilerplate route for teams that repeatedly run protein-level limma analyses with the same input/output conventions and that want the specific effect-separation and deregulogram outputs described below.

Input data

1. Quantitative matrix. A protein-by-sample matrix of quantitative measurements, as produced by DIA-NN (pg_matrix). Rows are protein groups and columns are samples.
2. Sample annotation table. A tab-delimited table describing the experimental design and linking each sample to its biological and technical metadata.

Both are shipped with the package as example files. The chunk below locates them with system.file() and reads the annotation table with base R (read.delim()); no external reader package is required. We also record whether the example files are available, so that the rest of the vignette degrades gracefully if they are not.

# Load package
library(ProteinBatcher)

# Input files shipped with the package
path_annotation <- system.file("extdata", "annotation_HaCaT.tsv",
                                package = "ProteinBatcher")
path_pgmatrix   <- system.file(
  "extdata", "2024MK017_HaCaT_Stimulation_1to24_Astral_report.pg_matrix.tsv",
  package = "ProteinBatcher")
path_output <- tempdir()

# Inspect annotation structure (base R, sin readr)
ann <- utils::read.delim(path_annotation, check.names = FALSE)

The annotation table read above has one row per sample. We can confirm its shape and look at the first rows:

dim(ann)
#> [1] 24  7
head(ann)
#>                                                                          file
#> 1 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_001_01_4pto.raw
#> 2 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_002_01_4pto.raw
#> 3 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_003_01_4pto.raw
#> 4 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_004_01_4pto.raw
#> 5 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_005_01_4pto.raw
#> 6 Y:\\data\\orbitrap_astral1\\Raw\\2502\\ELMO\\2024MK017_ELMO_006_01_4pto.raw
#>        sample sample_name condition replicate batch donor_id
#> 1 NoTreated_1 NoTreated_1 NoTreated         1  day1        1
#> 2     IL13 _1     IL13 _1      IL13         1  day1        1
#> 3      IL22_1      IL22_1      IL22         1  day1        1
#> 4    COMBO _1    COMBO _1     COMBO         1  day1        1
#> 5 NoTreated_2 NoTreated_2 NoTreated         1  day1        1
#> 6     IL13 _2     IL13 _2      IL13         1  day1        1

Sample annotation structure

For the high-level run_proteomics_pipeline() entry point, the annotation table is validated against a fixed set of column names. These names are required (not merely suggested): the input validator stops with an informative error if any of them is missing, because the pipeline locates each piece of metadata by name when it builds colData and when it assembles the limma design. The required columns are:

• file — raw-data file identifier for each sample (used as the matching label against the quantitative matrix columns).
• sample and sample_name — sample identifiers; sample_name must be unique and becomes the sample id in colData.
• condition — the main experimental factor of interest (the factor that the tests contrasts compare).
• replicate — repeated-measurement identifier, used by paired / correlation-aware designs.
• batch — a batch or sample-preparation factor (e.g. acquisition day), used as a covariate and as the interaction factor in the examples below.
• donor_id — biological replicate / donor identifier, used as the blocking factor in paired or repeated-measures models.

The set is intentionally fixed so that the same pipeline code can be reused across experiments without per-experiment column-mapping: the conventions are the contract. Which of these factors actually enters the model is then determined by the design formula and by the paired / block_effect parameters; columns that are present but not referenced by the formula simply remain available in colData for inspection.

Workflow parameters

In addition to the two inputs, the pipeline takes parameters that control filtering, modelling and output generation. The most important are:

• percent_missing — missingness threshold (in %) for feature filtering.
• tests — main condition contrasts to evaluate, written as "<level>_vs_<reference>".
• tests_interaction — interaction coefficients to evaluate, or "NA" to disable interaction testing.
• formula — the limma design formula (main effects and interactions).
• reference_condition — the reference level of condition.
• paired / block_effect — options for paired and correlation-aware models.
• ldv_source — strategy for left-censored (LDV) imputation ("global" or "per-condition").
• plots — whether to also write result tables and PDF figures to path_output.

Accessing the organised effects

effects is a list keyed by the main contrast, so a given contrast can be retrieved directly with $. Each element is itself a list of three SummarizedExperiment objects — all_common_effect, common_effect and interaction_effect — and the per-contrast statistics (log2FC, CI.L, CI.R, p.val, p.adj) are stored in their rowData() under the "<stat>_<contrast>" naming convention (for example log2FC_IL13_vs_NoTreated). This is the structure consumed by the downstream visualisation helpers.

# Contrasts available, and direct access with `$`
names(res$effects)
#> [1] "IL13_vs_NoTreated"  "IL22_vs_NoTreated"  "COMBO_vs_NoTreated"
il13 <- res$effects$IL13_vs_NoTreated
names(il13)
#> [1] "common_effect"      "interaction_effect" "all_common_effect"

# Each slot is a SummarizedExperiment; statistics live in rowData()
il13$all_common_effect
#> class: SummarizedExperiment 
#> dim: 9600 24 
#> metadata(5): log2transform exp lfq_type level imputation_map
#> assays(1): intensity
#> rownames(9600): A0A024RBG1 A0A096LP01 ... Q9Y6Y0 Q9Y6Y8
#> rowData names(11): Protein.Group Protein.Names ...
#>   CI.R_IL13_vs_NoTreated log2FC_IL13_vs_NoTreated
#> colnames(24): NoTreated_1 IL13 _1 ... IL22_6 COMBO _6
#> colData names(8): file sample ... donor_id label
head(colnames(SummarizedExperiment::rowData(il13$all_common_effect)))
#> [1] "Protein.Group" "Protein.Names" "Genes"         "name"         
#> [5] "Index"         "imputed"

Filtering and imputation summary tables

When plots = TRUE, the pipeline also writes plain-text (CSV) summaries of the preprocessing steps so that they are fully transparent and can be re-read with base R. The package ships precomputed copies of these tables as example outputs. The chunk below reads them with read.csv() if they are present:

• HaCaT_filtered_out_proteins.csv lists the proteins removed by the missingness filter, so users can verify exactly which features were excluded.
• HaCaT_Imputation_before.csv and HaCaT_Imputation_after.csv give the protein-by-sample matrix immediately before and after imputation, allowing a direct before/after comparison for quality control.

path_filtered <- system.file(
  "extdata", "HaCaT_filtered_out_proteins.csv",
  package = "ProteinBatcher"
)
path_before <- system.file(
  "extdata", "HaCaT_Imputation_before.csv",
  package = "ProteinBatcher"
)
path_after <- system.file(
  "extdata", "HaCaT_Imputation_after.csv",
  package = "ProteinBatcher"
)

if (nzchar(path_filtered)) {
  filtered_df <- read.csv(path_filtered, check.names = FALSE,
                          stringsAsFactors = FALSE)
  head(filtered_df)
}

if (nzchar(path_before) && nzchar(path_after)) {
  before_df <- read.csv(path_before, check.names = FALSE,
                        stringsAsFactors = FALSE)
  after_df  <- read.csv(path_after,  check.names = FALSE,
                        stringsAsFactors = FALSE)
  # First proteins, first few columns, before vs after imputation
  head(before_df[, seq_len(min(5, ncol(before_df)))])
  head(after_df[, seq_len(min(5, ncol(after_df)))])
}

The “before” table still contains missing values (NA); in the “after” table these have been replaced either by a within-condition mean (when a protein is only partially missing in a condition) or by a left-censored low-value draw (when it is essentially absent), which is why the imputed matrix has no missing entries.

Differential effects: IL13 vs NoTreated

From here we focus on the IL13_vs_NoTreated contrast to illustrate how the three effect types are read. The figures below are the PNG outputs shipped with the package; each chunk renders the figure only if it is present.