Make external likelihood terms proper Gaussians in the reported NLL

[lucalavezzo]

The problem

An external likelihood term is stored expanded, exactly as its docstring says:

-log L_ext = g^T x_sub + 0.5 x_sub^T H x_sub

For a Gaussian prior N(mu, C) with H = C^-1 and g = -H mu, the real thing is 0.5 (x-mu)^T H (x-mu) + lognorm, and expanding it gives

0.5 (x-mu)^T H (x-mu)  =  g^T x + 0.5 x^T H x  +  0.5 mu^T H mu

So the stored form is short by const = 0.5 mu^T H mu, and separately by the Gaussian normalisation. It evaluates to -0.5 mu^T H mu at the prior mean instead of 0.

How this surfaced. On a card with a 2-parameter lattice prior the offset is 62.4 units. That made a penalised fit appear to have a lower NLL (4373.45) than its own unpenalised version (4433.12) — impossible, since adding a penalty can only raise the minimum. We had been correcting it by hand.

Blast radius

No fit changes. Neither scalar depends on x, so the minimum, gradient, Hessian, covariance, EDM, impacts, scans and every NLL difference taken within a single fit are untouched. That is why fits have always been correct and only reported numbers were wrong.

What was wrong: nllvalreduced, nllvalfull, and chi2_val = 2.0 * nllvalreduced in rabbit_fit.py (the absolute saturated chi2), which was off by 2 * const.

Where each scalar goes

Following _compute_lc, which writes the native constraint term in centered form cw * 0.5 * (x - x0)^2 and gates only its normalisation on full_nll:

scalar	value	included in
const	0.5 mu^T H mu = 0.5 g^T H^-1 g	reduced and full — it is part of the quadratic, not a normalisation
lognorm	0.5 (k log 2pi - log det H)	full only

lognorm in one dimension is exactly the 0.9189 + log(sigma) that _compute_lc already adds per constrained parameter. So on a card mixing native constraints with an external term, nllvalfull was previously normalising some priors and not others.

Migration: none

Both scalars are derived at load time when H is dense, so existing cards are fixed with no regeneration and no re-fitting. Verified against a real card: the 62.408046 we had been adding by hand is now produced automatically, and the term evaluates to exactly 0 at mu.

A sparse H would need a sparse solve and log-determinant, neither of which belongs in a load path — those terms keep today's behaviour with a warning and must be given the scalars at write time. To that end add_external_likelihood_term() gains mean=, which derives g = -H mu and const = 0.5 mu^T H mu from a single matvec (no solve needed when mu is known) and is the recommended way to declare a Gaussian prior. const= / lognorm= overrides exist for the raw-sparse case.

Degenerate cases

The g^T x + 0.5 x^T H x form can express things that are not priors, so these are handled explicitly rather than raising:

• grad only (H absent) — a linear tilt with no minimum; nothing to centre, nothing added
• grad zero — already centred; const = 0 exactly
• singular H — pseudo-inverse for the constrained subspace; warns if g has a component outside range(H), where the term is unbounded below and no constant centres it
• non-positive-definite H — not a density; centred but not normalised, with a warning

Tests

New tests/test_external_nll_constant.py compares two cards identical but for the term and asserts the loss difference equals the analytic Gaussian at several points, not just the minimum — the existing test_external_term.py asserts only postfit parameter values, which by construction cannot see a constant. Also covers the mean= path, sparse storage, backward compatibility when the scalars are absent, the degenerate cases, and that gradient/Hessian are unchanged.

Both the new suite and the existing test_external_term.py pass. black / isort / flake8 run with the CI flags in the CI container image are clean repo-wide.

Note for reviewers

Absolute NLL values change for any card carrying an external term. That is the fix, but it will look like a regression to anyone diffing against previous output, so it may be worth a release note.

🤖 Generated with Claude Code

https://claude.ai/code/session_01CnJ9YKK8c1q1sCDouM6y1c