Quantities of interest

QoI(name, compute; reduce=mean)

A named quantity of interest: measure something per simulation, then combine the replicates.

Three parts of ModelManager need a number out of a group of simulations — sensitivity analysis, calibration, and the post-processing sink — and each used to ask in its own shape. A QoI is that measurement written once and passed to any of them.

Arguments

• name: identifies the quantity. It is the sink's column name and the key under which CalibrationProblem reports the value to its distance. It may not contain a ., which is reserved as the separator between a quantity and its components (see below).
• compute: called with one Simulation. It may return anything reduce understands — a scalar, a vector, a Dict — except when the QoI is used as a post_processor, where reduce is never called and compute's own return value is what gets stored.

Keywords

• reduce: collapses one parameter set's replicates, mean by default. It receives the vector of everything compute returned for that set — one entry per replicate — and returns that set's value. It is not Base.reduce, and the result need not be a scalar: what it reduces is the replicate dimension, so returning a Dict of name => value is a reduction over replicates that keeps several quantities, and every consumer here understands one.

What each consumer needs back

Neither compute nor reduce is constrained by QoI itself; the requirement comes from where the QoI is used, and it does not fall on the same function in each case:

The sink is the exception, and the reason is that it fires once per simulation: there is exactly one value and nothing to combine, so reduce is never called and the freedom to return a vector does not apply. Write richer per-simulation output to the simulation's own folder instead. Returning something a consumer cannot use is that consumer's error to raise, and the sink's names the QoI and the offending type.

A Dict becomes several quantities, not one

Two consumers spread a keyed value rather than demanding one number, and both name the pieces the same way — "<qoi name>.<key>" — so one measurement names its parts identically wherever it is used:

• Sensitivity analysis runs one analysis per key, labelled "<qoi name>.<key>" — so QoI("counts", …) reducing to Dict("tumor" => …, "immune" => …) gives counts.tumor and counts.immune. Every monad must reduce to the same keys; one that does not is refused, because a sensitivity index computed over a missing value is wrong rather than approximate. A Vector is not spread by index: only its length could be checked against the other monads', and equal length is not equal meaning.
• The sink names its columns the same way: "<qoi name>.<key>". Because those names are persisted, an anonymous compute is refused outright when it spreads — its derived anon_9 would prefix every column and vary between sessions. Name the QoI, or pass a named function.

Because reduce sees every replicate's value, a measurement that needs the replicates jointly rather than as summarised numbers is expressed by having compute return the raw material — a time series, say — and letting reduce do the pooled work.

reduce is the monad-level step, not merely an average

This matters whenever the quantity involves a nonlinearity applied after the replicates are combined. Sensitivity analysis on a discrepancy-to-data score is the common case: you want each measured value averaged across replicates and then compared to data, because averaging squared errors is not the same number as squaring the averaged error. A per-simulation compute cannot do it — it has no access to the mean — but reduce can, because it receives every replicate:

observed = Dict("tumor" => 2.0, "immune" => 3.0)

QoI("mse", endpointCounts;
    reduce = per_sim -> sum((mean(getindex.(per_sim, k)) - observed[k])^2 for k in keys(observed)))

compute returns each simulation's raw values, and reduce averages per quantity, takes the squared differences, and sums them into the one number sensitivity analysis needs. Reporting spread alongside the mean is a second QoI over the same compute, with reduce comparing std to the observed spread instead.

Reading a value the sink stored earlier

Post-processing runs while a simulation's output folder still exists; post-simulation cleanup may then delete what the quantity was computed from. A QoI whose compute reads the sink instead of the output folder therefore still works afterwards, and needs nothing new — the sink is keyed by simulation ID and by the QoI's own name:

tumor = QoI("tumor", s -> finalPopulationCount(s)["tumor"])
run(trial; post_processor=tumor)                      # stores it while the output exists

stored = QoI("tumor", s -> postProcessingTable([s.id]).tumor[1])
run(MOAT(), spec; functions=[stored])                 # reads it back, output folder or not

stored automates the lookup, and defaults to :never:

QoI("tumor", computeFromOutput; stored=:prefer)    # stored value if present, else compute
QoI("tumor", computeFromOutput; stored=:require)   # stored value or an error

It defaults off because nothing records which compute produced a stored value, and no fingerprint can: redefining a function's body in place leaves both hash and nameof unchanged, so a changed compute is undetectable, while two textually identical anonymous functions hash differently, so an unchanged one is equally unrecognisable. The sink stores no provenance either.

What is robust is recomputation. verifyStoredValues recomputes where a simulation's output survives and reports agreements, mismatches, and how many could not be checked because the output is gone — which is precisely the case stored exists for. Run it before trusting stored values for anything that matters.

Examples

tumor = QoI("tumor", s -> finalPopulationCount(s)["tumor"])

# A different way of combining replicates
QoI("tumor_median", s -> finalPopulationCount(s)["tumor"]; reduce=median)

# Spread across replicates, rather than their centre
QoI("spread", s -> finalPopulationCount(s)["tumor"]; reduce=std)

# Pooling the replicates instead of reducing per-simulation numbers
QoI("slope", timeSeries; reduce=series -> fitSlope(reduce(vcat, series)))

# One QoI, three consumers
run(MOAT(), spec; functions=[tumor])
CalibrationProblem(spec, observed, tumor, mseDistance)
run(trial; post_processor=tumor)

# Several quantities from one measurement, named the same way in both places: GSA spreads
# `reduce`'s keys into `counts.tumor` / `counts.immune`, and the sink spreads `compute`'s keys
# into columns `counts.tumor` / `counts.immune`.
counts = QoI("counts", finalPopulationCount;
             reduce = per_sim -> Dict(k => mean(getindex.(per_sim, k)) for k in ("tumor", "immune")))

verifyStoredValues(q::QoI, T; rtol=1e-8, limit=nothing) → NamedTuple

Check a QoI's stored values against freshly computed ones, for the simulations of T.

Returns (; n_checked, n_agreed, n_mismatched, n_unverifiable, n_missing, mismatches). A simulation is unverifiable when its output folder is gone, which is exactly the situation stored exists for — the value may be perfectly good, but nothing here can confirm it.

Use this before trusting stored=:prefer or stored=:require on results you care about. Nothing about a stored value records which compute produced it, so recomputation is the only real check.

Example

report = verifyStoredValues(tumor, my_sampling)
# `n_mismatched == 0` alone is NOT a pass: it is also what you get when every simulation was
# skipped. Require that something was actually compared.
report.n_agreed > 0 || error("nothing was verified: $(report.n_missing) had no stored value " *
                             "and $(report.n_unverifiable) had no output folder to recompute from")
report.n_mismatched == 0 || error("stored values disagree with a fresh computation")