R/tabdiag.R

In R4VN: Health Data Analysis and Publication-Ready Reporting

#' Comprehensive Diagnostic Accuracy and ROC Analysis
#'
#' @description
#' `tabdiag()` is the comprehensive diagnostic-accuracy command in R4VN. It is
#' designed for binary diagnostic tests, continuous or ordinal biomarkers,
#' prediction scores, and comparisons of paired ROC curves. The first variable
#' is the binary reference standard (gold standard); diagnostic tests can be
#' supplied directly in `...` and/or through `vars()`.
#'
#' The default output follows the R4VN philosophy: a compact publication-ready
#' Viewer table is shown, while the returned object retains the complete set of
#' diagnostic measures and ROC information for later inspection or export.
#' Interpretation is opt-in (`interpretation = FALSE` by default).
#'
#' @param outcome Binary reference-standard outcome. Supply an unquoted
#'   variable name or a single character variable name.
#' @param ... One or more diagnostic tests/markers. Binary variables are
#'   analyzed as 2 x 2 tests. Numeric variables, dates/times, and ordered
#'   factors are analyzed by ROC.
#' @param data Optional data frame. If omitted, the active R4VN data selected
#'   by `usedf()` or `opendata(..., active = TRUE)` is used.
#' @param vars Optional diagnostic-marker selection created by `vars()`. This
#'   may be used instead of, or together with, markers supplied in `...`.
#' @param event Positive level of the reference-standard outcome. If omitted,
#'   R4VN recognizes common positive encodings such as 1, TRUE, Yes, Positive,
#'   Case, or Co; otherwise it uses the second factor/observed level and reports
#'   the selected positive outcome explicitly.
#' @param positive Positive level for binary diagnostic tests. It may be a
#'   single value used for all binary tests, an unnamed vector in test order,
#'   or a named vector such as
#'   `positive = c(rapid = "Positive", ct = "Abnormal")`.
#' @param direction ROC direction: `"auto"` (default), `"<"`, or `">"`.
#'   With R4VN conventions, `"<"` means higher marker values indicate the
#'   positive outcome and R4VN reports cutoffs as `>=`; `">"` means lower
#'   values indicate the positive outcome and R4VN reports cutoffs as `<=`.
#' @param best Optimal-cutoff methods. Default `"youden"`. May contain one or
#'   more of `"youden"`, `"closest"`, `"ruleout"`, and `"rulein"`.
#'   `"ruleout"` finds the threshold with the highest specificity while
#'   sensitivity is at least `target`; `"rulein"` finds the threshold with the
#'   highest sensitivity while specificity is at least `target`. Set
#'   `best = NULL` to suppress automatically selected cutoffs.
#' @param cuts Optional numeric cutoff(s) requested by the user. Every cutoff
#'   becomes a separate diagnostic-performance row for each quantitative
#'   marker. This is useful for clinically predefined thresholds.
#' @param target Target sensitivity/specificity for `"ruleout"`/`"rulein"`,
#'   default 0.95.
#' @param ci `TRUE` (default), `FALSE`, or a character vector of measures for
#'   which confidence intervals are required. Examples include
#'   `c("sens", "spec")`, `c("sens","spec","ppv","npv")`, and
#'   `c("auc","lr+","lr-","dor")`. `TRUE` requests every currently
#'   implemented stable interval.
#' @param ci_level Confidence level, default 0.95.
#' @param ci_method Binomial confidence-interval method for 2 x 2 proportions:
#'   `"auto"`/`"wilson"` or `"exact"`. AUC uses an internally implemented DeLong CI for an ordinary
#'   empirical ROC and stratified bootstrap CI for partial AUC.
#' @param cut_ci Logical; calculate a bootstrap confidence interval for an
#'   automatically selected Youden/closest cutoff. Default FALSE because this
#'   can be computationally expensive.
#' @param boot Number of bootstrap replicates for cutoff CI, partial-AUC CI,
#'   and bootstrap ROC comparison. Default 2000.
#' @param prevalence Optional target disease prevalence, strictly between 0 and
#'   1. When supplied, PPV and NPV are recalculated for that prevalence.
#' @param roc Logical; calculate ROC analysis for eligible non-binary tests.
#'   Default TRUE.
#' @param partial_auc Optional numeric vector of length two defining a partial
#'   AUC range, for example `c(1, 0.80)` for specificity 100%-80% when
#'   `partial_focus = "specificity"`.
#' @param partial_focus `"specificity"` (default) or `"sensitivity"`.
#' @param partial_correct Logical; request corrected partial AUC where
#'   supported by the R4VN internal ROC engine.
#' @param compare Logical; for two or more ROC-eligible markers, calculate
#'   pairwise comparisons of AUC. Default FALSE.
#' @param compare_method ROC comparison method: `"delong"` (default) or
#'   `"bootstrap"`. Both are implemented internally by R4VN.
#' @param adjust Multiplicity adjustment for pairwise comparison p-values,
#'   passed to `stats::p.adjust()`. Examples: `"none"`, `"holm"`, `"BH"`.
#' @param all_cuts Logical; if TRUE, calculate full 2 x 2 diagnostic properties
#'   for every finite empirical ROC threshold and store them in
#'   `$all_cutoffs`. FALSE by default because this table can be very large.
#' @param missing Logical; include the marker-specific analysis sample table in
#'   printed/Viewer output. The sample table is always retained in
#'   `$descriptive`. Default FALSE.
#' @param show Logical; open the formatted result in the Viewer. Default TRUE.
#'   For backward compatibility, a character vector supplied to `show` is
#'   interpreted as the old diagnostic-column selection syntax.
#' @param measures Diagnostic measures displayed in the selected-threshold
#'   table. Default `"core"` gives N, sensitivity, specificity, PPV, NPV,
#'   LR+, LR-, and accuracy. Other convenient profiles are `"minimal"`,
#'   `"publication"`, `"counts"`, and `"all"`; or supply an explicit vector
#'   such as `c("sens","spec","ppv","npv","lr+","lr-","dor")`.
#'   Cell counts are presented as a conventional 2 x 2 diagnostic
#'   classification table rather than as separate TP/FP/FN/TN measures.
#'   The returned object always retains all calculated measures.
#' @param hide Optional diagnostic columns to remove from the displayed cutoff
#'   table, for example `hide = c("n","accuracy")`.
#' @param plot `FALSE` (default), `TRUE`, `"roc"`, `"cutoff"`, or `"both"`.
#'   `TRUE` is equivalent to `"roc"`. Requested figures are saved internally,
#'   embedded in the Viewer report, and also sent to the RStudio Plot pane when
#'   running interactively.
#' @param plot_args Named list of graphical arguments passed to
#'   `plot.r4vn_diag()`. Useful options include `color`, `lty`, `line_width`,
#'   `legend`, `legend_position`, `auc`, `diagonal`, `grid`, `xlim`, `ylim`,
#'   `xlab`, `ylab`, `main`, and cutoff-plot controls. ROC axes are constrained
#'   to 0-1 by default.
#' @param interpretation Logical; add a cautious deterministic interpretation
#'   table. Default FALSE. Interpretation describes discrimination and selected
#'   thresholds but does not claim clinical utility or causality.
#' @param digit Number of digits displayed for estimates. Default 2.
#' @param p_digit Number of digits displayed for p-values. Default 3.
#' @param zero_correction Continuity correction used only for a requested DOR
#'   confidence interval when a zero cell is present. Default 0.5.
#' @param title Optional title printed above the output.
#' @param console Logical; also print the traditional Console result. Default
#'   FALSE.
#' @param print Deprecated compatibility argument. When supplied, it overrides
#'   `console`.
#' @param export Optional export format accepted by `tabexport()`, such as
#'   `"html"`, `"docx"`, `"xlsx"`, `"pdf"`, or `"png"`. Multiple formats
#'   may be requested in one call.
#' @param file Optional export filename. Its extension may also determine the
#'   export format.
#' @param open Logical; open the exported file when supported.
#'
#' @details
#' ## Binary diagnostic tests
#'
#' A test with exactly two observed non-missing values is analyzed directly as
#' a 2 x 2 table. `tabdiag()` and `tabdiagi()` use the same internal engine, so
#' identical TP/FP/FN/TN counts give identical sensitivity, specificity,
#' predictive values, likelihood ratios, DOR, accuracy, balanced accuracy,
#' Youden index, F1 score, MCC, kappa, FPR, FNR, FDR, FOR, prevalence,
#' detection rate, and detection prevalence.
#'
#' ## Quantitative and ordered diagnostic markers
#'
#' ROC analysis is implemented inside R4VN using base R. `tabdiag()` therefore
#' does not require `pROC`, `ggplot2`, or another ROC package for routine ROC,
#' AUC, cutoff selection, DeLong inference, partial AUC, or paired AUC comparison.
#' This keeps installation light while preserving a complete analysis object.
#'
#' Numeric markers, Date/POSIX values, and ordered factors with more than two
#' levels are analyzed by ROC. Unordered factors with more than two levels are
#' rejected because a clinically meaningful ordering cannot safely be inferred.
#' A marker that becomes constant after removal of missing values is retained
#' in the sample description and omitted from ROC analysis with an explanatory
#' note instead of crashing the whole report.
#'
#' ## Cutoffs
#'
#' `best = "youden"` maximizes sensitivity + specificity. `"closest"`
#' minimizes distance to the upper-left ROC corner. `"ruleout"` and
#' `"rulein"` choose thresholds meeting a target sensitivity or specificity.
#' Tied optimal thresholds are deliberately retained rather than silently choosing one.
#' User-defined `cuts` are added as separate rows rather than replacing the
#' automatic cutoff.
#'
#' ## Confidence intervals
#'
#' Sensitivity, specificity, PPV, NPV, and accuracy use Wilson or exact
#' binomial intervals. LR+/LR- and DOR use conventional log-scale intervals.
#' Ordinary AUC uses DeLong CI. Partial-AUC CI and optimal-cutoff CI use
#' stratified bootstrap. If prevalence-adjusted PPV/NPV are requested, simple
#' binomial CIs are not attached to those adjusted predictive values because
#' they would not correspond to the target-prevalence estimates.
#'
#' ## AUC inference and paired ROC comparison
#'
#' For each ordinary empirical ROC, R4VN stores DeLong AUC variance, standard
#' error, z statistic, and a two-sided large-sample p-value for H0: AUC = 0.5.
#' When `compare = TRUE`, every pair of eligible markers is compared on its own
#' pairwise complete-case sample. This keeps the comparison paired and makes the
#' reported N explicit when marker missingness differs.
#'
#' ## Missing values
#'
#' Each marker is analyzed on its own complete cases with the outcome. Pairwise
#' ROC comparison uses complete cases for both markers and the outcome. Therefore
#' N may differ between marker-specific AUCs and pairwise comparisons. Use
#' `missing = TRUE` to show the sample accounting table; it is always available
#' as `$descriptive`.
#'
#' ## Returned result contract
#'
#' Backward-compatible components (`$summary`, `$thresholds`, `$coordinates`,
#' `$comparison`, `$all_cutoffs`, `$roc`) are retained. In addition, the object
#' follows the common R4VN reporting contract with `$descriptive`, `$estimates`,
#' `$tests`, `$diagnostics`, `$interpretation`, `$tables`, `$plots`, `$models`,
#' `$metadata`, and `$call`.
#'
#' ## Publication-oriented display
#'
#' The default `measures = "core"` keeps the Viewer readable. Use
#' `measures = "publication"` to add DOR, `measures = "all"` for the full
#' diagnostic panel, or `measures = "counts"` for the 2 x 2 diagnostic
#' classification table only. Cell counts are never shown as a vertical
#' TP/FP/FN/TN measure list. Calculations are never discarded by this display
#' choice. When a full performance table contains only a few rows but many
#' columns, the Viewer automatically presents measures vertically for easier
#' reading.
#'
#' ## ROC plotting
#'
#' ROC plots are drawn from false-positive rate (1-specificity) and sensitivity
#' coordinates using ordinary base graphics. The default axes are exactly 0 to
#' 1 (`xaxs = "i"`, `yaxs = "i"`), preventing the small negative/>1 axis
#' extensions that can otherwise appear in automatic plotting systems.
#'
#' @return
#' An object of class `r4vn_diag`. Important components include:
#' 
#' * `$summary`: marker-specific AUC results.
#' * `$thresholds`: full diagnostic metrics at binary tests and selected cutoffs.
#' * `$comparison`: pairwise AUC comparisons.
#' * `$descriptive`: marker-specific N/missing/case/control accounting.
#' * `$estimates`: structured AUC and threshold estimates.
#' * `$tests`: AUC-vs-0.5 and pairwise ROC tests.
#' * `$diagnostics`: ROC coordinates and optional all-cutoff table.
#' * `$tables`: publication-ready display tables used by Viewer/export, including
#'   `$tables$Confusion_matrix` for selected tests/cutoffs.
#' * `$interpretation`: NULL by default, or an interpretation table when requested.
#' * `$roc` / `$models$roc`: lightweight internal R4VN ROC objects.
#'
#' @section Examples - 1. Simplest binary test:
#' ```
#' d_bin <- data.frame(
#'   disease = c(1,1,1,1,0,0,0,0),
#'   rapid   = c(1,1,1,0,1,0,0,0)
#' )
#' tabdiag(disease, rapid, data = d_bin)
#' ```
#'
#' @section Examples - 2. Explicit positive outcome/test level:
#' ```
#' d_txt <- data.frame(
#'   truth = factor(c("No","Yes","Yes","No","Yes","No")),
#'   test  = factor(c("Neg","Pos","Pos","Neg","Neg","Pos"))
#' )
#' tabdiag(truth, test, data = d_txt, event = "Yes", positive = "Pos")
#' ```
#'
#' @section Examples - 3. Several binary tests with named positive levels:
#' ```
#' d_multi_bin <- data.frame(
#'   disease = c(1,1,1,0,0,0,1,0),
#'   rapid = c("Positive","Positive","Negative","Negative","Positive","Negative","Positive","Negative"),
#'   ct = c("Abnormal","Normal","Abnormal","Normal","Normal","Normal","Abnormal","Normal")
#' )
#' tabdiag(
#'   disease, rapid, ct, data = d_multi_bin,
#'   positive = c(rapid = "Positive", ct = "Abnormal")
#' )
#' ```
#'
#' @section Examples - 4. One continuous biomarker:
#' ```
#'   set.seed(101)
#'   d <- data.frame(
#'     disease = rep(c(0,1), each = 100),
#'     crp = c(rnorm(100, 5, 2), rnorm(100, 10, 3))
#'   )
#'   z <- tabdiag(disease, crp, data = d, show = FALSE)
#'   z$summary
#'   z$thresholds
#' ```
#'
#' @section Examples - 5. Select markers with vars():
#' ```
#'   tabdiag(disease, data = d, vars = vars(crp), show = FALSE)
#' ```
#'
#' @section Examples - 6. Automatic and clinically specified cutoffs:
#' ```
#'   tabdiag(
#'     disease, crp, data = d,
#'     best = "youden", cuts = c(5, 7.5, 10), show = FALSE
#'   )
#' ```
#'
#' @section Examples - 7. Several optimal-cutoff definitions:
#' ```
#'   tabdiag(
#'     disease, crp, data = d,
#'     best = c("youden", "closest", "ruleout", "rulein"),
#'     target = 0.90, show = FALSE
#'   )
#' ```
#'
#' @section Examples - 8. Disable automatic cutoffs and use only clinical cutoffs:
#' ```
#'   tabdiag(disease, crp, data = d, best = NULL, cuts = c(6, 8, 10), show = FALSE)
#' ```
#'
#' @section Examples - 9. Confidence-interval control:
#' ```
#'   tabdiag(disease, crp, data = d, ci = c("auc", "sens", "spec"), show = FALSE)
#'   tabdiag(disease, crp, data = d, ci = TRUE, show = FALSE)
#'   tabdiag(disease, crp, data = d, ci = FALSE, show = FALSE)
#' ```
#'
#' @section Examples - 10. Exact binomial CI for diagnostic proportions:
#' ```
#' tabdiag(disease, rapid, data = d_bin, ci = TRUE, ci_method = "exact", show = FALSE)
#' ```
#'
#' @section Examples - 11. Bootstrap CI for the optimal cutoff:
#' ```
#'   tabdiag(
#'     disease, crp, data = d, cut_ci = TRUE,
#'     boot = 200, show = FALSE
#'   )
#' ```
#'
#' @section Examples - 12. Predictive values at a target prevalence:
#' ```
#'   tabdiag(disease, crp, data = d, prevalence = 0.10, show = FALSE)
#' ```
#'
#' @section Examples - 13. Pairwise DeLong comparison:
#' ```
#'   set.seed(102)
#'   d$pct <- c(rnorm(100, 1, 0.8), rnorm(100, 3, 1.2))
#'   tabdiag(disease, crp, pct, data = d, compare = TRUE, show = FALSE)
#' ```
#'
#' @section Examples - 14. Multiple-comparison adjustment:
#' ```
#'   d$score3 <- d$crp + rnorm(nrow(d), 0, 2)
#'   tabdiag(
#'     disease, crp, pct, score3, data = d,
#'     compare = TRUE, adjust = "holm", show = FALSE
#'   )
#' ```
#'
#' @section Examples - 15. Bootstrap ROC comparison:
#' ```
#'   tabdiag(
#'     disease, crp, pct, data = d,
#'     compare = TRUE, compare_method = "bootstrap", boot = 200,
#'     show = FALSE
#'   )
#' ```
#'
#' @section Examples - 16. Partial AUC:
#' ```
#'   tabdiag(
#'     disease, crp, data = d,
#'     partial_auc = c(1, 0.80), partial_focus = "specificity",
#'     partial_correct = TRUE, show = FALSE
#'   )
#' ```
#'
#' @section Examples - 17. Ordered diagnostic score:
#' ```
#'   d_ord <- data.frame(
#'     disease = c(0,0,0,0,1,1,1,1,1,0),
#'     score = ordered(c("Low","Low","Medium","Low","Medium","High","High","Medium","High","Medium"),
#'                     levels = c("Low","Medium","High"))
#'   )
#'   tabdiag(disease, score, data = d_ord, show = FALSE)
#' ```
#'
#' @section Examples - 18. Lower values indicate disease:
#' ```
#'   d$low_marker <- -d$crp
#'   tabdiag(disease, low_marker, data = d, direction = ">", show = FALSE)
#' ```
#'
#' @section Examples - 19. Marker-specific missing values:
#' ```
#'   d_missing <- d
#'   d_missing$crp[1:12] <- NA
#'   d_missing$pct[21:35] <- NA
#'   zmis <- tabdiag(
#'     disease, crp, pct, data = d_missing,
#'     compare = TRUE, missing = TRUE, show = FALSE
#'   )
#'   zmis$descriptive
#'   zmis$comparison
#' ```
#'
#' @section Examples - 20. Compact, publication, counts, and full displays:
#' ```
#'   tabdiag(disease, crp, data = d, measures = "minimal", show = FALSE)
#'   tabdiag(disease, crp, data = d, measures = "publication", show = FALSE)
#'   tabdiag(disease, crp, data = d, measures = "counts", show = FALSE)
#'   tabdiag(disease, crp, data = d, measures = "all", show = FALSE)
#' ```
#'
#' @section Examples - 21. Explicit measure selection and hiding columns:
#' ```
#'   tabdiag(
#'     disease, crp, data = d,
#'     measures = c("sens","spec","ppv","npv","lr+","lr-","dor"),
#'     hide = "npv", show = FALSE
#'   )
#' ```
#'
#' @section Examples - 22. Full metrics at every empirical ROC threshold:
#' ```
#'   zall <- tabdiag(disease, crp, data = d, all_cuts = TRUE, show = FALSE)
#'   head(zall$all_cutoffs)
#' ```
#'
#' @section Examples - 23. Interpretation is explicit opt-in:
#' ```
#'   zint <- tabdiag(disease, crp, data = d, interpretation = TRUE, show = FALSE)
#'   zint$interpretation
#' ```
#'
#' @section Examples - 24. ROC graph with publication controls:
#' ```
#'   tabdiag(
#'     disease, crp, pct, data = d, plot = "roc",
#'     plot_args = list(
#'       line_width = 2, diagonal = TRUE, grid = TRUE,
#'       legend_position = "bottomright",
#'       xlab = "1 - Specificity", ylab = "Sensitivity",
#'       main = "ROC curves"
#'     )
#'   )
#' ```
#'
#' @section Examples - 25. Sensitivity/specificity across cutoffs:
#' ```
#'   tabdiag(
#'     disease, crp, data = d, plot = "cutoff",
#'     plot_args = list(cutoff_mark = TRUE)
#'   )
#' ```
#'
#' @section Examples - 26. Draw both graph families from a saved object:
#' ```
#'   z <- tabdiag(disease, crp, pct, data = d, plot = FALSE, show = FALSE)
#'   plot(z, what = "roc")
#'   plot(z, what = "cutoff")
#' ```
#'
#' @section Examples - 27. Active R4VN data:
#' ```
#' \donttest{
#'   usedf(d)
#'   tabdiag(disease, crp, show = FALSE)
#' }
#' ```
#'
#' @section Examples - 28. Export publication tables:
#' ```
#' \donttest{
#'   tabdiag(
#'     disease, crp, data = d, show = FALSE,
#'     export = "xlsx", file = tempfile(fileext = ".xlsx")
#'   )
#' }
#' ```
#'
#' @section Examples - 29. Inspect the stable R4VN result contract:
#' ```
#'   z <- tabdiag(disease, crp, data = d, show = FALSE)
#'   names(z$tables)
#'   z$estimates$auc
#'   z$tests$auc_vs_0.5
#'   z$diagnostics$coordinates
#'   z$tables$Confusion_matrix
#'   z$metadata
#' ```
#'
#' @section Examples - 30. Full confusion matrix and full performance panel:
#' ```
#'   zfull <- tabdiag(
#'     disease, crp, data = d,
#'     measures = "all", missing = TRUE, show = FALSE
#'   )
#'   zfull$tables$Diagnostic_performance
#'   zfull$tables$Confusion_matrix
#' ```
#'
#' @section Examples - 31. Save publication-ready ROC and cutoff plots:
#' ```
#' \donttest{
#'   zplot <- tabdiag(disease, crp, data = d, show = FALSE)
#'   plot(
#'     zplot, what = "roc",
#'     file = tempfile(fileext = ".png"),
#'     width = 7, height = 7, res = 300,
#'     line_width = 2, grid = TRUE
#'   )
#'   plot(
#'     zplot, what = "cutoff",
#'     file = tempfile(fileext = ".pdf"),
#'     width = 7, height = 7, cutoff_mark = TRUE
#'   )
#' }
#' ```
#'
#' @examples
#'   set.seed(123)
#'   dd <- data.frame(
#'     disease = rep(c(0, 1), each = 40),
#'     marker = c(rnorm(40, 0, 1), rnorm(40, 1.5, 1))
#'   )
#'   z <- tabdiag(disease, marker, data = dd, show = FALSE)
#'   z$summary
#'
#' @export
tabdiag <- function(
  outcome,
  ...,
  data = NULL,
  vars = NULL,
  event = NULL,
  positive = NULL,
  direction = c("auto", "<", ">"),
  best = "youden",
  cuts = NULL,
  target = 0.95,
  ci = TRUE,
  ci_level = 0.95,
  ci_method = c("auto", "wilson", "exact"),
  cut_ci = FALSE,
  boot = 2000,
  prevalence = NULL,
  roc = TRUE,
  partial_auc = NULL,
  partial_focus = c("specificity", "sensitivity"),
  partial_correct = FALSE,
  compare = FALSE,
  compare_method = c("delong", "bootstrap"),
  adjust = "none",
  all_cuts = FALSE,
  missing = FALSE,
  show = TRUE,
  measures = "core",
  hide = NULL,
  plot = FALSE,
  plot_args = list(),
  interpretation = FALSE,
  digit = 2,
  p_digit = 3,
  zero_correction = 0.5,
  title = NULL,
  console = FALSE,
  print = NULL,
  export = NULL,
  file = NULL,
  open = FALSE
) {
  call <- match.call()
  direction <- match.arg(direction)
  ci_method <- match.arg(ci_method)
  partial_focus <- match.arg(partial_focus)
  compare_method <- match.arg(compare_method)

  viewer_show <- TRUE
  if (is.logical(show) && length(show) == 1L && !is.na(show)) {
    viewer_show <- isTRUE(show)
  } else {
    measures <- show
    viewer_show <- TRUE
  }
  if (!is.null(print)) console <- isTRUE(print)

  if (!is.logical(roc) || length(roc) != 1L || is.na(roc)) {
    stop("`roc` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.logical(compare) || length(compare) != 1L || is.na(compare)) {
    stop("`compare` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.logical(all_cuts) || length(all_cuts) != 1L || is.na(all_cuts)) {
    stop("`all_cuts` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.logical(missing) || length(missing) != 1L || is.na(missing)) {
    stop("`missing` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.logical(interpretation) || length(interpretation) != 1L || is.na(interpretation)) {
    stop("`interpretation` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.logical(open) || length(open) != 1L || is.na(open)) {
    stop("`open` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.list(plot_args) ||
      (length(plot_args) && (is.null(names(plot_args)) || any(!nzchar(names(plot_args)))))) {
    stop("`plot_args` must be a named list.", call. = FALSE)
  }
  if (!is.numeric(ci_level) || length(ci_level) != 1L || is.na(ci_level) ||
      ci_level <= 0 || ci_level >= 1) {
    stop("`ci_level` must be strictly between 0 and 1.", call. = FALSE)
  }
  if (!is.numeric(target) || length(target) != 1L || is.na(target) || target <= 0 || target > 1) {
    stop("`target` must be in (0, 1].", call. = FALSE)
  }
  if (!is.numeric(boot) || length(boot) != 1L || is.na(boot) || boot < 10 || boot != floor(boot)) {
    stop("`boot` must be an integer value of at least 10.", call. = FALSE)
  }
  if (!is.numeric(digit) || length(digit) != 1L || is.na(digit) || digit < 0) {
    stop("`digit` must be a non-negative number.", call. = FALSE)
  }
  if (!is.numeric(p_digit) || length(p_digit) != 1L || is.na(p_digit) || p_digit < 1) {
    stop("`p_digit` must be a positive number.", call. = FALSE)
  }
  if (!is.null(cuts) && (!is.numeric(cuts) || any(!is.finite(cuts)))) {
    stop("`cuts` must be NULL or contain only finite numeric values.", call. = FALSE)
  }
  if (!is.null(partial_auc) &&
      (!is.numeric(partial_auc) || length(partial_auc) != 2L || any(!is.finite(partial_auc)) ||
       any(partial_auc < 0 | partial_auc > 1))) {
    stop("`partial_auc` must be NULL or a finite numeric vector of length two within 0-1.", call. = FALSE)
  }
  if (!is.numeric(zero_correction) || length(zero_correction) != 1L || is.na(zero_correction) ||
      zero_correction <= 0) {
    stop("`zero_correction` must be a positive number.", call. = FALSE)
  }

  dat <- .r4vn_diag_resolve_data(data)

  outcome_expr <- substitute(outcome)
  marker_exprs <- as.list(substitute(list(...)))[-1L]
  outcome_name <- .r4vn_diag_name(outcome_expr, "outcome")
  marker_names <- .r4vn_diag_names(marker_exprs)

  if (!is.null(vars)) {
    vr <- .r4vn_resolve_vars_input(vars, data = dat, arg = "vars", default_type = "auto", strict = TRUE)
    marker_names <- c(marker_names, vr$variable)
  }
  marker_names <- unique(marker_names)
  marker_names <- marker_names[marker_names != outcome_name]

  if (!length(marker_names)) {
    stop("At least one diagnostic test/marker must be supplied in `...` or `vars`.", call. = FALSE)
  }

  missing_names <- setdiff(c(outcome_name, marker_names), names(dat))
  if (length(missing_names)) {
    stop("Variables not found in `data`: ", paste(missing_names, collapse = ", "), call. = FALSE)
  }

  out_bin <- .r4vn_diag_binary(dat[[outcome_name]], positive = event, what = outcome_name)
  y_all <- out_bin$value
  outcome_label <- .r4vn_diag_var_label(dat[[outcome_name]], outcome_name)
  marker_labels <- stats::setNames(
    vapply(marker_names, function(nm) .r4vn_diag_var_label(dat[[nm]], nm), character(1)),
    marker_names
  )

  threshold_rows <- list()
  summary_rows <- list()
  descriptive_rows <- list()
  roc_objects <- list()
  coord_objects <- list()
  all_cut_rows <- list()
  notes <- character()

  tr_i <- sr_i <- ac_i <- dr_i <- 0L

  for (j in seq_along(marker_names)) {
    marker <- marker_names[j]
    marker_label <- unname(marker_labels[[marker]])
    x0 <- dat[[marker]]
    desc_ok <- !is.na(y_all) & !is.na(x0)
    binary_test <- .r4vn_diag_is_binary(x0)
    analysis_type <- if (binary_test) "Binary 2 x 2" else if (isTRUE(roc)) "ROC" else "Not analyzed"

    dr_i <- dr_i + 1L
    descriptive_rows[[dr_i]] <- data.frame(
      marker = marker,
      marker_label = marker_label,
      analysis = analysis_type,
      total = nrow(dat),
      complete = sum(desc_ok),
      missing = nrow(dat) - sum(desc_ok),
      cases = sum(y_all[desc_ok] == 1L),
      controls = sum(y_all[desc_ok] == 0L),
      event_prevalence = if (sum(desc_ok) > 0L) mean(y_all[desc_ok] == 1L) else NA_real_,
      stringsAsFactors = FALSE
    )

    if (binary_test) {
      pos_j <- .r4vn_diag_test_positive_for(positive, marker, j)
      tr_i <- tr_i + 1L
      rr <- .r4vn_diag_binary_row(
        marker = marker,
        y = y_all,
        test = x0,
        test_positive = pos_j,
        prevalence = prevalence,
        ci = ci,
        ci_level = ci_level,
        ci_method = ci_method,
        zero_correction = zero_correction
      )
      rr$marker_label <- marker_label
      rr$cutoff_low <- NA_real_
      rr$cutoff_high <- NA_real_
      threshold_rows[[tr_i]] <- rr
      next
    }

    if (!isTRUE(roc)) next

    x <- .r4vn_diag_prepare_predictor(x0, marker)
    ok <- !is.na(y_all) & !is.na(x)
    yy <- y_all[ok]
    xx <- x[ok]

    if (length(unique(yy)) != 2L) {
      notes <- c(notes, paste0(marker_label, ": ROC omitted because complete cases did not contain both outcome levels."))
      descriptive_rows[[dr_i]]$analysis <- "ROC omitted: one outcome level"
      next
    }
    if (length(unique(xx)) < 2L) {
      notes <- c(notes, paste0(marker_label, ": ROC omitted because the marker was constant among complete cases."))
      descriptive_rows[[dr_i]]$analysis <- "ROC omitted: constant marker"
      next
    }

    roc_obj <- tryCatch(
      .r4vn_diag_roc(yy, xx, direction = direction),
      error = function(e) e
    )
    if (inherits(roc_obj, "error")) {
      notes <- c(notes, paste0(marker_label, ": ROC omitted: ", conditionMessage(roc_obj)))
      descriptive_rows[[dr_i]]$analysis <- "ROC omitted"
      next
    }
    roc_objects[[marker]] <- roc_obj

    auc_val <- roc_obj$auc
    auc_var <- roc_obj$variance
    auc_se <- if (is.finite(auc_var) && auc_var >= 0) sqrt(auc_var) else NA_real_
    auc_z <- if (is.finite(auc_se) && auc_se > 0) (auc_val - 0.5) / auc_se else NA_real_
    auc_p <- if (is.finite(auc_z)) 2 * stats::pnorm(-abs(auc_z)) else NA_real_

    auc_low <- auc_high <- NA_real_
    if (.r4vn_diag_ci_requested(ci, "auc")) {
      auc_ci <- .r4vn_diag_auc_ci(roc_obj, level = ci_level)
      auc_low <- unname(auc_ci[1L])
      auc_high <- unname(auc_ci[2L])
    }

    p_auc <- p_auc_low <- p_auc_high <- NA_real_
    if (!is.null(partial_auc)) {
      p_auc <- .r4vn_diag_partial_auc(
        roc_obj,
        range = partial_auc,
        focus = partial_focus,
        correct = partial_correct
      )
      if (.r4vn_diag_ci_requested(ci, "auc")) {
        p_ci <- .r4vn_diag_partial_auc_ci(
          roc_obj,
          range = partial_auc,
          focus = partial_focus,
          correct = partial_correct,
          level = ci_level,
          boot = boot
        )
        p_auc_low <- unname(p_ci[1L])
        p_auc_high <- unname(p_ci[2L])
      }
    }

    sr_i <- sr_i + 1L
    summary_rows[[sr_i]] <- data.frame(
      marker = marker,
      marker_label = marker_label,
      n = length(yy),
      cases = sum(yy == 1L),
      controls = sum(yy == 0L),
      auc = auc_val,
      auc_se = auc_se,
      auc_z = auc_z,
      auc_p = auc_p,
      auc_low = auc_low,
      auc_high = auc_high,
      partial_auc = p_auc,
      partial_auc_low = p_auc_low,
      partial_auc_high = p_auc_high,
      direction = if (roc_obj$direction == "<") "higher = positive" else "lower = positive",
      stringsAsFactors = FALSE
    )

    coords <- roc_obj$coordinates
    coords$marker <- marker
    coords$marker_label <- marker_label
    coord_objects[[marker]] <- coords[, c("marker", "marker_label", "threshold", "sensitivity", "specificity"), drop = FALSE]

    bt <- .r4vn_diag_best_thresholds(roc_obj, best = best, target = target)

    if (nrow(bt)) {
      for (r in seq_len(nrow(bt))) {
        tr_i <- tr_i + 1L
        rr <- .r4vn_diag_threshold_row(
          marker = marker,
          criterion = bt$criterion[r],
          cutoff = bt$cutoff[r],
          direction = roc_obj$direction,
          y = yy,
          x = xx,
          prevalence = prevalence,
          ci = ci,
          ci_level = ci_level,
          ci_method = ci_method,
          zero_correction = zero_correction
        )
        rr$marker_label <- marker_label
        rr$cutoff_low <- NA_real_
        rr$cutoff_high <- NA_real_

        if (isTRUE(cut_ci)) {
          crit_low <- tolower(bt$criterion[r])
          method_key <- if (grepl("^youden", crit_low)) {
            "youden"
          } else if (grepl("^closest", crit_low)) {
            "closest"
          } else {
            NA_character_
          }
          if (!is.na(method_key)) {
            cci <- .r4vn_diag_cut_ci(
              roc_obj, method = method_key,
              level = ci_level, boot = boot
            )
            rr$cutoff_low <- unname(cci[1L])
            rr$cutoff_high <- unname(cci[2L])
          }
        }
        threshold_rows[[tr_i]] <- rr
      }
    }

    if (!is.null(cuts)) {
      for (cu in cuts) {
        if (!is.finite(cu)) next
        tr_i <- tr_i + 1L
        rr <- .r4vn_diag_threshold_row(
          marker = marker,
          criterion = "User specified",
          cutoff = cu,
          direction = roc_obj$direction,
          y = yy,
          x = xx,
          prevalence = prevalence,
          ci = ci,
          ci_level = ci_level,
          ci_method = ci_method,
          zero_correction = zero_correction
        )
        rr$marker_label <- marker_label
        rr$cutoff_low <- NA_real_
        rr$cutoff_high <- NA_real_
        threshold_rows[[tr_i]] <- rr
      }
    }

    if (isTRUE(all_cuts)) {
      finite_cuts <- unique(coords$threshold[is.finite(coords$threshold)])
      if (length(finite_cuts)) {
        rows <- vector("list", length(finite_cuts))
        for (r in seq_along(finite_cuts)) {
          rows[[r]] <- .r4vn_diag_threshold_row(
            marker = marker,
            criterion = "All ROC thresholds",
            cutoff = finite_cuts[r],
            direction = roc_obj$direction,
            y = yy,
            x = xx,
            prevalence = prevalence,
            ci = ci,
            ci_level = ci_level,
            ci_method = ci_method,
            zero_correction = zero_correction
          )
          rows[[r]]$marker_label <- marker_label
        }
        ac_i <- ac_i + 1L
        all_cut_rows[[ac_i]] <- do.call(rbind, rows)
      }
    }
  }

  thresholds <- if (length(threshold_rows)) do.call(rbind, threshold_rows) else data.frame()
  summary_df <- if (length(summary_rows)) do.call(rbind, summary_rows) else data.frame()
  descriptive_df <- if (length(descriptive_rows)) do.call(rbind, descriptive_rows) else data.frame()
  coordinates <- if (length(coord_objects)) do.call(rbind, coord_objects) else data.frame()
  all_cutoffs <- if (length(all_cut_rows)) do.call(rbind, all_cut_rows) else data.frame()

  comparison <- data.frame()
  roc_marker_names <- names(roc_objects)
  if (isTRUE(compare) && length(roc_marker_names) >= 2L) {
    cmb <- utils::combn(roc_marker_names, 2L, simplify = FALSE)
    comp_rows <- list()
    ci_comp <- 0L

    for (i in seq_along(cmb)) {
      m1 <- cmb[[i]][1L]
      m2 <- cmb[[i]][2L]
      x1 <- .r4vn_diag_prepare_predictor(dat[[m1]], m1)
      x2 <- .r4vn_diag_prepare_predictor(dat[[m2]], m2)
      ok <- !is.na(y_all) & !is.na(x1) & !is.na(x2)
      yy <- y_all[ok]
      xx1 <- x1[ok]
      xx2 <- x2[ok]

      if (length(unique(yy)) != 2L || length(unique(xx1)) < 2L || length(unique(xx2)) < 2L) {
        notes <- c(
          notes,
          paste0(marker_labels[[m1]], " vs ", marker_labels[[m2]],
                 ": ROC comparison omitted because the paired complete-case sample was not estimable.")
        )
        next
      }

      dir1 <- if (identical(direction, "auto")) roc_objects[[m1]]$direction else direction
      dir2 <- if (identical(direction, "auto")) roc_objects[[m2]]$direction else direction
      score1 <- if (dir1 == "<") xx1 else -xx1
      score2 <- if (dir2 == "<") xx2 else -xx2

      pair <- tryCatch(
        if (identical(compare_method, "bootstrap")) {
          .r4vn_diag_compare_bootstrap(
            yy, score1, score2,
            level = ci_level, boot = boot
          )
        } else {
          .r4vn_diag_compare_delong(
            yy, score1, score2,
            level = ci_level
          )
        },
        error = function(e) e
      )

      if (inherits(pair, "error")) {
        notes <- c(
          notes,
          paste0(marker_labels[[m1]], " vs ", marker_labels[[m2]],
                 ": ROC comparison omitted: ", conditionMessage(pair))
        )
        next
      }

      a1 <- pair$auc1
      a2 <- pair$auc2
      lo <- pair$lower
      hi <- pair$upper

      ci_comp <- ci_comp + 1L
      comp_rows[[ci_comp]] <- data.frame(
        marker1 = m1,
        marker1_label = unname(marker_labels[[m1]]),
        marker2 = m2,
        marker2_label = unname(marker_labels[[m2]]),
        n = sum(ok),
        auc1 = a1,
        auc2 = a2,
        difference = a1 - a2,
        lower = lo,
        upper = hi,
        p = as.numeric(pair$p),
        stringsAsFactors = FALSE
      )
    }

    if (length(comp_rows)) {
      comparison <- do.call(rbind, comp_rows)
      comparison$p_adjusted <- stats::p.adjust(comparison$p, method = adjust)
    }
  } else if (isTRUE(compare) && length(roc_marker_names) < 2L) {
    notes <- c(notes, "Pairwise ROC comparison was requested, but fewer than two estimable ROC markers were available.")
  }

  if (!is.null(prevalence)) {
    notes <- c(
      notes,
      paste0(
        "PPV and NPV were adjusted to a disease prevalence of ",
        formatC(prevalence * 100, format = "f", digits = 1),
        "%; observed predictive values remain available in the internal 2 x 2 calculations."
      )
    )
  }
  if (isTRUE(cut_ci)) {
    notes <- c(
      notes,
      paste0("Optimal-cutoff confidence intervals used stratified bootstrap with ", boot, " replicates.")
    )
  }
  if (isTRUE(compare) && nrow(comparison)) {
    notes <- c(
      notes,
      paste0(
        "Pairwise ROC comparisons used ", compare_method,
        " tests on pair-specific complete cases; adjusted p-values use `", adjust, "`."
      )
    )
  }

  settings <- list(
    outcome = outcome_name,
    outcome_label = outcome_label,
    event = out_bin$positive,
    negative = out_bin$negative,
    markers = marker_names,
    marker_labels = marker_labels,
    direction = direction,
    best = best,
    cuts = cuts,
    target = target,
    ci = ci,
    ci_level = ci_level,
    ci_method = ci_method,
    cut_ci = cut_ci,
    boot = boot,
    prevalence = prevalence,
    roc = roc,
    partial_auc = partial_auc,
    partial_focus = partial_focus,
    partial_correct = partial_correct,
    compare = compare,
    compare_method = compare_method,
    adjust = adjust,
    all_cuts = all_cuts,
    missing = missing,
    show = measures,
    viewer = viewer_show,
    console = console,
    hide = hide,
    plot = plot,
    plot_args = plot_args,
    interpretation = interpretation,
    digit = as.integer(digit),
    p_digit = as.integer(p_digit),
    title = title
  )

  out <- list(
    summary = summary_df,
    thresholds = thresholds,
    coordinates = coordinates,
    comparison = comparison,
    all_cutoffs = all_cutoffs,
    roc = roc_objects,
    descriptive = descriptive_df,
    settings = settings,
    notes = unique(notes),
    call = call
  )
  class(out) <- "r4vn_diag"
  out <- .r4vn_diag_contract(out, interpretation = interpretation)

  out$plot_titles <- list()
  if (!identical(plot, FALSE)) {
    what <- if (isTRUE(plot)) "roc" else as.character(plot)[1L]
    plot_types <- if (identical(what, "both")) c("roc", "cutoff") else what
    plot_types <- intersect(plot_types, c("roc", "cutoff"))

    if (length(out$roc) && length(plot_types)) {
      plot_dir <- tempfile(pattern = "r4vn-tabdiag-plots-")
      dir.create(plot_dir, recursive = TRUE, showWarnings = FALSE)
      save_args <- plot_args
      save_args[c("x", "what", "file", "width", "height", "res")] <- NULL

      for (tp in plot_types) {
        path <- file.path(plot_dir, paste0(tp, ".png"))
        ok_plot <- tryCatch({
          do.call(
            plot.r4vn_diag,
            c(
              list(x = out, what = tp, file = path, width = 7, height = 7, res = 144),
              save_args
            )
          )
          file.exists(path) && isTRUE(file.info(path)$size > 0)
        }, error = function(e) {
          warning(
            "Diagnostic plot `", tp, "` could not be added to the Viewer: ",
            conditionMessage(e), call. = FALSE
          )
          FALSE
        })
        if (isTRUE(ok_plot)) {
          out$plots[[tp]] <- path
          out$plot_titles[[tp]] <- if (tp == "roc") {
            "ROC curve"
          } else {
            "Sensitivity and specificity across thresholds"
          }
        }
      }

      if (interactive()) {
        pane_args <- plot_args
        pane_args[c("x", "what", "file")] <- NULL
        tryCatch(
          do.call(plot.r4vn_diag, c(list(x = out, what = what), pane_args)),
          error = function(e) warning(
            "Diagnostic Plot-pane figure was unavailable: ", conditionMessage(e),
            call. = FALSE
          )
        )
      }
    } else if (!length(out$roc)) {
      out$notes <- unique(c(out$notes, "Plot requested, but no estimable ROC curve was available."))
    }
  }

  if (!is.null(export) || !is.null(file)) {
    if (!exists("tabexport", mode = "function", inherits = TRUE)) {
      stop("`tabexport()` is required for direct export from `tabdiag()`.", call. = FALSE)
    }
    export_tables <- out$tables[vapply(out$tables, function(z) is.data.frame(z) && nrow(z) > 0L, logical(1))]
    if (!length(export_tables)) {
      stop("No diagnostic tables are available to export.", call. = FALSE)
    }
    out$export <- tabexport(
      export_tables, export = export, file = file, open = open,
      title = title %||% "R4VN diagnostic accuracy analysis"
    )
  }

  .r4vn_show(out, show = viewer_show, console = console)
}

.r4vn_diag_var_label <- function(x, fallback) {
  lab <- attr(x, "label", exact = TRUE)
  if (!is.null(lab) && length(lab)) {
    lab <- as.character(lab[1L])
    if (!is.na(lab) && nzchar(trimws(lab))) return(lab)
  }
  fallback
}

.r4vn_diag_ci_label <- function(level = 0.95) {
  paste0(formatC(100 * level, format = "fg", digits = 6), "% CI")
}

.r4vn_diag_comparison_display <- function(z, digit = 2, p_digit = 3, ci_level = 0.95) {
  if (is.null(z) || !nrow(z)) return(data.frame())
  lab1 <- if ("marker1_label" %in% names(z)) z$marker1_label else z$marker1
  lab2 <- if ("marker2_label" %in% names(z)) z$marker2_label else z$marker2
  out <- data.frame(
    Comparison = paste(lab1, "vs", lab2),
    N = z$n,
    `AUC 1` = .r4vn_diag_format_num(z$auc1, digit),
    `AUC 2` = .r4vn_diag_format_num(z$auc2, digit),
    `AUC difference` = .r4vn_diag_format_est_ci(
      z$difference, z$lower, z$upper, digit = digit, percent = FALSE
    ),
    `p-value` = .r4vn_diag_format_p(z$p, p_digit),
    `Adjusted p` = .r4vn_diag_format_p(z$p_adjusted, p_digit),
    stringsAsFactors = FALSE,
    check.names = FALSE
  )
  attr(out, "ci_label") <- .r4vn_diag_ci_label(ci_level)
  out
}

.r4vn_diag_confusion_display <- function(x) {
  if (is.null(x$thresholds) || !nrow(x$thresholds)) return(data.frame())
  z <- x$thresholds
  rows <- vector("list", 3L * nrow(z))
  k <- 0L
  for (i in seq_len(nrow(z))) {
    marker <- if ("marker_label" %in% names(z)) z$marker_label[i] else z$marker[i]
    cutoff <- if (is.na(z$cutoff[i])) "" else paste0(z$direction[i], .r4vn_diag_format_num(z$cutoff[i], x$settings$digit))
    tp <- as.integer(z$tp[i]); fp <- as.integer(z$fp[i])
    fn <- as.integer(z$fn[i]); tn <- as.integer(z$tn[i])

    k <- k + 1L
    rows[[k]] <- data.frame(
      Marker = marker,
      Criterion = z$criterion[i],
      Cutoff = cutoff,
      `Test result` = "Positive",
      `Reference positive` = tp,
      `Reference negative` = fp,
      Total = tp + fp,
      stringsAsFactors = FALSE, check.names = FALSE
    )
    k <- k + 1L
    rows[[k]] <- data.frame(
      Marker = marker,
      Criterion = z$criterion[i],
      Cutoff = cutoff,
      `Test result` = "Negative",
      `Reference positive` = fn,
      `Reference negative` = tn,
      Total = fn + tn,
      stringsAsFactors = FALSE, check.names = FALSE
    )
    k <- k + 1L
    rows[[k]] <- data.frame(
      Marker = marker,
      Criterion = z$criterion[i],
      Cutoff = cutoff,
      `Test result` = "Total",
      `Reference positive` = tp + fn,
      `Reference negative` = fp + tn,
      Total = tp + fp + fn + tn,
      stringsAsFactors = FALSE, check.names = FALSE
    )
  }
  do.call(rbind, rows)
}

.r4vn_diag_confusion_view_sections <- function(x, tab) {
  if (!is.data.frame(tab) || !nrow(tab)) return(character())
  ids <- intersect(c("Marker", "Criterion", "Cutoff"), names(tab))
  if (!length(ids)) return(.r4vn_view_section("2 x 2 diagnostic classification", tab))

  key <- interaction(tab[ids], drop = TRUE, lex.order = TRUE)
  groups <- split(seq_len(nrow(tab)), key)
  blocks <- character()

  for (ii in groups) {
    z <- tab[ii, , drop = FALSE]
    marker <- as.character(z$Marker[1L])
    criterion <- as.character(z$Criterion[1L])
    cutoff <- as.character(z$Cutoff[1L])

    body <- z[, c("Test result", "Reference positive", "Reference negative", "Total"), drop = FALSE]
    title <- paste0("2 x 2 diagnostic classification - ", marker)

    details <- character()
    if (nzchar(criterion)) {
      if (grepl("^Positive:", criterion, ignore.case = TRUE)) {
        details <- c(details, paste0("Test positive = ", trimws(sub("^Positive:", "", criterion, ignore.case = TRUE))))
      } else {
        details <- c(details, paste0("Criterion: ", criterion))
      }
    }
    if (nzchar(cutoff)) details <- c(details, paste0("Cutoff: ", cutoff))
    description <- if (length(details)) paste(details, collapse = "; ") else NULL

    blocks <- c(blocks, .r4vn_view_section(title, body, description = description))
  }
  blocks
}

.r4vn_diag_performance_long <- function(tab) {
  if (!is.data.frame(tab) || !nrow(tab)) return(data.frame())
  id <- intersect(c("Marker", "Criterion", "Cutoff"), names(tab))
  metrics <- setdiff(names(tab), id)
  if (!length(metrics)) return(tab)

  out <- list()
  k <- 0L
  for (i in seq_len(nrow(tab))) {
    vals <- as.character(tab[i, id, drop = TRUE])
    vals <- vals[!is.na(vals) & nzchar(trimws(vals))]
    label <- paste(vals, collapse = " | ")
    for (m in metrics) {
      k <- k + 1L
      out[[k]] <- data.frame(
        Result = label,
        Measure = m,
        Value = as.character(tab[[m]][i]),
        stringsAsFactors = FALSE
      )
    }
  }
  do.call(rbind, out)
}

.r4vn_diag_performance_display <- function(x) {
  s <- x$settings
  if (!nrow(x$thresholds)) return(data.frame())
  disp <- .r4vn_diag_threshold_display(
    x$thresholds,
    ci = s$ci,
    show = s$show,
    hide = s$hide,
    digit = s$digit
  )
  if (isTRUE(s$cut_ci) && any(!is.na(x$thresholds$cutoff_low))) {
    cutci <- rep("", nrow(x$thresholds))
    ok <- !is.na(x$thresholds$cutoff_low) & !is.na(x$thresholds$cutoff_high)
    cutci[ok] <- paste0(
      .r4vn_diag_format_num(x$thresholds$cutoff_low[ok], s$digit),
      "\u2013",
      .r4vn_diag_format_num(x$thresholds$cutoff_high[ok], s$digit)
    )
    loc <- match("Cutoff", names(disp))
    nm <- paste0("Cutoff ", .r4vn_diag_ci_label(s$ci_level))
    tmp <- list(cutci)
    names(tmp) <- nm
    disp <- append(disp, tmp, after = loc)
    disp <- as.data.frame(disp, stringsAsFactors = FALSE, check.names = FALSE)
  }
  disp
}

.r4vn_diag_interpretation <- function(x) {
  out <- list()
  k <- 0L
  add <- function(section, marker, text) {
    k <<- k + 1L
    out[[k]] <<- data.frame(
      section = section,
      marker = marker,
      interpretation = text,
      stringsAsFactors = FALSE
    )
  }

  if (nrow(x$summary)) {
    for (i in seq_len(nrow(x$summary))) {
      z <- x$summary[i, , drop = FALSE]
      auc <- z$auc[1L]
      quality <- if (!is.finite(auc)) {
        "could not be classified"
      } else if (auc < 0.60) {
        "showed little discrimination"
      } else if (auc < 0.70) {
        "showed limited discrimination"
      } else if (auc < 0.80) {
        "showed acceptable discrimination"
      } else if (auc < 0.90) {
        "showed good discrimination"
      } else {
        "showed excellent discrimination"
      }
      ci_txt <- ""
      if (is.finite(z$auc_low[1L]) && is.finite(z$auc_high[1L])) {
        ci_txt <- paste0(
          " (", .r4vn_diag_ci_label(x$settings$ci_level), " ",
          .r4vn_diag_format_num(z$auc_low[1L], x$settings$digit), "\u2013",
          .r4vn_diag_format_num(z$auc_high[1L], x$settings$digit), ")"
        )
      }
      add(
        "Discrimination", z$marker_label[1L],
        paste0(
          z$marker_label[1L], " ", quality, " with AUC ",
          .r4vn_diag_format_num(auc, x$settings$digit), ci_txt,
          ". This describes discrimination in the analyzed sample, not clinical utility by itself."
        )
      )
    }
  }

  if (nrow(x$thresholds)) {
    for (m in unique(x$thresholds$marker)) {
      zz_all <- x$thresholds[x$thresholds$marker == m, , drop = FALSE]
      zz <- zz_all[!is.na(zz_all$cutoff), , drop = FALSE]
      if (!nrow(zz)) {
        z <- zz_all[1L, , drop = FALSE]
        add(
          "Binary test", z$marker_label[1L],
          paste0(
            z$marker_label[1L], " had sensitivity ",
            .r4vn_diag_format_num(100 * z$sensitivity[1L], x$settings$digit),
            "% and specificity ",
            .r4vn_diag_format_num(100 * z$specificity[1L], x$settings$digit),
            "%. Predictive values depend on disease prevalence, and these results require validation in the intended population."
          )
        )
        next
      }
      pick <- which(grepl("^Youden", zz$criterion))
      if (!length(pick)) pick <- which(grepl("^Closest", zz$criterion))
      if (!length(pick)) pick <- 1L
      z <- zz[pick[1L], , drop = FALSE]
      add(
        "Selected cutoff", z$marker_label[1L],
        paste0(
          "At the ", z$criterion[1L], " cutoff (", z$direction[1L],
          .r4vn_diag_format_num(z$cutoff[1L], x$settings$digit), "), sensitivity was ",
          .r4vn_diag_format_num(100 * z$sensitivity[1L], x$settings$digit),
          "% and specificity was ",
          .r4vn_diag_format_num(100 * z$specificity[1L], x$settings$digit),
          "%. Cutoff performance should be externally validated before clinical adoption."
        )
      )
    }
  }

  if (nrow(x$comparison)) {
    for (i in seq_len(nrow(x$comparison))) {
      z <- x$comparison[i, , drop = FALSE]
      p <- z$p_adjusted[1L]
      conclusion <- if (is.finite(p) && p < 0.05) {
        "The paired AUCs differed statistically after the requested multiplicity adjustment."
      } else {
        "The paired analysis did not show a statistically detectable AUC difference after the requested multiplicity adjustment."
      }
      add(
        "ROC comparison",
        paste(z$marker1_label[1L], "vs", z$marker2_label[1L]),
        paste0(
          z$marker1_label[1L], " vs ", z$marker2_label[1L], ": ", conclusion,
          " Adjusted p = ", .r4vn_diag_format_p(p, x$settings$p_digit), "."
        )
      )
    }
  }

  if (!length(out)) return(data.frame())
  do.call(rbind, out)
}

.r4vn_diag_contract <- function(x, interpretation = FALSE) {
  s <- x$settings
  x$estimates <- list(
    auc = x$summary,
    performance = x$thresholds
  )
  x$tests <- list(
    auc_vs_0.5 = if (nrow(x$summary)) {
      x$summary[, intersect(
        c("marker", "marker_label", "n", "auc", "auc_se", "auc_z", "auc_p", "auc_low", "auc_high"),
        names(x$summary)
      ), drop = FALSE]
    } else data.frame(),
    roc_comparison = x$comparison
  )
  x$diagnostics <- list(
    coordinates = x$coordinates,
    all_cutoffs = x$all_cutoffs
  )
  x$models <- list(roc = x$roc)
  x["interpretation"] <- list(if (isTRUE(interpretation)) .r4vn_diag_interpretation(x) else NULL)

  sample_disp <- x$descriptive
  if (nrow(sample_disp)) {
    sample_disp <- data.frame(
      Marker = sample_disp$marker_label,
      Analysis = sample_disp$analysis,
      Total = sample_disp$total,
      Complete = sample_disp$complete,
      Missing = sample_disp$missing,
      Cases = sample_disp$cases,
      Controls = sample_disp$controls,
      `Outcome prevalence` = ifelse(
        is.na(sample_disp$event_prevalence), "",
        paste0(.r4vn_diag_format_num(100 * sample_disp$event_prevalence, s$digit), "%")
      ),
      stringsAsFactors = FALSE,
      check.names = FALSE
    )
  }

  roc_disp <- if (nrow(x$summary)) {
    .r4vn_diag_summary_display(x$summary, ci = s$ci, digit = s$digit, p_digit = s$p_digit)
  } else data.frame()
  perf_disp <- .r4vn_diag_performance_display(x)
  confusion_disp <- .r4vn_diag_confusion_display(x)
  comp_disp <- .r4vn_diag_comparison_display(
    x$comparison, digit = s$digit, p_digit = s$p_digit, ci_level = s$ci_level
  )

  x$tables <- list(
    Sample = sample_disp,
    ROC_summary = roc_disp,
    Diagnostic_performance = perf_disp,
    Confusion_matrix = confusion_disp,
    ROC_comparison = comp_disp
  )
  if (!is.null(x$interpretation) && nrow(x$interpretation)) {
    x$tables$Interpretation <- x$interpretation
  }
  x$plots <- list()
  x$metadata <- list(
    outcome = s$outcome,
    outcome_label = s$outcome_label,
    positive = s$event,
    negative = s$negative,
    markers = s$markers,
    marker_labels = s$marker_labels,
    ci_level = s$ci_level,
    prevalence = s$prevalence,
    pairwise_complete_case_comparison = isTRUE(s$compare),
    missing_table_shown = isTRUE(s$missing),
    interpretation = isTRUE(interpretation)
  )
  x
}

#' @export
print.r4vn_diag <- function(x, ...) {
  s <- x$settings
  if (!is.null(s$title) && nzchar(s$title)) {
    cat("\n", s$title, "\n", sep = "")
  } else {
    cat("\nR4VN diagnostic accuracy analysis\n")
  }
  cat(strrep("-", 58), "\n", sep = "")
  cat("Outcome: ", s$outcome_label %||% s$outcome, "\n", sep = "")
  cat("Positive outcome: ", s$event, "\n", sep = "")

  if (isTRUE(s$missing) && !is.null(x$tables$Sample) && nrow(x$tables$Sample)) {
    cat("\nAnalysis sample by diagnostic marker\n")
    .r4vn_diag_print_df(x$tables$Sample)
  }

  if (!is.null(x$tables$ROC_summary) && nrow(x$tables$ROC_summary)) {
    cat("\nROC summary\n")
    .r4vn_diag_print_df(x$tables$ROC_summary)
  }

  binary_markers <- if (!is.null(x$descriptive) && nrow(x$descriptive)) {
    x$descriptive$marker[x$descriptive$analysis == "Binary 2 x 2"]
  } else character()
  counts_profile <- length(s$show) == 1L && is.character(s$show) &&
    tolower(s$show) %in% c("counts", "cells")
  show_all_cells <- length(s$show) == 1L && is.character(s$show) &&
    tolower(s$show) %in% c("all", "counts", "cells")

  if (!is.null(x$tables$Confusion_matrix) && nrow(x$tables$Confusion_matrix) &&
      (isTRUE(show_all_cells) || length(binary_markers))) {
    cm <- x$tables$Confusion_matrix
    if (!isTRUE(show_all_cells) && length(binary_markers)) {
      cm <- cm[cm$Marker %in% unname(x$settings$marker_labels[binary_markers]), , drop = FALSE]
    }
    if (nrow(cm)) {
      cat("\n2 x 2 diagnostic classification\n")
      .r4vn_diag_print_df(cm)
    }
  }

  if (!isTRUE(counts_profile) && !is.null(x$tables$Diagnostic_performance) &&
      nrow(x$tables$Diagnostic_performance)) {
    cat("\nDiagnostic performance at selected thresholds/tests\n")
    .r4vn_diag_print_df(x$tables$Diagnostic_performance)
  }

  if (!is.null(x$tables$ROC_comparison) && nrow(x$tables$ROC_comparison)) {
    cat("\nPairwise ROC comparison\n")
    .r4vn_diag_print_df(x$tables$ROC_comparison)
  }

  if (isTRUE(s$all_cuts) && nrow(x$all_cutoffs)) {
    cat(
      "\nFull metrics for every empirical ROC threshold are stored in ",
      "`$all_cutoffs` (", nrow(x$all_cutoffs), " rows).\n",
      sep = ""
    )
  }

  if (!is.null(x$interpretation) && nrow(x$interpretation)) {
    cat("\nInterpretation\n")
    for (i in seq_len(nrow(x$interpretation))) {
      cat("- ", x$interpretation$section[i], " [", x$interpretation$marker[i], "]: ",
          x$interpretation$interpretation[i], "\n", sep = "")
    }
  }

  if (length(x$notes)) {
    cat("\nNotes:\n")
    for (z in x$notes) cat("- ", z, "\n", sep = "")
  }

  invisible(x)
}

#' Plot an R4VN diagnostic analysis
#'
#' Draw ROC curves and/or sensitivity-specificity curves over empirical
#' thresholds from an object returned by `tabdiag()`.
#'
#' @param x An object returned by `tabdiag()`.
#' @param what `"roc"`, `"cutoff"`, or `"both"`.
#' @param color Optional vector of line colors. Defaults to a distinct base-R
#'   qualitative palette.
#' @param lty ROC line type(s).
#' @param line_width ROC/cutoff line width(s).
#' @param legend Logical; show the ROC legend.
#' @param legend_position Base-graphics legend position, default
#'   `"bottomright"`.
#' @param auc Logical; append AUC to ROC legend labels.
#' @param diagonal Logical; draw the no-discrimination diagonal.
#' @param diagonal_lty Line type for the no-discrimination diagonal.
#' @param grid Logical; draw a light reference grid.
#' @param xlim,ylim ROC axis limits. Values are clamped to `c(0, 1)`. Defaults are
#'   exactly `c(0, 1)` so the ROC axes never extend below 0 or above 1.
#' @param xlab,ylab ROC axis labels.
#' @param main Optional ROC title.
#' @param cutoff_mark Logical; mark selected cutoffs with vertical reference
#'   lines on cutoff plots.
#' @param cutoff_legend Logical; show sensitivity/specificity legend on cutoff
#'   plots.
#' @param cutoff_xlab,cutoff_ylab Cutoff-plot axis labels.
#' @param cutoff_main Optional cutoff-plot title. When multiple markers are
#'   plotted, the marker label is appended automatically.
#' @param file Optional graphics filename. Supported extensions are `.png`,
#'   `.jpg`/`.jpeg`, `.tif`/`.tiff`, `.pdf`, and `.svg`. For `what = "both"`,
#'   save the ROC and cutoff plots separately with two calls to `plot()`.
#' @param width,height Figure width and height in inches when `file` is used.
#'   Defaults are 7 and 7.
#' @param res Raster resolution in dpi for PNG/JPEG/TIFF output. Default 300.
#' @param bg Graphics-device background color. Default `"white"`.
#' @param font_family Optional base-graphics font family, for example
#'   `"Arial"`, when available on the current system.
#' @param cex_axis,cex_lab,cex_main Text-size controls for axes, labels, title.
#' @param bty Box type passed to base graphics.
#' @param ... Additional arguments passed to the initial base `plot()` call.
#'
#' @return Invisibly returns `x`.
#' @export
plot.r4vn_diag <- function(
  x,
  what = c("roc", "cutoff", "both"),
  color = NULL,
  lty = 1,
  line_width = 2,
  legend = TRUE,
  legend_position = "bottomright",
  auc = TRUE,
  diagonal = TRUE,
  diagonal_lty = 2,
  grid = FALSE,
  xlim = c(0, 1),
  ylim = c(0, 1),
  xlab = "1 - Specificity",
  ylab = "Sensitivity",
  main = NULL,
  cutoff_mark = TRUE,
  cutoff_legend = TRUE,
  cutoff_xlab = "Threshold",
  cutoff_ylab = "Probability",
  cutoff_main = NULL,
  file = NULL,
  width = 7,
  height = 7,
  res = 300,
  bg = "white",
  font_family = NULL,
  cex_axis = 1,
  cex_lab = 1,
  cex_main = 1,
  bty = "l",
  ...
) {
  what <- match.arg(what)

  if (!is.null(file)) {
    if (!is.character(file) || length(file) != 1L || is.na(file) || !nzchar(file)) {
      stop("`file` must be NULL or a single non-empty filename.", call. = FALSE)
    }
    if (identical(what, "both")) {
      stop(
        "When `file` is supplied, use separate plot() calls for `what = 'roc'` and `what = 'cutoff'`.",
        call. = FALSE
      )
    }
    if (!is.numeric(width) || length(width) != 1L || !is.finite(width) || width <= 0 ||
        !is.numeric(height) || length(height) != 1L || !is.finite(height) || height <= 0) {
      stop("`width` and `height` must be positive finite numbers in inches.", call. = FALSE)
    }
    if (!is.numeric(res) || length(res) != 1L || !is.finite(res) || res <= 0) {
      stop("`res` must be a positive finite number.", call. = FALSE)
    }
  }

  if (!length(x$roc)) {
    stop("No ROC curves are available in this object.", call. = FALSE)
  }
  if (!is.numeric(xlim) || length(xlim) != 2L || any(!is.finite(xlim)) ||
      !is.numeric(ylim) || length(ylim) != 2L || any(!is.finite(ylim))) {
    stop("`xlim` and `ylim` must each contain two finite numeric values.", call. = FALSE)
  }
  xlim <- pmax(0, pmin(1, xlim))
  ylim <- pmax(0, pmin(1, ylim))
  if (diff(range(xlim)) <= 0 || diff(range(ylim)) <= 0) {
    stop("ROC `xlim` and `ylim` must span a positive range within 0-1.", call. = FALSE)
  }
  xlim <- sort(xlim)
  ylim <- sort(ylim)

  if (!is.null(file)) {
    file <- path.expand(file)
    out_dir <- dirname(file)
    if (!dir.exists(out_dir)) {
      ok_dir <- dir.create(out_dir, recursive = TRUE, showWarnings = FALSE)
      if (!ok_dir && !dir.exists(out_dir)) {
        stop("Could not create graphics output directory: ", out_dir, call. = FALSE)
      }
    }
    ext <- tolower(tools::file_ext(file))
    if (!nzchar(ext)) {
      stop("`file` must include a supported graphics extension.", call. = FALSE)
    }
    if (ext == "png") {
      grDevices::png(file, width = width, height = height, units = "in", res = res, bg = bg)
    } else if (ext %in% c("jpg", "jpeg")) {
      grDevices::jpeg(file, width = width, height = height, units = "in", res = res, bg = bg)
    } else if (ext %in% c("tif", "tiff")) {
      grDevices::tiff(file, width = width, height = height, units = "in", res = res, bg = bg, compression = "lzw")
    } else if (ext == "pdf") {
      grDevices::pdf(file, width = width, height = height, bg = bg, onefile = TRUE)
    } else if (ext == "svg") {
      grDevices::svg(file, width = width, height = height, bg = bg)
    } else {
      stop("Unsupported graphics extension: .", ext, call. = FALSE)
    }
    on.exit(grDevices::dev.off(), add = TRUE)
  }

  if (!is.null(font_family)) {
    if (!is.character(font_family) || length(font_family) != 1L || is.na(font_family)) {
      stop("`font_family` must be NULL or a single character value.", call. = FALSE)
    }
    old_family <- graphics::par("family")
    graphics::par(family = font_family)
    on.exit(graphics::par(family = old_family), add = TRUE)
  }

  nm <- names(x$roc)
  labels <- x$settings$marker_labels[nm]
  labels[is.na(labels) | !nzchar(labels)] <- nm[is.na(labels) | !nzchar(labels)]
  n_curve <- length(nm)
  if (is.null(color)) color <- grDevices::hcl.colors(max(1L, n_curve), "Dark 3")
  color <- rep(color, length.out = n_curve)
  lty <- rep(lty, length.out = n_curve)
  line_width <- rep(line_width, length.out = n_curve)

  draw_roc <- function() {
    graphics::plot(
      NA_real_, NA_real_, type = "n",
      xlim = xlim, ylim = ylim, xaxs = "i", yaxs = "i",
      xlab = xlab, ylab = ylab,
      main = main %||% "ROC curve",
      cex.axis = cex_axis, cex.lab = cex_lab, cex.main = cex_main,
      bty = bty, ...
    )
    if (isTRUE(grid)) graphics::grid()
    if (isTRUE(diagonal)) graphics::abline(0, 1, lty = diagonal_lty)

    for (i in seq_along(nm)) {
      z <- x$coordinates[x$coordinates$marker == nm[i], , drop = FALSE]
      if (!nrow(z)) next
      fpr <- 1 - z$specificity
      ord <- order(fpr, z$sensitivity)
      graphics::lines(
        fpr[ord], z$sensitivity[ord],
        col = color[i], lty = lty[i], lwd = line_width[i]
      )
    }

    if (isTRUE(legend)) {
      leg <- unname(labels)
      if (isTRUE(auc)) {
        auc_lookup <- stats::setNames(x$summary$auc, x$summary$marker)
        leg <- paste0(
          leg, " (AUC=",
          vapply(nm, function(m) {
            a <- auc_lookup[[m]]
            if (is.null(a) || !is.finite(a)) "NA" else .r4vn_diag_format_num(a, x$settings$digit)
          }, character(1)),
          ")"
        )
      }
      graphics::legend(
        legend_position, legend = leg,
        col = color, lty = lty, lwd = line_width,
        bty = "n"
      )
    }
  }

  draw_cutoff <- function() {
    old <- NULL
    if (n_curve > 1L) {
      old <- graphics::par(mfrow = c(n_curve, 1L))
      on.exit(graphics::par(old), add = TRUE)
    }

    for (i in seq_along(nm)) {
      m <- nm[i]
      z <- x$coordinates[x$coordinates$marker == m & is.finite(x$coordinates$threshold), , drop = FALSE]
      if (!nrow(z)) next
      z <- z[order(z$threshold), , drop = FALSE]
      ttl <- cutoff_main
      if (is.null(ttl)) ttl <- paste("Sensitivity and specificity -", labels[[i]])
      else if (n_curve > 1L) ttl <- paste(ttl, "-", labels[[i]])

      graphics::plot(
        z$threshold, z$sensitivity,
        type = "l", ylim = c(0, 1), yaxs = "i",
        xlab = cutoff_xlab, ylab = cutoff_ylab, main = ttl,
        lwd = line_width[i], col = color[i],
        cex.axis = cex_axis, cex.lab = cex_lab, cex.main = cex_main,
        bty = bty, ...
      )
      if (isTRUE(grid)) graphics::grid()
      graphics::lines(
        z$threshold, z$specificity,
        lty = 2, lwd = line_width[i], col = color[i]
      )

      if (isTRUE(cutoff_mark) && nrow(x$thresholds)) {
        sel <- x$thresholds[
          x$thresholds$marker == m & !is.na(x$thresholds$cutoff) &
            x$thresholds$criterion != "All ROC thresholds",
          , drop = FALSE
        ]
        if (nrow(sel)) {
          for (v in unique(sel$cutoff[is.finite(sel$cutoff)])) {
            graphics::abline(v = v, lty = 3)
          }
        }
      }
      if (isTRUE(cutoff_legend)) {
        graphics::legend(
          "bottomright",
          legend = c("Sensitivity", "Specificity"),
          col = rep(color[i], 2L), lty = c(1, 2),
          lwd = rep(line_width[i], 2L), bty = "n"
        )
      }
    }
  }

  if (what == "roc") draw_roc()
  if (what == "cutoff") draw_cutoff()
  if (what == "both") {
    draw_roc()
    draw_cutoff()
  }

  invisible(x)
}

#' @export
as.data.frame.r4vn_diag <- function(x, row.names = NULL, optional = FALSE, ...) {
  if (nrow(x$thresholds)) return(x$thresholds)
  x$summary
}


# ============================================================================
# Viewer renderer for tabdiag()
# ============================================================================
.r4vn_viewer_tabdiag <- function(x) {
  s <- x$settings
  blocks <- character()

  meta <- paste0(
    '<div class="r4vn-meta">',
    '<span class="r4vn-chip">Outcome: ', .r4vn_view_escape(s$outcome_label %||% s$outcome), '</span>',
    '<span class="r4vn-chip">Positive: ', .r4vn_view_escape(s$event), '</span>',
    '<span class="r4vn-chip">CI: ', .r4vn_view_escape(.r4vn_diag_ci_label(s$ci_level)), '</span>',
    if (!is.null(s$prevalence)) paste0(
      '<span class="r4vn-chip">Target prevalence: ',
      formatC(100 * s$prevalence, format = "f", digits = 1), '%</span>'
    ) else '',
    '</div>'
  )

  if (isTRUE(s$missing) && !is.null(x$tables$Sample) && nrow(x$tables$Sample)) {
    blocks <- c(blocks, .r4vn_view_section("Analysis sample by diagnostic marker", x$tables$Sample))
  }
  if (!is.null(x$tables$ROC_summary) && nrow(x$tables$ROC_summary)) {
    blocks <- c(blocks, .r4vn_view_section("ROC summary", x$tables$ROC_summary))
  }
  binary_markers <- if (!is.null(x$descriptive) && nrow(x$descriptive)) {
    x$descriptive$marker[x$descriptive$analysis == "Binary 2 x 2"]
  } else character()
  counts_profile <- length(s$show) == 1L && is.character(s$show) &&
    tolower(s$show) %in% c("counts", "cells")
  show_all_cells <- length(s$show) == 1L && is.character(s$show) &&
    tolower(s$show) %in% c("all", "counts", "cells")

  # A conventional 2 x 2 table is easier to interpret than listing TP/FP/FN/TN
  # as separate measures. Binary tests receive this table by default. For
  # continuous markers it is added when the user requests all/counts.
  if (!is.null(x$tables$Confusion_matrix) && nrow(x$tables$Confusion_matrix) &&
      (isTRUE(show_all_cells) || length(binary_markers))) {
    cm <- x$tables$Confusion_matrix
    if (!isTRUE(show_all_cells) && length(binary_markers)) {
      binary_labels <- unname(x$settings$marker_labels[binary_markers])
      cm <- cm[cm$Marker %in% binary_labels, , drop = FALSE]
    }
    if (nrow(cm)) blocks <- c(blocks, .r4vn_diag_confusion_view_sections(x, cm))
  }

  if (!isTRUE(counts_profile) && !is.null(x$tables$Diagnostic_performance) &&
      nrow(x$tables$Diagnostic_performance)) {
    perf <- x$tables$Diagnostic_performance
    if (nrow(perf) <= 4L && ncol(perf) > 12L) {
      long_perf <- .r4vn_diag_performance_long(perf)
      for (lab in unique(long_perf$Result)) {
        zz <- long_perf[long_perf$Result == lab, c("Measure", "Value"), drop = FALSE]
        blocks <- c(
          blocks,
          .r4vn_view_section(paste0("Diagnostic performance - ", lab), zz)
        )
      }
    } else {
      blocks <- c(
        blocks,
        .r4vn_view_section(
          "Diagnostic performance at selected thresholds/tests",
          perf
        )
      )
    }
  }
  if (!is.null(x$tables$ROC_comparison) && nrow(x$tables$ROC_comparison)) {
    blocks <- c(blocks, .r4vn_view_section("Pairwise ROC comparison", x$tables$ROC_comparison))
  }
  if (!is.null(x$interpretation) && nrow(x$interpretation)) {
    blocks <- c(blocks, .r4vn_view_section("Interpretation", x$interpretation))
  }

  notes <- x$notes
  if (isTRUE(s$all_cuts) && nrow(x$all_cutoffs)) {
    notes <- c(
      notes,
      paste0(
        "Full metrics for every empirical ROC threshold are stored in `$all_cutoffs` (",
        nrow(x$all_cutoffs), " rows)."
      )
    )
  }

  # Keep every figure requested in tabdiag() inside the same Viewer report.
  # The image files are copied beside index.html, so the Viewer does not depend
  # on browser access to an unrelated temporary directory or on extra packages.
  report_dir <- tempfile(pattern = "r4vn-tabdiag-viewer-")
  dir.create(report_dir, recursive = TRUE, showWarnings = FALSE)
  if (length(x$plots)) {
    for (nm in names(x$plots)) {
      src <- x$plots[[nm]]
      if (!is.character(src) || length(src) != 1L || !file.exists(src)) next
      ext <- tools::file_ext(src)
      dest_name <- paste0("plot-", nm, if (nzchar(ext)) paste0(".", ext) else "")
      dest <- file.path(report_dir, dest_name)
      copied <- isTRUE(file.copy(src, dest, overwrite = TRUE))
      if (!copied || !file.exists(dest)) next
      ttl <- x$plot_titles[[nm]] %||% if (nm == "roc") "ROC curve" else "Diagnostic plot"
      blocks <- c(
        blocks,
        paste0(
          '<section class="r4vn-section"><h2>', .r4vn_view_escape(ttl), '</h2>',
          '<div class="r4vn-diag-figure"><img src="', .r4vn_view_escape(dest_name),
          '" alt="', .r4vn_view_escape(ttl), '"></div></section>'
        )
      )
    }
  }

  title <- s$title %||% "Diagnostic accuracy analysis"
  subtitle <- "ROC and diagnostic performance"
  html <- paste0(
    '<!DOCTYPE html><html><head><meta charset="UTF-8">',
    '<meta name="viewport" content="width=device-width,initial-scale=1">',
    '<style>', .r4vn_view_css(),
    '.r4vn-diag-figure{text-align:center;margin:8px 0 2px;}',
    '.r4vn-diag-figure img{max-width:100%;height:auto;display:inline-block;}',
    '</style></head><body><main class="r4vn-page">',
    '<div class="r4vn-brand">R4VN</div><h1>', .r4vn_view_escape(title), '</h1>',
    '<div class="r4vn-subtitle">', .r4vn_view_escape(subtitle), '</div>',
    meta, paste0(blocks, collapse = ""), .r4vn_view_notes(notes),
    '</main></body></html>'
  )
  file <- file.path(report_dir, "index.html")
  writeLines(enc2utf8(html), file, useBytes = TRUE)
  list(html = html, file = file)
}