stanassay

Bayesian 5-Parameter Logistic Curve Fitting for Luminex and ELISA Assays

stanassay is an R package that fits hierarchical Bayesian 5PL (5-Parameter Logistic) standard curves to multi-plate Luminex (xMAP) and ELISA immunoassay data. It wraps pre-compiled Stan models inside an ergonomic R6 class, providing a single, consistent API for fitting, back-calculation, LOQ estimation, and visualisation.

---

Installation

The package is not yet on CRAN. Install directly from the repository:

# Install dependencies first
install.packages(c("rstan", "R6", "dplyr", "plotly", "ggplot2"))

# Install stanassay (replace path with your local clone)
devtools::install("path/to/stanassay")

Note: rstan requires a working C++ toolchain. See the RStan Getting Started guide if you have not used Stan before.

---

Why stanassay?

The Exact Model (error_model = "exact")

The standard hierarchical 5PL model treats the nominal standard concentrations as exact known values. It uses a heteroscedastic Student-T likelihood and pools information across plates through a shared hyperprior, allowing robust parameter estimation even on plates with sparse data.

Use this when:

• Concentrations come from a certified reference standard
• Dilutions were made robotically with negligible pipetting error
• Note: Can be used for both ELISA and Luminex (xMAP)

The Error-in-Variables (EiV) Model (error_model = "eiv")

Serial dilution introduces compounding pipetting uncertainty: each transfer step adds ~2 % CV, so the lowest standards have the most uncertain true concentration. The EiV model treats nominal concentrations as measured with noise and infers the latent true concentrations jointly with the 5PL parameters.

Use this when:

• Standards are made by manual serial dilution
• You want uncertainty in standard concentrations to propagate into back-calculated sample concentrations
• Assay range spans ≥ 3 log₁₀ decades (compounding error is non-trivial)
• Note: Can be used for both Luminex (xMAP) and ELISA

---

Quick Start

library(stanassay)

# ── 1. Prepare standard curve data ─────────────────────────────────────────
xmap_std <- data.frame(
  conc     = c(50000, 10000, 3333, 1111, 370, 123, 41, 14, 5, 2,
               50000, 10000, 3333, 1111, 370, 123, 41, 14, 5, 2),
  mfi      = c(19129, 17880,   49,   22, 289,  13,  47, 41, 144, 43,
                1126,  1107,  214,  171,  50,  51,  42, 44,  39,  9),
  plate_id = c(rep("Plate_1", 10), rep("Plate_2", 10))
)

# ── 2. Initialise the R6 object ─────────────────────────────────────────────
assay <- StanAssay$new(
  std_data          = xmap_std,
  concentration_col = "conc",
  response_col      = "mfi",
  plate_col         = "plate_id",
  assay_type        = "xmap",     # "xmap" or "elisa"
  error_model       = "eiv"       # "eiv" (pipetting error) or "exact" (perfect stds)
)

# ── 3. Fit the Bayesian model ───────────────────────────────────────────────
assay$fit(
  n_iter   = 2000,   # total iterations per chain (including warmup)
  n_chains = 4,      # parallel chains
  cores    = 4       # auto-caps at (available physical cores - 1) for safety
)

# ── 4. Extract raw curve geometry for custom ggplot2 figures ───────────────
pd <- assay$extract_plot_data("Plate_1")

# pd$curve_df    — 100-row data frame: concentration, mfi_median, lower/upper 95% CI
# pd$asymptotes  — list: a_mean (lower), d_mean (upper)
# pd$raw_data    — observed standard points for this plate

library(ggplot2)
ggplot(pd$curve_df, aes(x = concentration)) +
  geom_ribbon(aes(ymin = mfi_lower_95, ymax = mfi_upper_95),
              fill = "steelblue", alpha = 0.25) +
  geom_line(aes(y = mfi_median), colour = "steelblue", linewidth = 1) +
  geom_point(data = pd$raw_data, aes(y = mfi), size = 2) +
  scale_x_log10() + scale_y_log10() +
  labs(title = "Plate_1 — Bayesian 5PL Fit",
       x = "Concentration (log scale)", y = "MFI (log scale)") +
  theme_bw()

# ── 5. Back-calculate sample concentrations ─────────────────────────────────
samples <- data.frame(
  mfi      = c(500, 2000, 8000, 15000),
  plate_id = "Plate_1"
)

results <- assay$predict_samples(samples)
print(results[, c("mfi", "predicted_conc_mean",
                  "predicted_conc_lower", "predicted_conc_upper",
                  "gate_class")])

# ── 6. Interactive plotly curve ─────────────────────────────────────────────
p <- assay$plot("Plate_1", sample_data = results)
print(p)   # renders in RStudio Viewer / browser

# ── 7. CDAN precision profile & LOQ bounds ─────────────────────────────────
cdan <- assay$compute_cdan("Plate_1", n_grid = 200)

cat("LLOQ (15% CV):", cdan$lloq_15, "\n")
cat("ULOQ (15% CV):", cdan$uloq_15, "\n")

---

Package Structure

File	Contents
R/StanAssay.R	R6 class: $new(), $fit(), $predict_samples(), $plot(), $compute_cdan(), $extract_plot_data()
R/predict.R	backcalc_samples_bayesian() — posterior inversion engine
R/cdan.R	compute_cdan_precision_profile(), compute_precision_profile_loq()
R/plot.R	build_plot_data(), plot_bayesian_plate()
R/utils_priors.R	compute_dynamic_priors()
inst/stan/hierarchical_5pl.stan	Standard 5PL Stan model
inst/stan/hierarchical_5pl_eiv.stan	Error-in-Variables 5PL Stan model

---

License

MIT © 2026