docs: record that the reconstruction noise map is the unconstrained uncertainty - #472
Conversation
…ncertainty
reconstruction_covariance_matrix is [F + reg_coeff*H]^-1, the posterior
covariance of the positive-negative (unconstrained) Warren & Dye solve. But
Settings.use_positive_only_solver defaults to True, so the reconstruction is
normally NNLS. Constraining s >= 0 truncates the posterior, so this overstates
per-pixel uncertainty -- more so near the s = 0 boundary, and not meaningfully
at all for pixels pinned at exactly zero.
Documentation only; no behaviour change. The maths is left alone deliberately,
because measurement showed the practical effect is small at the operating point
and because the obvious "fix" is wrong.
Measured on real ray-traced fits (Isothermal + shear,
RectangularBilinearAdaptDensity, Constant regularization), varying the
regularization coefficient and locating the Bayesian-evidence optimum that a
model-fit actually converges on:
at the evidence-optimal coefficient (lambda* = 10 in every case tested)
r_eff 0.05: 96.6% of mesh pinned, noise overstated x1.263 median, flux 0.0%
r_eff 0.10: 87.1% pinned, x1.055 median, flux -1.3%
r_eff 0.30: 42.3% pinned, x1.007 median, flux 0.0%
well below the evidence optimum (under-regularized): up to x2.8 median and
x10 on individual pixels, with source flux through a S/N >= 5 cut moving by
as much as -49%.
The Bayesian evidence selects away from the regime where this matters, so for
most fits the bias is small and one-directional (conservative).
Also recorded: restricting the covariance to the solver's free set is NOT the
correction. That treats the active set as known and so understates. The two
bracket the true truncated-Gaussian posterior, and swapping one bound for the
other would not be more correct.
Caveat on the measurement, not repeated in the docstring: the lens mass was
fixed at truth, so a real fit with a free and imperfect mass model may need a
lower coefficient to absorb residuals -- which is the regime where the gap
opens. Tracked in PyAutoMind
draft/bug/autoarray/reconstruction_noise_map_solver_mismatch.md.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_0133X4XhMV91SFjzV2mK4Ejh
Correction — one of the "defects left open" in this PR's description is not a defectThis PR's body closes by naming two things left open, one of which was:
That is intended behaviour, confirmed by the author. The control-flow description was accurate — The same claim appears in the "Out of scope" section of #468. It is wrong there too. Corrected in The other item named there stands unchanged: the covariance ignores Only residue worth anything, and it is cosmetic: the coupling is undocumented. Neither No code change; this PR's diff is unaffected. Generated by Claude Code |
…ange
Measuring this change on a real ray-traced fit (PyAutoLens + PyAutoGalaxy +
PyAutoFit at main, against this branch) turned up an overclaim in the previous
commit and supplied the magnitude it said was missing.
The invariant was stated too strongly
-------------------------------------
The docstring and the end-to-end test claimed a biconditional: a pixel reads
exactly 0.0 in the reconstruction exactly where the noise map is NaN. That is
false, and the difference is not marginal. The non-negative solver ALSO pins
pixels it did solve for at exactly 0.0, wherever the fit wants no flux there,
and those are fitted values carrying real error bars.
On a 28x28 RectangularBilinearAdaptDensity fit, 784 parameters:
reconstruction == 0.0 : 603
NaN in noise map : 108 <- only these were never estimated
solved, but pinned at 0 by NNLS : 495
So `NaN => reconstruction == 0.0` holds, the converse does not, and
`np.isnan(reconstruction_noise_map)` is the only way to identify the
never-estimated pixels. The test passed only because the 3x3 fixture it runs on
solves a single pixel, which makes the two sets coincide by accident. It now
asserts the NaN set equals the excluded set exactly, and the one-way
implication, with a comment saying why the converse is deliberately absent.
This is also precisely the Defect 1 / Defect 2 boundary: 495 pixels pinned by
the constraint are the NNLS active-set question (documented in #472, deferred),
and 108 structurally zeroed ones are this change.
Measured magnitude, at the evidence-optimal coefficient
-------------------------------------------------------
The previous commit said the magnitude had not been measured on a real fit and
that the structural proxy standing in for it should not be trusted. It should
not have been: the proxy predicted a uniform 1-2% shift and the real answer is
a much smaller median with a long tail.
Rectangular 28x28 (108 of 784 zeroed), lambda* = 1 by log evidence on a
9-point grid, interior to the grid in every case:
old/new p50 1.0002 p90 1.069 p99 1.44 max 1.81
21% of solved pixels move >1%, 9% move >10%
Delaunay, set up as the workspace's own delaunay.py does (Overlay image mesh
plus a ring of edge points, zeroed_pixels=30) -- which answers the "Delaunay
untested" caveat this change inherited:
old/new p50 1.0045 p99 1.90 max 2.37
With the default zeroed_pixels=0 the property is unchanged to the bit, and
solve_ids_to_keep is None, confirming the no-op path.
Downstream source flux and magnification through the workspace's S/N >= 5 cut
moved by 0.00% in every case. The pixels whose noise changes most sit at the
mesh edge, where the reconstruction has no flux to move: the counts of lit
pixels surviving the cut are identical old and new (7/7, 24/24, 121/121).
So the release note should say: a systematic, one-directional correction
concentrated at the mesh edge, negligible in the median, up to ~2x on
individual edge pixels, with no effect on published source flux or
magnification.
Full suite: 1177 passed, 55 skipped.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Lkq5ww6eLEvJgPFGMgMU1C
Documentation only — no behaviour change, one docstring, +28 lines.
Follow-up to #469. That PR fixed how the covariance is computed; this one records what it means.
The gap
reconstruction_covariance_matrixis[F + reg_coeff*H]^-1— the posterior covariance of the positive-negative (unconstrained) Warren & Dye (2003) eq. 12 solve. ButSettings.use_positive_only_solverdefaults toTrue, so the reconstruction is normally NNLS. Constrainings >= 0truncates the posterior, and a truncated Gaussian's covariance is not the untruncated one.So
reconstruction_noise_mapoverstates per-pixel uncertainty — more near thes = 0boundary, and not meaningfully at all for pixels the solver pinned at exactly zero.Why documentation and not a code fix
Two measured reasons.
1. The Bayesian evidence selects away from the regime where it matters. The regularization coefficient is a free parameter (
LogUniform(1e-6, 1e6)), and a pixelized fit maximises the Bayesian evidence when choosing it. Locating that optimum on real ray-traced fits (Isothermal+ shear,RectangularBilinearAdaptDensity,Constant):λ* = 10 in every case, well inside the scanned grid. Well below the evidence optimum — an under-regularized fit — the factor grows to ~2.8 median, ~10× on individual pixels, and source flux through an
S/N >= 5cut moves by as much as −49%. But that is a regime the evidence penalises.2. The obvious fix is wrong. Restricting the covariance to the solver's free set treats the active set as known, so it understates. The two quantities bracket the true truncated-Gaussian posterior — swapping one bound for the other would not be more correct. A real fix means computing the truncated posterior, which is not worth it for a bias this size.
So the honest move is to say what the number is, and let anyone quoting per-pixel error bars on a very compact source know it is good to a few tens of percent rather than exact.
Caveat on the measurement
Deliberately kept out of the docstring but recorded here and in the tracking prompt: the lens mass was fixed at truth. A real fit has it free, and a poor mass model may need a lower coefficient to absorb residuals — exactly the regime where the gap opens. Also: Nautilus samples a posterior over λ, so some mass sits below λ*. And
RectangularBilinearAdaptDensityonly; Delaunay untested.Tracked in
PyAutoMind/draft/bug/autoarray/reconstruction_noise_map_solver_mismatch.md, which also holds two related defects left open here: the noise map ignoringzeroed_ids_to_keepunderuse_edge_zeroed_pixels, and that setting being silently ignored when the positive-only solver is off.Tests
test_autoarray/inversion/— 346 passed. No code changed.🤖 Generated with Claude Code
https://claude.ai/code/session_0133X4XhMV91SFjzV2mK4Ejh
Generated by Claude Code