fix: form the reconstruction covariance on the parameters the solve solved for

[Jammy2211]

Summary

AbstractInversion.reconstruction does not always solve the full linear system. Under use_edge_zeroed_pixels it subsets curvature_reg_matrix to zeroed_ids_to_keep, solves the reduced problem and scatters the answer back with exact zeros at the excluded pixels — the mesh's poorly-constrained boundary vertices, zeroed precisely to keep the inversion stable.

reconstruction_covariance_matrix inverted the full matrix regardless. So it re-admitted into an explicit inverse the very rows the solve dropped to stay stable, and reconstruction_noise_map reported a finite noise value for a pixel whose reconstruction reads exactly 0.0 because it was never solved for.

Not opt-in, as the issue originally scoped it. Delaunay.zeroed_pixels is a count defaulting to 0, but RectangularRTUAdaptDensity.zeroed_pixels returns the whole edge ring unconditionally — 108 of 784 parameters on a 28×28, 156 of 1600 on a 40×40. And the workspace's own delaunay.py example uses zeroed_pixels=30, so the documented Delaunay idiom hits this path too.

Closes #492.

The structural fix, not just the local one

The predicate deciding whether the solve subsets the system lived inline in reconstruction and nowhere else, so the covariance had no way to know. It is now one property, solve_ids_to_keep, which both read. It returns None (not "every index") when the solve is full, so that path never enters the indexing and scatter-back code and stays byte-identical. It reproduces the deliberate nesting — use_edge_zeroed_pixels is consulted only when use_positive_only_solver is on — and the 84b9ed4 comment forbidding its "fix" is preserved and extended.

Excluded entries are NaN, not zero

Zero is a legitimate covariance value ("known exactly"), indistinguishable from a real result. These parameters were held at zero by construction and never estimated. Both consumers handle NaN and neither needed changing: norm_from derives colour limits with np.nanmax, and save_reconstruction_csv already writes nan in this column when the covariance cannot be computed. Verified by rendering the mapper subplot in both linear and log10 scales and writing the CSV.

A signal-to-noise map formed as reconstruction / reconstruction_noise_map gives 0.0 / NaN = NaN silently, and NaN >= threshold is False, so these pixels fall outside a S/N cut rather than counting as significant. (Reporting 0 instead would have made that 0/0 — the same NaN, plus a RuntimeWarning on every call.)

API Changes

No signature changes and no removals. One property added, and the values of two existing properties change whenever the solve subsets the system:

• Added AbstractInversion.solve_ids_to_keep — the single answer to "which parameters did the solve include?"
• Changed values reconstruction_covariance_matrix — formed on the solved index set, scattered back to full [total_params, total_params] shape with NaN outside it. Shape is unchanged, so callers never branch on the settings.
• Changed values reconstruction_noise_map — inherits the above; NaN at pixels the solve excluded.
• Unchanged to the bit when use_edge_zeroed_pixels=False, when use_positive_only_solver=False, when there is no mapper, or on a mesh with zeroed_pixels=0.

The value change on the pixels that ARE solved — for the release note

Inverting the submatrix is not the corresponding block of the full inverse. For symmetric positive-definite A, [A⁻¹]_keep ⪰ (A_keep)⁻¹ in the PSD ordering, so every kept pixel's variance is lower than before. Measured on real ray-traced fits (Isothermal einstein_radius=1.6 + shear, compact Sérsic source, r=3.0" mask, over_sample_size_pixelization=4, PSF and Poisson noise) at the coefficient the Bayesian evidence selects (λ*=1, interior to a 9-point grid):

mesh	zeroed	p50	p90	p99	max	>1%	>10%
RectangularBilinearAdaptDensity(28,28)	108/784	1.0002	1.069	1.44	1.81	21%	9%
Delaunay, workspace idiom, zeroed_pixels=30	30/562	1.0045	—	1.90	2.37	—	—
Delaunay, default zeroed_pixels=0	0	1.000	—	1.000	1.000	0%	0%

Downstream: 0.00% in every case. Source flux and magnification through the workspace's S/N >= 5 cut are unmoved, and lit-pixel counts are identical old vs new (7/7, 24/24, 121/121) — the pixels whose noise changes most sit at the mesh edge, where the reconstruction has no flux to move.

So: 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.

This is unrelated to the NNLS active-set caveat documented on reconstruction_noise_map (#472), where the shipped and free-set-restricted values bracket the truth. Here the excluded parameters are fixed at zero before the solver runs rather than chosen by the data, so conditioning on them is exact and there is no bracket.

The invariant, stated correctly

An earlier commit in this branch claimed reconstruction == 0.0 ⟺ noise_map is NaN. That is false, and the real-fit measurement caught it. The non-negative solver also pins pixels it did solve for at exactly 0.0. On the 28×28 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

NaN ⟹ reconstruction == 0.0 holds; the converse does not. np.isnan(reconstruction_noise_map) is the only way to identify the never-estimated pixels. Fixed in the docstring, the test, and the workspace prose.

Test Plan

• Full suite: 1177 passed, 55 skipped (pytest test_autoarray/ -x -q -n auto)
• 9 unit cases over a 4×4 mesh (4 interior parameters kept, 12 zeroed): shape, exact-NaN placement, kept block equals inv(submatrix), NaN-not-zero, the one-directional noise bound, all three full-system paths unchanged, and LinAlgError raised only when the solve actually used the non-finite entry
• 1 end-to-end case through the real solver and mesh, asserting the NaN set equals the excluded set exactly plus the one-way implication
• Consumers smoke-tested directly under a NaN-heavy noise map: CSV writer (writes nan, no RuntimeWarning), norm_from in both scales, and the mapper subplot rendered in linear and log10
• Measured on real ray-traced fits with PyAutoLens/PyAutoGalaxy/PyAutoFit at main

Full API Changes (for automation & release notes)

Added

• AbstractInversion.solve_ids_to_keep → Optional[np.ndarray] — the global parameter indices the reconstruction actually solved for, or None when it solved the full system. Read by both reconstruction and reconstruction_covariance_matrix.
• MockMapper.mesh — accessor for the mesh the mock already stored but never exposed.

Removed

• Nothing.

Changed (values only, no signatures)

• AbstractInversion.reconstruction_covariance_matrix — formed on solve_ids_to_keep and scattered back with NaN at excluded rows/columns. Shape unchanged. The finiteness guard now runs on the submatrix, so a non-finite entry in a row the solve excluded no longer raises — matching the reconstruction, which does not fail on it either.
• AbstractInversion.reconstruction_noise_map — NaN at pixels the solve excluded.
• AbstractInversion.reconstruction_noise_map_with_covariance — deprecated alias, inherits both.

Migration

• To find never-estimated pixels: np.isnan(inversion.reconstruction_noise_map).
• Not inversion.reconstruction == 0.0 — that also selects the pixels NNLS pinned, which were solved for and carry real error bars (495 vs 108 on a representative fit).
• Code that assumed the noise map was everywhere finite under use_edge_zeroed_pixels needs a np.nan* reduction or an explicit mask. Both in-repo consumers already did.

Also

• Corrected three comments describing mapper_indices as "ids of values which are on edge so zero-d and not solved for". It is not edge zeroing — it drops the linear objects carrying no regularization. Conflating the two made the reduced matrices and the zeroed solve look like one inconsistent treatment of a single index set when they are unrelated concerns.

Known, not fixed here

Pre-existing and unreachable from this change: the np.errstate(all="ignore") guard in plot/utils.py:norm_from does not suppress np.nanmax's All-NaN RuntimeWarning, which is issued through warnings rather than the floating-point error state. A per-mapper noise map cannot be all-NaN — a rectangular mesh is at least 3×3 and keeps its interior, Delaunay zeroes a count strictly below its pixel total — so this change cannot trigger it.

---

Companion workspace PR (documentation only, merges after this one): autolens_workspace feature/reconstruction-noise-map-zeroed-pixels.

Generated by the PyAutoLabs agent workflow.

---

Generated by Claude Code