R/epi.R

In R4VN: Health Data Analysis and Publication-Ready Reporting

# R4VN epidemiological commands using individual-level data
#
# epi() is the main 2 x 2 and stratified command. cc(), cs(), mcc(), and ir()
# remain separate user-facing functions/aliases for familiar workflows.

# ----------------------------------------------------------------------------
# Viewer helpers
# ----------------------------------------------------------------------------

# Render the 2 x 2 table with explicit exposure and outcome variable headings.
# Example:
#
#                         Low birth weight [nhecan]
# Sex [gioi]              No              Yes          Total
# Male                    ...             ...          ...
# Female                  ...             ...          ...
# Total                   ...             ...          ...
#
# The first header identifies the row/exposure variable. The spanning header
# identifies the column/outcome variable. This is used only for epi()/cc()/cs()
# when variable metadata are available; immediate epii() keeps the generic
# table renderer.
.r4vn_epi_view_2x2 <- function(x, meta = NULL) {
  if (is.null(x)) return("")

  if (is.table(x)) x <- as.matrix(x)

  if (is.matrix(x)) {
    rn <- rownames(x)
    cn <- colnames(x)
    values <- x
  } else if (is.data.frame(x)) {
    rn <- rownames(x)
    cn <- names(x)
    values <- as.matrix(x)
  } else {
    return(.r4vn_view_table(x))
  }

  if (is.null(rn)) rn <- rep("", nrow(values))
  if (is.null(cn)) cn <- paste0("Column ", seq_len(ncol(values)))

  # A standard epii() display has two outcome columns plus Total.
  # If the structure is unexpected, fall back to the shared renderer.
  total_column <- length(cn) >= 1L && identical(tolower(tail(cn, 1L)), "total")
  outcome_n <- if (total_column) length(cn) - 1L else length(cn)

  if (outcome_n < 1L) return(.r4vn_view_table(x))

  exposure_heading <- if (!is.null(meta$exposure_display) && nzchar(meta$exposure_display)) {
    meta$exposure_display
  } else {
    "Exposure"
  }

  outcome_heading <- if (!is.null(meta$outcome_display) && nzchar(meta$outcome_display)) {
    meta$outcome_display
  } else {
    "Outcome"
  }

  first_header <- paste0(
    '<tr>',
    '<th class="stub" rowspan="2" style="vertical-align:bottom;text-align:left;">',
    .r4vn_view_escape(exposure_heading),
    '</th>',
    '<th class="value" colspan="', outcome_n,
    '" style="text-align:center;">',
    .r4vn_view_escape(outcome_heading),
    '</th>',
    if (total_column) '<th class="value" rowspan="2" style="vertical-align:bottom;text-align:right;">Total</th>' else '',
    '</tr>'
  )

  second_header <- paste0(
    '<tr>',
    paste0(
      '<th class="value">',
      .r4vn_view_escape(cn[seq_len(outcome_n)]),
      '</th>',
      collapse = ""
    ),
    '</tr>'
  )

  rows <- character(nrow(values))

  for (i in seq_len(nrow(values))) {
    cells <- character(ncol(values))

    for (j in seq_len(ncol(values))) {
      val <- values[i, j]
      val <- if (is.na(val)) "" else as.character(val)
      cells[j] <- paste0(
        '<td class="value">',
        .r4vn_view_escape(val),
        '</td>'
      )
    }

    rows[i] <- paste0(
      '<tr>',
      '<td class="stub">', .r4vn_view_escape(rn[i]), '</td>',
      paste0(cells, collapse = ""),
      '</tr>'
    )
  }

  paste0(
    '<div class="r4vn-table-scroll">',
    '<table class="r4vn-result-table">',
    '<thead>', first_header, second_header, '</thead>',
    '<tbody>', paste0(rows, collapse = ""), '</tbody>',
    '</table>',
    '</div>'
  )
}

.r4vn_epi_view_section_2x2 <- function(title, object, meta = NULL) {
  ttl <- if (is.null(title) || !nzchar(title)) {
    ""
  } else {
    paste0('<h2>', .r4vn_view_escape(title), '</h2>')
  }

  paste0(
    '<section class="r4vn-section epi-2x2-section">',
    ttl,
    .r4vn_epi_view_2x2(object, meta),
    '</section>'
  )
}

.r4vn_epi_viewer <- function(x, subtitle = "Epidemiological analysis using individual-level data") {
  blocks <- character()

  preferred <- c(
    "2 x 2 table", "Epidemiological measures", "Association tests",
    "Stratum-specific estimates", "Adjusted and interaction analysis",
    "Matched-pair table", "Matched estimate", "Rates", "Rate measures"
  )

  nms <- c(
    intersect(preferred, names(x$sections)),
    setdiff(names(x$sections), preferred)
  )

  for (nm in nms) {
    if (identical(nm, "2 x 2 table") && !is.null(x$epi_meta)) {
      blocks <- c(
        blocks,
        .r4vn_epi_view_section_2x2(
          nm,
          x$sections[[nm]],
          x$epi_meta
        )
      )
    } else {
      blocks <- c(
        blocks,
        .r4vn_view_section(nm, x$sections[[nm]])
      )
    }
  }

  .r4vn_view_document(
    x$title,
    paste0(blocks, collapse = ""),
    notes = x$notes,
    subtitle = subtitle,
    prefix = "r4vn-epi-"
  )
}

.r4vn_epi_display <- function(x, call, show = TRUE, console = FALSE) {
  x$call <- call
  .r4vn_show(x, show = show, console = console)
}

.r4vn_viewer_epi <- function(x) .r4vn_epi_viewer(x, "2 x 2 and stratified epidemiological analysis")
.r4vn_viewer_cc <- function(x) .r4vn_epi_viewer(x, "Case-control analysis")
.r4vn_viewer_cs <- function(x) .r4vn_epi_viewer(x, "Cross-sectional / cohort 2 x 2 analysis")
.r4vn_viewer_mcc <- function(x) .r4vn_epi_viewer(x, "Matched case-control analysis")
.r4vn_viewer_ir <- function(x) .r4vn_epi_viewer(x, "Incidence-rate comparison")

# ----------------------------------------------------------------------------
# Helpers for variable names and variable labels
# ----------------------------------------------------------------------------

.r4vn_epi_expression_name <- function(expr, fallback = "Variable") {
  if (is.symbol(expr)) return(as.character(expr))

  out <- paste(deparse(expr, width.cutoff = 500L), collapse = "")
  if (!length(out) || is.na(out) || !nzchar(out)) fallback else out
}

.r4vn_epi_variable_label <- function(x) {
  label <- attr(x, "label", exact = TRUE)

  if (
    is.null(label) ||
      !length(label) ||
      is.na(label[1L]) ||
      !nzchar(trimws(as.character(label[1L])))
  ) {
    return(NULL)
  }

  trimws(as.character(label[1L]))
}

.r4vn_epi_variable_display <- function(x, variable_name) {
  label <- .r4vn_epi_variable_label(x)

  if (is.null(label)) return(variable_name)

  if (identical(label, variable_name)) return(variable_name)

  paste0(label, " [", variable_name, "]")
}

# ----------------------------------------------------------------------------
# Helpers for value labels in epidemiological tables
# ----------------------------------------------------------------------------

# Return a named map: original code -> display label.
#
# R4VN preserves original codes for factor variables in the `r4vn_values`
# attribute, for example c("0" = "No", "1" = "Yes"). Haven-labelled
# variables use the `labels` attribute and are converted to the same map.
.r4vn_epi_value_labels <- function(x) {
  # Return a normalized named character vector in the form:
  #   original code -> display label
  #
  # R4VN data may carry `r4vn_values` in either direction depending on the
  # source/older package version. Therefore this helper accepts both:
  #   c("0" = "No",  "1" = "Yes")
  # and
  #   c("No" = "0", "Yes" = "1")
  # and normalizes both to code -> label.
  map <- attr(x, "r4vn_values", exact = TRUE)

  if (!is.null(map) && length(map) && !is.null(names(map))) {
    nm <- as.character(names(map))
    val <- as.character(unname(map))

    valid <- !is.na(nm) & nzchar(nm) & !is.na(val)
    nm <- nm[valid]
    val <- val[valid]

    if (length(nm)) {
      # Prefer the orientation that clearly contains the usual binary codes
      # in the code position. This also handles numeric codes other than 0/1.
      name_has_binary_code <- any(nm %in% c("0", "1"))
      value_has_binary_code <- any(val %in% c("0", "1"))

      if (name_has_binary_code || !value_has_binary_code) {
        out <- val
        names(out) <- nm
      } else {
        out <- nm
        names(out) <- val
      }

      # Keep the last declaration when duplicated codes exist, but never use
      # [[name]] later because unusual/empty names can otherwise trigger
      # "subscript out of bounds" in some imported metadata.
      out <- out[!duplicated(names(out), fromLast = TRUE)]
      return(out)
    }
  }

  # haven_labelled variables use labels in the form label -> numeric code,
  # e.g. c("No" = 0, "Yes" = 1). Normalize to code -> label.
  labs <- attr(x, "labels", exact = TRUE)

  if (!is.null(labs) && length(labs) && !is.null(names(labs))) {
    nm <- as.character(names(labs))
    code <- as.character(unname(labs))
    valid <- !is.na(nm) & nzchar(nm) & !is.na(code) & nzchar(code)

    nm <- nm[valid]
    code <- code[valid]

    if (length(code)) {
      out <- nm
      names(out) <- code
      out <- out[!duplicated(names(out), fromLast = TRUE)]
      return(out)
    }
  }

  NULL
}

# Convert the currently stored values to display labels where value-label
# metadata are available. Factors already contain their display labels in
# as.character(), so they are returned unchanged.
.r4vn_epi_display_values <- function(x) {
  raw <- as.character(x)

  if (is.factor(x)) return(raw)

  map <- .r4vn_epi_value_labels(x)
  if (is.null(map)) return(raw)

  shown <- unname(map[raw])
  missing_label <- is.na(shown) & !is.na(raw)
  shown[missing_label] <- raw[missing_label]
  shown[is.na(raw)] <- NA_character_
  shown
}

# Return observed analysis levels in a stable order. For factor variables the
# displayed factor-level order is respected. Numeric values are sorted, while
# character variables retain order of first appearance.
.r4vn_epi_observed_levels <- function(x) {
  raw <- as.character(x)
  observed <- raw[!is.na(raw)]

  if (!length(observed)) return(character())

  if (is.factor(x)) {
    lev <- levels(x)
    return(lev[lev %in% observed])
  }

  if (is.logical(x)) {
    lev <- c("FALSE", "TRUE")
    return(lev[lev %in% observed])
  }

  if (is.numeric(x) || is.integer(x)) {
    z <- sort(unique(x[!is.na(x)]))
    return(as.character(z))
  }

  unique(observed)
}

# Normalize category labels for conservative recognition of common positive
# and negative binary categories.
.r4vn_epi_normalize_level <- function(x) {
  z <- trimws(tolower(as.character(x)))

  ascii <- suppressWarnings(iconv(z, from = "", to = "ASCII//TRANSLIT"))
  use_ascii <- !is.na(ascii) & nzchar(ascii)
  z[use_ascii] <- ascii[use_ascii]

  z <- gsub("[^a-z0-9]+", " ", z)
  z <- gsub("\\s+", " ", z)
  trimws(z)
}

.r4vn_epi_matches_any <- function(x, patterns) {
  vapply(
    x,
    function(z) {
      any(vapply(patterns, function(p) grepl(p, z, perl = TRUE), logical(1)))
    },
    logical(1)
  )
}

# Build a map from original value codes to the analysis levels that are
# actually present in the vector. This is the key step that lets a factor such
# as levels c("Không", "Có") retain the knowledge that "Có" originated from
# code 1 when R4VN imported or labelled the data.
.r4vn_epi_code_level_map <- function(x, raw_levels, display_levels) {
  map <- .r4vn_epi_value_labels(x)
  out <- character()

  if (!is.null(map) && length(map)) {
    for (code in names(map)) {
      # Use single-bracket lookup deliberately. Imported value-label metadata
      # can contain unusual names; map[[code]] may throw "subscript out of
      # bounds", while map[code] safely returns NA when no match exists.
      label <- unname(map[code])[1L]
      if (is.na(label) || !nzchar(label)) next

      hit <- which(raw_levels == code)

      if (!length(hit)) {
        hit <- which(raw_levels == label)
      }

      if (!length(hit)) {
        hit <- which(display_levels == label)
      }

      if (!length(hit)) {
        hit <- which(
          tolower(display_levels) == tolower(label)
        )
      }

      if (length(hit) == 1L) {
        out[code] <- raw_levels[hit]
      }
    }
  }

  # If the currently stored levels themselves are codes (for example 0/1),
  # they are valid code mappings even when no metadata attribute exists.
  for (lev in raw_levels) {
    if (!lev %in% names(out)) out[lev] <- lev
  }

  out
}

.r4vn_epi_level_from_code <- function(code, x, raw_levels, display_levels) {
  code <- as.character(code)[1L]
  map <- .r4vn_epi_code_level_map(x, raw_levels, display_levels)

  if (!is.na(code) && nzchar(code) && code %in% names(map)) {
    level <- unname(map[code])[1L]
    if (!is.na(level) && level %in% raw_levels) return(level)
  }

  NULL
}

.r4vn_epi_code_for_level <- function(level, x, raw_levels, display_levels) {
  map <- .r4vn_epi_code_level_map(x, raw_levels, display_levels)
  hit <- names(map)[unname(map) == as.character(level)[1L]]

  if (length(hit)) return(hit[1L])
  NULL
}

# Default epidemiological orientation.
#
# Priority when event/exposed is omitted:
#   1. original value code 1, when identifiable;
#   2. TRUE for logical binary variables;
#   3. common affirmative labels such as Yes/Có/Positive;
#   4. the category opposite a recognized negative label;
#   5. the second observed level as a deterministic fallback.
#
# The same rule is used for outcome and exposure: code 1 means Case for the
# outcome and Exposed for the exposure. Users can always override this with
# `event =` and `exposed =`.
.r4vn_epi_default_binary_level <- function(x, raw_levels, display_levels,
                                           role = c("event", "exposed")) {
  role <- match.arg(role)
  if (!length(raw_levels)) return(NULL)

  # First choice: original code 1. This works for numeric 0/1 variables,
  # factors created/imported by R4VN with r4vn_values metadata, and labelled
  # data that preserve original value labels.
  code_one <- .r4vn_epi_level_from_code(
    "1",
    x,
    raw_levels,
    display_levels
  )

  if (!is.null(code_one)) return(code_one)

  raw_norm <- .r4vn_epi_normalize_level(raw_levels)
  shown_norm <- .r4vn_epi_normalize_level(display_levels)

  if ("true" %in% raw_norm) {
    return(raw_levels[match("true", raw_norm)])
  }

  positive_patterns <- if (identical(role, "event")) {
    c(
      "^yes$", "^y$", "^co$", "^co\\b",
      "^case$", "^event$", "^positive$", "^duong tinh$",
      "^present$", "^true$", "^mac$", "^mac\\b",
      "^benh$", "^benh\\b"
    )
  } else {
    c(
      "^yes$", "^y$", "^co$", "^co\\b",
      "^exposed$", "^treated$", "^treatment$", "^positive$",
      "^present$", "^true$"
    )
  }

  negative_patterns <- if (identical(role, "event")) {
    c(
      "^no$", "^n$", "^khong$", "^khong\\b",
      "^non ?case$", "^nonevent$", "^negative$", "^am tinh$",
      "^absent$", "^false$"
    )
  } else {
    c(
      "^no$", "^n$", "^khong$", "^khong\\b",
      "^unexposed$", "^untreated$", "^control$", "^negative$",
      "^absent$", "^false$"
    )
  }

  pos <- which(.r4vn_epi_matches_any(shown_norm, positive_patterns))
  if (length(pos) == 1L) return(raw_levels[pos])

  neg <- which(.r4vn_epi_matches_any(shown_norm, negative_patterns))
  if (length(neg) == 1L && length(raw_levels) == 2L) {
    return(raw_levels[setdiff(seq_along(raw_levels), neg)])
  }

  if (length(raw_levels) >= 2L) return(raw_levels[2L])
  raw_levels[1L]
}

# Resolve an explicit event/exposed argument against:
#   * an original code, e.g. event = 1;
#   * the currently stored analysis level;
#   * a displayed value label, e.g. event = "Có".
.r4vn_epi_resolve_level <- function(value, x, raw_levels, display_levels, arg,
                                    default = NULL) {
  if (is.null(value)) {
    if (!is.null(default) && length(default)) {
      return(as.character(default)[1L])
    }
    return(raw_levels[2L])
  }

  candidate <- as.character(value)[1L]

  # Prefer original-code interpretation when metadata are available. Thus
  # event = 1 works even when the factor currently displays "Có".
  from_code <- .r4vn_epi_level_from_code(
    candidate,
    x,
    raw_levels,
    display_levels
  )
  if (!is.null(from_code)) return(from_code)

  if (candidate %in% raw_levels) return(candidate)

  hit <- which(display_levels == candidate)
  if (length(hit) == 1L) return(raw_levels[hit])

  hit <- which(
    tolower(display_levels) == tolower(candidate)
  )
  if (length(hit) == 1L) return(raw_levels[hit])

  available <- vapply(
    seq_along(raw_levels),
    function(i) {
      code <- .r4vn_epi_code_for_level(
        raw_levels[i],
        x,
        raw_levels,
        display_levels
      )
      if (!is.null(code) && !identical(code, raw_levels[i])) {
        paste0(display_levels[i], " [code ", code, "]")
      } else {
        display_levels[i]
      }
    },
    character(1)
  )

  stop(
    sprintf(
      "`%s` was not found in the variable. Available levels are: %s.",
      arg,
      paste(available, collapse = ", ")
    ),
    call. = FALSE
  )
}

.r4vn_epi_label_for_raw <- function(raw_value, raw_levels, display_levels) {
  hit <- match(as.character(raw_value), raw_levels)
  if (is.na(hit)) as.character(raw_value) else display_levels[hit]
}

.r4vn_epi_label_with_code <- function(label, code = NULL) {
  if (is.null(code) || !length(code) || is.na(code) || !nzchar(as.character(code))) {
    return(as.character(label))
  }

  paste0(as.character(label), " (code ", as.character(code)[1L], ")")
}

# Make the generic immediate epii() output self-explanatory when epi()
# knows the actual variable names, labels, and category labels.
.r4vn_epi_apply_metadata <- function(out, meta) {
  if (!inherits(out, "r4vn_stat")) return(out)

  if (!is.null(out$sections[["Epidemiological measures"]])) {
    measures <- out$sections[["Epidemiological measures"]]

    if (
      is.data.frame(measures) &&
        "Measure" %in% names(measures) &&
        nrow(measures) >= 2L
    ) {
      measures$Measure[1L] <- paste0(
        "Risk/prevalence: ",
        meta$exposed_label
      )

      measures$Measure[2L] <- paste0(
        "Risk/prevalence: ",
        meta$unexposed_label
      )

      out$sections[["Epidemiological measures"]] <- measures
    }
  }

  analysis_note <- paste0(
    "Outcome: ", meta$outcome_display,
    ". Case/event = ",
    .r4vn_epi_label_with_code(meta$event_label, meta$event_code),
    "; Noncase = ",
    .r4vn_epi_label_with_code(meta$nonevent_label, meta$nonevent_code),
    ". Exposure: ", meta$exposure_display,
    ". Exposed = ",
    .r4vn_epi_label_with_code(meta$exposed_label, meta$exposed_code),
    "; Unexposed = ",
    .r4vn_epi_label_with_code(meta$unexposed_label, meta$unexposed_code),
    "."
  )

  old_notes <- out$notes
  if (is.null(old_notes)) old_notes <- character()

  out$notes <- unique(c(old_notes, analysis_note))
  out$epi_meta <- meta
  out
}

# ============================================================================
# Function source: epi.R
# ============================================================================
#' @rdname epi
#' @export
epi <- function(outcome, exposure, by = NULL, data = NULL, event = NULL, exposed = NULL,
                level = 0.95, correction = 0.5, digits = 3, p_digits = 3,
                show = TRUE, console = FALSE) {
  # `event` determines the Case column (always displayed first).
  # `exposed` determines the Exposed row (always displayed first).
  # When omitted, original value code 1 is used whenever it can be
  # identified. Otherwise common affirmative coding (TRUE/Yes/Co/Positive)
  # is detected, followed by the second observed level as a fallback.
  call <- match.call()
  env <- parent.frame()
  data <- .r4vn_stat_data(data)

  outcome_expr <- substitute(outcome)
  exposure_expr <- substitute(exposure)

  y <- .r4vn_eval_var(outcome_expr, data, env, "outcome")
  x <- .r4vn_eval_var(exposure_expr, data, env, "exposure")

  outcome_name <- .r4vn_epi_expression_name(outcome_expr, "outcome")
  exposure_name <- .r4vn_epi_expression_name(exposure_expr, "exposure")

  outcome_display <- .r4vn_epi_variable_display(y, outcome_name)
  exposure_display <- .r4vn_epi_variable_display(x, exposure_name)

  y_raw <- as.character(y)
  x_raw <- as.character(x)
  y_display <- .r4vn_epi_display_values(y)
  x_display <- .r4vn_epi_display_values(x)

  y_ok <- !is.na(y_raw)
  x_ok <- !is.na(x_raw)

  # Stable observed-level order is important because the selected event and
  # exposed levels define both table orientation and effect direction.
  yl <- .r4vn_epi_observed_levels(y)
  xl <- .r4vn_epi_observed_levels(x)

  if (length(yl) != 2L || length(xl) != 2L) {
    stop(
      "Outcome and exposure must each have two observed levels.",
      call. = FALSE
    )
  }

  y_display_levels <- vapply(
    yl,
    function(z) {
      i <- which(y_ok & y_raw == z)[1L]
      y_display[i]
    },
    character(1)
  )

  x_display_levels <- vapply(
    xl,
    function(z) {
      i <- which(x_ok & x_raw == z)[1L]
      x_display[i]
    },
    character(1)
  )

  default_event <- .r4vn_epi_default_binary_level(
    y,
    yl,
    y_display_levels,
    role = "event"
  )

  default_exposed <- .r4vn_epi_default_binary_level(
    x,
    xl,
    x_display_levels,
    role = "exposed"
  )

  # Epidemiological orientation is controlled by event/exposed, not by a
  # cosmetic row/column reversal:
  #   first outcome column = Case/event
  #   second outcome column = Noncase/non-event
  #   first exposure row   = Exposed
  #   second exposure row  = Unexposed
  ev <- .r4vn_epi_resolve_level(
    event,
    y,
    yl,
    y_display_levels,
    "event",
    default = default_event
  )

  ex <- .r4vn_epi_resolve_level(
    exposed,
    x,
    xl,
    x_display_levels,
    "exposed",
    default = default_exposed
  )

  non_event <- setdiff(yl, ev)
  non_exposed <- setdiff(xl, ex)

  if (length(non_event) != 1L || length(non_exposed) != 1L) {
    stop(
      "Outcome and exposure must each contain exactly two observed levels.",
      call. = FALSE
    )
  }

  event_label <- .r4vn_epi_label_for_raw(
    ev,
    yl,
    y_display_levels
  )

  nonevent_label <- .r4vn_epi_label_for_raw(
    non_event,
    yl,
    y_display_levels
  )

  exposed_label <- .r4vn_epi_label_for_raw(
    ex,
    xl,
    x_display_levels
  )

  unexposed_label <- .r4vn_epi_label_for_raw(
    non_exposed,
    xl,
    x_display_levels
  )

  event_code <- .r4vn_epi_code_for_level(
    ev,
    y,
    yl,
    y_display_levels
  )

  nonevent_code <- .r4vn_epi_code_for_level(
    non_event,
    y,
    yl,
    y_display_levels
  )

  exposed_code <- .r4vn_epi_code_for_level(
    ex,
    x,
    xl,
    x_display_levels
  )

  unexposed_code <- .r4vn_epi_code_for_level(
    non_exposed,
    x,
    xl,
    x_display_levels
  )

  meta <- list(
    outcome_name = outcome_name,
    outcome_label = .r4vn_epi_variable_label(y),
    outcome_display = outcome_display,
    exposure_name = exposure_name,
    exposure_label = .r4vn_epi_variable_label(x),
    exposure_display = exposure_display,
    event = ev,
    event_label = event_label,
    event_code = event_code,
    nonevent = non_event,
    nonevent_label = nonevent_label,
    nonevent_code = nonevent_code,
    exposed = ex,
    exposed_label = exposed_label,
    exposed_code = exposed_code,
    unexposed = non_exposed,
    unexposed_label = unexposed_label,
    unexposed_code = unexposed_code
  )

  make <- function(idx) {
    yy <- y_raw[idx]
    xx <- x_raw[idx]

    m <- matrix(
      c(
        sum(xx == ex & yy == ev),
        sum(xx == ex & yy == non_event),
        sum(xx == non_exposed & yy == ev),
        sum(xx == non_exposed & yy == non_event)
      ),
      nrow = 2L,
      ncol = 2L,
      byrow = TRUE
    )

    # Keep both the category labels and the variable identification.
    # statistics-utils.R preserves these dimnames in .r4vn_as_2x2().
    dimnames(m) <- stats::setNames(
      list(
        c(
          paste0(exposed_label, " (Exposed)"),
          paste0(unexposed_label, " (Unexposed)")
        ),
        c(
          paste0(event_label, " (Case)"),
          paste0(nonevent_label, " (Noncase)")
        )
      ),
      c(
        exposure_display,
        outcome_display
      )
    )

    m
  }

  ok <- stats::complete.cases(y, x)

  if (missing(by) || identical(substitute(by), quote(NULL))) {
    out <- epii(
      make(ok),
      level = level,
      correction = correction,
      digits = digits,
      p_digits = p_digits,
      show = FALSE,
      console = FALSE
    )
  } else {
    by_expr <- substitute(by)
    g <- .r4vn_eval_var(by_expr, data, env, "by")
    ok <- stats::complete.cases(y, x, g)
    gf <- droplevels(factor(g[ok]))
    inds <- split(
      which(ok),
      factor(g[ok], levels = levels(gf))
    )
    tabs <- lapply(inds, make)

    by_name <- .r4vn_epi_expression_name(by_expr, "by")
    by_display <- .r4vn_epi_variable_display(g, by_name)
    meta$by_name <- by_name
    meta$by_label <- .r4vn_epi_variable_label(g)
    meta$by_display <- by_display

    out <- epii(
      by = tabs,
      level = level,
      correction = correction,
      digits = digits,
      p_digits = p_digits,
      show = FALSE,
      console = FALSE
    )
  }

  out <- .r4vn_epi_apply_metadata(out, meta)

  .r4vn_epi_display(
    out,
    call,
    show,
    console
  )
}

# ============================================================================
# Function source: cc.R
# ============================================================================
#' @rdname epi
#' @export
cc <- function(outcome, exposure, by = NULL, data = NULL, event = NULL, exposed = NULL,
               level = 0.95, correction = 0.5, digits = 3, p_digits = 3,
               show = TRUE, console = FALSE) {
  call <- match.call()
  mc <- call
  mc[[1L]] <- epi
  mc$show <- FALSE
  mc$console <- FALSE
  out <- eval(mc, envir = parent.frame())
  .r4vn_epi_display(out, call, show, console)
}

# ============================================================================
# Function source: cs.R
# ============================================================================
#' @rdname epi
#' @export
cs <- function(outcome, exposure, by = NULL, data = NULL, event = NULL, exposed = NULL,
               level = 0.95, correction = 0.5, digits = 3, p_digits = 3,
               show = TRUE, console = FALSE) {
  call <- match.call()
  mc <- call
  mc[[1L]] <- epi
  mc$show <- FALSE
  mc$console <- FALSE
  out <- eval(mc, envir = parent.frame())
  .r4vn_epi_display(out, call, show, console)
}

# ============================================================================
# Function source: mcc.R
# ============================================================================
#' @rdname mcc
#' @export
mcc <- function(case, control, data = NULL, exposed = NULL, level = 0.95, digits = 3,
                p_digits = 3, show = TRUE, console = FALSE) {
  call <- match.call()
  env <- parent.frame()
  data <- .r4vn_stat_data(data)
  ca <- .r4vn_eval_var(substitute(case), data, env, "case")
  co <- .r4vn_eval_var(substitute(control), data, env, "control")
  ok <- stats::complete.cases(ca, co)
  lev <- unique(c(as.character(ca[ok]), as.character(co[ok])))
  ex <- as.character(exposed %||% tail(lev, 1))
  a <- sum(as.character(ca[ok]) == ex & as.character(co[ok]) == ex)
  b <- sum(as.character(ca[ok]) == ex & as.character(co[ok]) != ex)
  c <- sum(as.character(ca[ok]) != ex & as.character(co[ok]) == ex)
  d <- sum(as.character(ca[ok]) != ex & as.character(co[ok]) != ex)

  out <- mcci(
    a, b, c, d,
    level = level,
    digits = digits,
    p_digits = p_digits,
    show = FALSE,
    console = FALSE
  )

  .r4vn_epi_display(out, call, show, console)
}

# ============================================================================
# Function source: ir.R
# ============================================================================
#' @rdname ir
#' @export
ir <- function(cases, exposure, time, data = NULL, exposed = NULL, digits = 4, p_digits = 3,
               level = 0.95, show = TRUE, console = FALSE) {
  call <- match.call()
  env <- parent.frame()
  data <- .r4vn_stat_data(data)
  ca <- .r4vn_eval_var(substitute(cases), data, env, "cases")
  exv <- .r4vn_eval_var(substitute(exposure), data, env, "exposure")
  tm <- .r4vn_eval_var(substitute(time), data, env, "time")
  ok <- stats::complete.cases(ca, exv, tm)
  lev <- unique(as.character(exv[ok]))
  ex <- as.character(exposed %||% tail(lev, 1))

  out <- iri(
    sum(ca[ok & as.character(exv) == ex]),
    sum(ca[ok & as.character(exv) != ex]),
    sum(tm[ok & as.character(exv) == ex]),
    sum(tm[ok & as.character(exv) != ex]),
    level = level,
    digits = digits,
    p_digits = p_digits,
    show = FALSE,
    console = FALSE
  )

  .r4vn_epi_display(out, call, show, console)
}