@@ -305,6 +305,33 @@ not exclude samples from this analysis.
   multiplicative (CZM) replacement; `replacement_delta` sets the replacement
   fraction and optional `bayesian_alpha` overrides GBM's concentration. The
   offset is not used. See [zero-handling notes](docs/statistics/zero-handling.md).
+- `ancom_bc`: ANCOM-BC sampling-fraction bias-corrected log-linear model
+  (Lin & Peddada 2020, *Nat. Commun.* 11:3514). Estimates within-group sample
+  offsets and between-group sampling-fraction bias via an E-M Gaussian mixture
+  on `log(counts + 1)`, propagating bias-estimation variance into two-sided Wald
+  tests. The effect is a bias-corrected log2 fold change; useful when asymmetric
+  shifts or load changes would bias standard CLR centering (needs at least 3
+  tested taxa; assumes fewer than half are differential).
+- `aldex2`: Dirichlet(`count + 0.5`) Monte-Carlo CLR model (Fernandes et al.
+  2013, 2014). Draws `mc_samples` (default 128, seeded by the study's `seed`)
+  instances per sample and reports standardized CLR effect sizes with
+  consistent-sign posterior predictive Welch p-values; suited to sparse or
+  variable-depth cohorts where deterministic zero imputation is undesirable.
+- `songbird`: compositional baseline-category multinomial regression with
+  reference-frame invariant CLR-centred differentials summing to zero across
+  fitted taxa (Morton et al. 2019, *Nat. Commun.* 10:2719; Zhang et al. 2017).
+  Conditions on read totals without extra-multinomial overdispersion; taxa
+  separated in one group are marked `failed` rather than diverging.
+- `multinomial`: one multinomial logistic model for the whole composition
+  (Songbird's model in pure Julia), which needs no size factor and no zero
+  replacement. The effect is a log2 fold change in relative abundance against
+  `mn_reference_taxon`; `mn_penalty` adds an L2 ridge penalty, and a taxon with
+  no reads in one group is refused rather than reported.
+- `dirichlet_multinomial`: the same joint model plus a Dirichlet prior on the
+  composition (La Rosa et al. 2012), which absorbs between-sample
+  overdispersion; `dm_overdispersion_method` picks `mle` (profiled) or `mom`.
+  Both multinomial models bound their Newton fit with `max_iter`, `tolerance`
+  and `max_seconds`, and a fit that does not converge is refused.
 
 ZINB and hurdle models can overfit sparse data; do not select ZINB solely
 because counts are sparse. Compare against NB diagnostics. The methods report
@@ -316,6 +343,14 @@ histories and is not selectable in the standard analysis configuration or UI.
 It currently has no detection covariates. Do not use it as a substitute for
 the full workflow tracked in [issue #205](https://github.com/hyperpolymath/MetaManifold-WebUI/issues/205).
 See [occupancy model notes](docs/occupancy-models.md).
+
+The multinomial models are the count-based, compositionally coherent
+alternative to `nb_glm` (which fits each taxon independently and so ignores the
+constraint the composition imposes) and to `clr_welch` (which is Gaussian on
+log-ratios and ignores the count likelihood). They are fitted in Julia without R
+and without a pseudocount; the dense observed information limits a fit to 1000
+taxa, above which raise `min_prevalence` or analyse a coarser rank.
+
 Across the differential methods, `min_prevalence` (default 0) filters taxa by the fraction of
 samples with nonzero counts. Benjamini–Hochberg adjustment includes only taxa
 that produced a p-value. Failed or filtered taxa have no substitute p-value;
@@ -326,6 +361,119 @@ Implementation and reference tests: `src/analysis/differential.jl`,
 `src/server/routes/analysis.jl`, `bench/compositional_differential/benchmark.jl`,
 and `test/unit/test_differential.jl`.
 
+#### Why a count model rather than a t-test
+
+Sequencing measures a sample of the library, not an absolute quantity: a
+sample of twice the depth has twice every count, and one taxon can only increase
+at another's expense. Comparisons of raw or rarefied abundances therefore carry
+the compositionality problem (Gloor et al. 2017). The multinomial models fit the
+whole composition at once and condition on each sample's total, so each
+coefficient is the log ratio of two taxa. That ratio is unchanged when every
+count in a sample is multiplied by a constant, which is why no size factor and no
+pseudocount is needed. The intercept is the log-odds of the taxon against the
+reference at the reference group; the group coefficient is the change in that
+log-odds, i.e. a log fold change in relative abundance. Because a log ratio is
+the sum of log ratios through any third taxon, every reported effect is coherent
+under compositionality and can be re-expressed against any other taxon.
+
+#### Choosing the reference taxon
+
+Every reported effect is anchored on the reference taxon, so the choice matters
+for interpretability and for stability, but not for the fit: the likelihood and
+the log ratio between two non-reference taxa are the same whichever taxon is the
+reference. An empty `mn_reference_taxon` takes the most abundant taxon that
+passes the prevalence filter. An explicit one must be a taxon at the analysed
+rank with reads in at least half the samples, and reads in both groups; a taxon
+seen in fewer samples makes every ratio it anchors unstable, and a taxon absent
+from one group has no finite log fold change, so both are refused with the
+reason. Taxa with no reads in one group are removed from the composition and
+listed as `failed` with their note, rather than fitted: in a joint model a single
+diverging coefficient would otherwise pull the whole fit with it.
+
+#### Overdispersion
+
+The multinomial treats every read as independent, so its variance is the sampling
+variance and nothing else. Real samples vary more than that: the same community,
+sequenced twice, differs by more than resampling reads from one library. The
+Dirichlet-multinomial adds that between-sample variation as a concentration `c`,
+with `theta = 1/(1+c)` the model's intra-sample correlation; the count variance
+becomes `N p (1 - p) (1 + (N - 1) theta)`. When the data are overdispersed the
+multinomial's standard errors are too small, and the Dirichlet-multinomial
+corrects them; when they are not, `c` runs to its upper bound, the fit collapses
+to the multinomial exactly, and the affected rows are marked `boundary` with a
+note saying so. The two models are otherwise the same model on the same scale,
+so their coefficients are directly comparable.
+
+#### Which model to use
+
+| Data and question | Method |
+| --- | --- |
+| Counts, two groups, few samples, expected overdispersion | `nb_glm` (per-taxon negative binomial; the default) |
+| Counts, many taxa and samples, a joint composition and a single reference | `multinomial` (Songbird's model) |
+| As above, with visible between-sample overdispersion | `dirichlet_multinomial` |
+| Composition already transformed, or a quick robust rank-based check | `clr_welch` |
+| Excess zeros beyond negative-binomial overdispersion, with evidence for them | `zinb` or `hurdle_nb` (two-part; see the caution above) |
+| Asymmetric shifts or a community-wide load change | `ancom_bc` (bias-corrected) |
+| Sparse or variable-depth cohorts, Monte-Carlo CLR | `aldex2` |
+| Many taxa, a joint multinomial fit, but through R and with CLR-centred effects | `songbird` |
+
+`nb_glm` fits a negative-binomial GLM per taxon with a library-size offset: it
+tolerates overdispersion through its own dispersion parameter and is the safer
+choice with very few samples per group, but each taxon is fitted independently
+and the offsets are estimated, not modelled. The multinomial models borrow
+strength across taxa by fitting one model, and their effect has a direct
+compositional reading, at the cost of a joint fit that a single separated taxon
+can compromise (hence the separation check and the optional ridge penalty), and
+of a dense information matrix that caps a fit at 1000 taxa. Songbird's L1
+sparsity path is not implemented: a lasso estimate has no valid Wald standard
+error, so this implements the L2-penalised model and Wald inference instead.
+
+The multinomial models are fitted in Julia and need no R package.
+`max_iter`, `tolerance` and `max_seconds` bound one fit; a fit that runs out of
+time fails with an explicit error instead of holding the request open. Runtime
+and memory at 100, 1000 and 10000 taxa are measured by
+`bench/compositional_scaling/scaling.jl`, which writes `results/scaling.csv`.
+
+Implementation and reference tests: `src/analysis/compositional.jl` (the two
+models), `src/analysis/differential.jl`, `src/server/routes/analysis.jl`, and
+`test/unit/test_differential.jl`. The reference values in the tests come from an
+independent implementation and are quoted to 12 digits; the special functions
+are cross-checked against mpmath at 40 digits.
+
+### Evidence Mode
+
+The **Evidence** tab turns a run's outputs into a citable epistemic record. A
+claim is a sentence a reader could disagree with, a **scope** naming the part of
+the data it is about (its *fibre*: the rows whose value at a rank, or whose
+category in a set, is one of the named labels), and a **warrant** saying why it
+follows. Claims are stored in `projects/{study}/[<group>/]<run>/evidence/claims.yml`
+and can be downloaded as that exact YAML.
+
+Three views share one record:
+
+- **Claims** — the epistemic editor (avec fibre). A claim cannot be saved
+  without a scope and at least one value, and cannot be marked `supported`
+  without a measurement, citation or inference warrant: a definition alone
+  supports nothing.
+- **Fibre & echo** — the fibre's share of every sample, its largest members, and
+  the same fibre seen through a second projection: the other classifier at the
+  same rank, and the parent rank. Agreement is reported with a verdict;
+  disagreement is a limit on how strongly the claim may be stated, and an echo
+  that cannot be computed says so rather than vanishing.
+- **Residuals** — where the fibre departs from a null you declare, either
+  `proportional` (the fibre distributed across samples as the library is) or a
+  stated `claim_share`, with Pearson residuals, a log2 ratio and a χ² statistic
+  printed without a p-value. The worst departure is shown first.
+
+The tab's Provenance card states whether `attestation.yml` is absent, degraded
+or internally inconsistent before any claim is read. Longer notes on the
+definitions, the schema and the limits are kept in `docs/evidence-mode.md`
+(documentation under `docs/` is kept local by policy — see `.gitignore`).
+
+Implementation: `src/core/evidence.jl` (schema, residuals, verdicts),
+`src/server/routes/evidence.jl` (the SQL), `frontend/src/evidence/`
+(editor, fibre visualiser, residual explorer) and `test/unit/test_evidence.jl`.
+
 ### Figures and report
 
 The Figures tab lays analysis charts out on a page (A4, letter or a journal
@@ -662,13 +810,23 @@ analysis:
     normalisation: hellinger     # hellinger | none | rarefy | srs
     depth: 0                     # rarefy/srs reads per sample; 0 = smallest non-empty sample
   differential:
-    method: nb_glm               # nb_glm | clr_welch | ancom_bc | aldex2 | songbird
+    method: nb_glm               # nb_glm | clr_welch | zinb | hurdle_nb | ancom_bc
+                                 #   | aldex2 | songbird | multinomial | dirichlet_multinomial
+    offset: tss                  # count models: tss | rle
     dispersion: ml               # nb_glm: ml | glmGamPoi
+    min_prevalence: 0.0          # fraction of samples a taxon must be seen in
     replacement_method: bayesian_multiplicative
                                  # clr_welch: bayesian_multiplicative | multiplicative_replacement
     replacement_delta: 0.65      # clr_welch: strictly between 0 and 1
     bayesian_alpha: ~            # clr_welch: null (GBM estimate) or a number > 0
     mc_samples: 128              # aldex2: Dirichlet Monte-Carlo instances in [16, 4096]
+    mn_reference_taxon: ~        # multinomial models: null picks the most abundant taxon
+    mn_penalty: 0.0              # multinomial models: L2 ridge penalty, >= 0
+    dm_overdispersion_method: mle
+                                 # dirichlet_multinomial: mle | mom
+    max_iter: 500                # multinomial models: Newton iteration cap
+    tolerance: 1.0e-8            # multinomial models: convergence tolerance
+    max_seconds: 300.0           # multinomial models: wall-clock budget
   nmds:
     max_stress: 0.2              # warn if NMDS stress exceeds this value
 ```