Add per-model bias c update for multi-model AMICA (#27)

[neuromechanist]

Closes the concrete, fixable slice of #27: the omitted per-model bias c update.

What changed

Ported Fortran's update_c (amica17.f90:1423-1429 / 1899-1901) into both AMICATorchNG and the NumPy oracle pyAMICA.py:

• c[i,h] = sum_t v_h(t)*x(i,t) / sum_t v_h(t) — the per-model, per-channel responsibility-weighted data-space mean.
• The E-step now centers each model's data before unmixing (b = W(x - c)); transform() in both backends does the same.
• Replaces the old dc = sum(g) accumulator, which was accumulated but never applied (c was frozen at 0). Renamed the update-dict key dc -> dc_numer.
• Guarded to a no-op for n_models=1: with v == 1 the update collapses to the (zero) mean of mean-removed data, and skipping it keeps single-model parity bit-exact (issue Verify NG init-basin vs Newton bit-parity via Fortran load_* init-matching #24).

Why the gap persists

A controlled 2-model A/B against the Fortran binary (same config/seed, only c toggled):

run	final LL	mean cross-corr
NG, c OFF (pre)	-3.3754	0.631
NG, c ON (fix)	-3.3762	0.642
Fortran	-3.3596	1.000

The c omission was a genuine but minor contributor (+0.011 cross-corr, LL unchanged). The dominant residual gap is intrinsic partition ambiguity (mixture-of-ICA has many near-degenerate partitions; NG is self-consistent, cross-corr 1.0 across block sizes), so >0.95 is not reachable via c alone. Details: .context/issue-27/multimodel_c_update.md.

Tests

• New: single-model c stays exactly zero after fit (issue Verify NG init-basin vs Newton bit-parity via Fortran load_* init-matching #24 regression guard); multi-model c equals the responsibility-weighted data mean and the two models center differently.
• Existing NG<->NumPy sufficient-stat parity holds (key renamed dc -> dc_numer).
• Full suite: 48 passed, 3 xfailed (pre-existing), including the single-model Fortran-parity end-to-end test. Real sample EEG only (no synthetic/mock).

[neuromechanist]

PR review (pr-review-toolkit, reviewers on Sonnet)

Ran four reviewers on the diff: code-reviewer, pr-test-analyzer, silent-failure-hunter, comment-analyzer.

Addressed

• Silent failure (important): the new c = dc_numer/dgm division could be 0/0 = NaN for a dead model (dgm[h]==0). Unlike log(gm[h])=-inf (which softmax tolerates), a NaN c poisons the next iteration's cross-model softmax for every model. Added a containment guard in both backends that keeps a zero-responsibility model's prior c, mirroring the existing mu/beta/rho non-finite guards. This also resolves the reviewer's NumPy restart-preserves-NaN-c concern (with the guard, c can no longer go NaN from a dead model).
• Test coverage: added the multi-model coverage the change opened up but the first commit left unexercised: NumPy backend c update on real data, NG↔NumPy finalized-c parity, transform() with nonzero c (verified against W(x-c) by hand), the dead-model containment guard, multi-model dc_numer blocking invariance, and do_reject + multi-model c finiteness. Strengthened the Newton multi-model test to assert finite c and full iteration count.
• Comments: fixed two pre-existing docstrings this change made stale (transform()/get_weights() said X^T @ W, now (X-c)^T @ W); corrected the Fortran citation for wc = W@c (:2178, get_unmixing_matrices); clarified update_c is a flag, not a routine; AGENTS.md "LL unchanged" → "LL comparable".

Not changed (with reason)

• code-reviewer: no issues at confidence ≥80; confirmed formula/backend parity, broadcasting, the n_models==1 guard, do_reject consistency, and persistence.
• AMICA.fit() doesn't check stop_reason before is_fitted_=True (silent-failure, MEDIUM): pre-existing wrapper behavior, not introduced here, and broader than Multi-model NG partition matching vs Fortran (>0.95 cross-corr) #27 (affects every degenerate fit). The containment guard makes the c-specific NaN path benign, lowering the practical urgency. Deferring to a separate issue rather than expanding this PR's scope.
• NumPy _reject_outliers may not filter self.data (code-reviewer FYI): pre-existing, affects all accumulators equally, unrelated to c. Out of scope.

Fast suites green: 48 passed (32 NG + 16 NumPy). Single-model Fortran parity unaffected (guard is n_models>1).

[neuromechanist]

Multi-model validation added (the confirmation that makes this mergeable)

Since mixture-of-ICA is not partition-identifiable, exact partition parity with Fortran is the wrong acceptance bar (the >0.95 in #27's title asks the algorithm to be more identifiable than it is). The well-posed test is whether the two implementations sample the same distribution over solutions.

Ran N=20 fits per implementation on the real sample EEG (n_models=2, 100 iters) and compared three distributions of the stacked 2×32 Hungarian cross-corr:

distribution	mean	sd
within-Fortran	0.634	0.042
within-NG	0.644	0.046
between (NG↔Fortran)	0.638	0.047

• Mann-Whitney (H1: between < within-Fortran): p = 0.97 — no evidence NG-vs-Fortran agreement is worse than Fortran-vs-itself.
• TOST equivalence within ±0.05: p ≈ 1e-32 → EQUIVALENT.

So NG's multi-model ensemble is statistically equivalent to Fortran's. The single-run ~0.64 cross-corr was never a defect — Fortran agrees with itself at 0.63.

Full method / figure / acceptance criteria: .context/issue-27/multimodel_distributional_equivalence.md. Reproduction harness committed alongside.

One residual, tracked as #51 (not a blocker): NG's LL distribution is ~0.02 lower and more variable than Fortran's — optimizer quality, not a correctness bug (one M-step is bit-exact vs Fortran).

This is the validation #27 needed. Ready to merge.