Fit a model of triage error, with its between-team variance

Models the called sighting category conditional on the true one, so the object being fitted is a row of the confusion matrix rather than a marginal rate. Where sdt_fit_behaviour() answers "how often does a team call red", this answers "given that this casualty is red, what does a team call, and how far do teams differ in getting it wrong".

Usage

sdt_fit_triage_error(
  data,
  truth = "triage",
  estimate = "called_triage",
  group = "group_id",
  block = "session_id",
  card = "card",
  engine = "lme4",
  verbose = TRUE,
  card_effect = c(hit = "fixed", split = "none"),
  frame = c("all", "matched"),
  min_coverage = NULL,
  joint = FALSE
)

Arguments

data: A data frame of decision episodes, one row per (group, card), carrying both the true and the called category – the shape sdt_apply_gold_standard() returns.
truth, estimate: Columns holding the true and the called sighting category. Both must hold sdt_triage_levels().
group, block, card: Columns holding the team, the session the team played in, and the patient card. block may be NULL to fit without a session term. Group labels that repeat across blocks are nested inside their block, with a warning.
engine: "lme4" (default) or "glmmTMB".
verbose: If TRUE, report progress via cli.
card_effect: How the card enters each part: one or two of "fixed", "random" and "none". See the sections above.
frame: "all" (default) or "matched" for the composition-matched frame of sdt_fit_behaviour(). The default differs from that function on purpose: there dispersion is measured off the same frame that was fitted, while here the card term is in the model and the between-team dispersion is measured separately by sdt_error_fidelity(), which does its own matching. Matching here would cost a category its cards – on this corpus the matched frame keeps 2 of the 5 yellow cards.
min_coverage: Passed to the matched frame; NULL chooses it.
joint: If TRUE, also fit the full multinomial mixed model per category as a diagnostic. Slow, and not used by the simulator.

Value

An S3 object of class sdt_triage_error_model:

cells: One row per (category, part) fitted, with n_obs, prevalence, the card effect actually used, sd_group, sd_block, sd_card, the null fit's sd_block_null and aic_null, delta_aic, chisq, df, p_value, and the singular, boundary and separated flags.
fits,fits_null: Named lists of fitted models, keyed "category.part", with and without the team term.
joint: NULL, or one row per category with the joint fit's team correlation, per-part SDs, AIC and why it was rejected.
terms: Column names, card levels, the card-to-truth lookup and the wrong-answer ordering the split refers to.
data: The frame actually fitted on.
meta: Counts, frame, engine, versions, runtime and the observed confusion matrix.
skipped: One row per cell not modelled, with the reason.

The parameterisation

For each true category k the response has three levels and therefore two free logits. They are fitted as a sequential pair:

hit – P(called = k | true = k), over every row with truth k.
split – P(called = w1 | called != k, true = k), over the error rows only, where w1 is the more severe of the two remaining categories.

Each is a binomial mixed model

y ~ card + (1 | block) + (1 | group)

This is an exact reparameterisation of the three-category multinomial, not an approximation to it: the multinomial likelihood factorises into the two binomial likelihoods, so fitting them separately fits the multinomial with the team random effects of the two logits held independent. The single parameter the decomposition gives up is the correlation between a team's accuracy and its error direction, and joint = TRUE estimates exactly that (see below).

Why the card term is load-bearing

Cards inside one true category are not alike. On the validation corpus red card 14 was called correctly by 26% of teams and 15 of the 30 red cards by every team. Removing the card term does not shrink the team effect, it erases it: the fitted team SD goes to exactly 0 for red and for yellow, against 0.544 and 0.639 with the card term, and the AIC is 473 and 137 worse. Every group saw nearly the same deck, so the spread of raw hit rates across teams is narrower than binomial noise around a single pooled probability, and a model with one probability per category can only account for that by setting the team variance to zero.

card_effect therefore defaults to "fixed" for the hit part. It takes one value for both parts, or two – c("fixed", "none") sets the hit part then the split part – or a named vector such as c(hit = "fixed", split = "random").

The split part defaults to "none", and that is measured rather than assumed. Its cells hold few error rows with a small minority (176, 124 and 529 rows with 10, 24 and 38 minority calls on this corpus); a card-fixed split runs the yellow team SD up to 41 on the logit scale, is worse by AIC than no card term for red, and leaves any card that no team ever miscalled with no coefficient at all. It also costs little, because for red both wrong answers are under-triage and for green both are over-triage: the split moves cells inside the confusion matrix without moving the over- and under-triage rates.

The ladder, and what gets skipped

A card term is attempted only where there are at least as many minority events as card levels, and a fit is rejected if it is rank deficient, fails to converge, or lands a random-effect SD above 5 on the logit scale – which is separation, not a tendency. Each rejection demotes the cell one rung down "fixed", "random", "none" and is recorded in the cell's demotions. A cell that runs out of rungs, or that is constant, or that has too few events or groups, appears in $skipped with its reason and is never dropped in silence. A fit that lands a zero team variance is kept and flagged singular: teams not differing is a finding, not a failure.

The full multinomial

joint = TRUE stacks a category's two parts and fits them together with an unstructured team covariance, which is the full multinomial mixed model for that category. It is a diagnostic: it measures the one parameter the default gives up, and sdt_simulate_triage() never reads it. Run on the validation corpus it converged for green only – correlation -0.050, and a joint AIC equal to the sum of the two separate ones – while red came back rank deficient with a team correlation of 0.99 and yellow ran its split SD to 33-40. The full multinomial does not fit reliably at this n, and where it does fit it says the constraint costs nothing.

The team term

Every cell is fitted twice: once as above and once with (1 | group) removed from exactly the same rows. $cells carries sd_group, the null fit's AIC, delta_aic and a likelihood-ratio test, so the necessity of the team term is a measured claim per cell rather than a stated one. The null fits are kept, and sdt_simulate_triage() can simulate from them with team_effects = FALSE.

Examples

set.seed(1)
cards <- sprintf("c%02d", seq_len(12))
truth <- rep(c("red", "yellow", "green"), each = 4)
d <- expand.grid(card = cards, group_id = sprintf("g%02d", seq_len(20)),
                 stringsAsFactors = FALSE)
d$session_id <- ifelse(d$group_id < "g11", "s1", "s2")
d$triage <- truth[match(d$card, cards)]
b <- stats::rnorm(20, sd = 1)
p <- stats::plogis(1 + b[match(d$group_id, unique(d$group_id))])
hit <- stats::runif(nrow(d)) < p
wrong <- vapply(seq_len(nrow(d)), function(i) {
  sample(setdiff(sdt_triage_levels(), d$triage[[i]]), 1)
}, character(1))
d$called_triage <- ifelse(hit, d$triage, wrong)
if (requireNamespace("lme4", quietly = TRUE)) {
  sdt_fit_triage_error(d, verbose = FALSE)
}
#> Warning: ! 1 cell fitted a zero team variance: "green.split".
#> ✖ Simulated teams will be identical there.
#> ℹ That is a finding about the cell, not a failure of the fit.
#> <sdt_triage_error_model> lme4 | all frame (all rows)
#> 12 cards x 20 groups in 2 blocks = 240 rows (of 240)
#> 6 cells fitted in 0.4s:
#>   cell               n    prev sd_group sd_block     card    dAIC
#>   red.hit           80   0.750    1.284    0.000    fixed    -2.0
#>   red.split         20   0.250    0.000    0.000     none     2.0
#>   yellow.hit        80   0.800    0.969    0.000    fixed     0.5
#>   yellow.split      16   0.375    0.453    0.000     none     2.0
#>   green.hit         80   0.750    1.009    0.000    fixed    -0.4
#>   green.split       20   0.600    0.000    0.000     none     2.0 singular