Changelog

0.9.5 (2026-05-29)

New modules — simulation pipeline and GLASS mocks

  • sys_mapping.simulation — end-to-end validation framework for systematic contamination recovery. Provides inject_systematics() (per-galaxy contamination weights from the forward model, Eq. 11-13 of Berlfein et al. 2024), ContaminationConfig (9-point grid of level × scenario), FITS I/O helpers, and run_wtheta_recovery() (runs all decontamination methods and compares recovered \(w(\theta)\) to truth).

  • sys_mapping.glass_mocks — full-sky lognormal mock generation via GLASS (Tessore et al. 2023): generate_glass_fullsky_mock() and measure_nz(). Advantage over the existing mocks.make_mock_catalog: proper full-sky catalogs with physically realistic clustering, no galactic-cut artefacts.

Bug fixes (five pipeline bugs, all in regression / simulation / scripts)

  • Bug A — w(θ) correction applied to wrong sample. run_wtheta_recovery was multiplying the correction weight by the clean galaxy catalog instead of the contaminated one. The effective weight is now correctly weight_cont × weight_recovered. Before this fix, the recovered \(w(\theta)\) was anti-contaminated — biased in the opposite direction from truth — producing residuals of order 10 000× the contamination bias.

  • Bug B — bias metric blow-up at zero-crossings. The fractional-bias metric used a per-bin denominator |w_true(θ)| + 1e-6. For shot-noise-dominated mocks (e.g. GLASS full-sky) where \(w_{\rm true} \sim 10^{-4}\), this produced artifacts of order \(10^6\). Fixed to scalar max_θ|w_true(θ)| throughout run_simulation_tests.py, plot_simulation_tests.py, and run_paper_validation.py.

  • Bug C — unclipped OLS / MCMC correction weights. _compute_weights returned values up to \(10^6\) for pixels where \(1 + \hat{a}\cdot t \to 0\). OLS and MCMC methods now clip to \([1/20,\, 20]\), matching the existing ISD and ElasticNet behaviour. Same fix applied in method_comparison(), run_mock_analysis_real_templates.py, and run_mock_analysis_diagnostic.py.

  • Bug D — MCMC-comb correction weight was a flawed approximation. The combined-model weight used \(1/(1+(a+b)\cdot t)\), which overcorrects by ~2× when the inferred \(\hat{b}\neq 0\) in a purely additive scenario (where \(b_{\rm true}=0\)). Replaced with the exact pixel-level inverse of the contamination model: \(w(p) = (1+\hat\delta_g^{\rm clean}(p))/(1+\delta_g^{\rm obs}(p))\) computed via invert_contamination(). This cancels weight_cont exactly when parameters equal their true values.

  • Bug E — MCMC walker initialisation too narrow; chain too short. Walkers were initialised at \(\mathcal{N}(0,\,10^{-3})\) — 20–100× narrower than the typical contamination amplitude (0.02–0.10). Changed to \(\mathcal{N}(0,\,0.05)\). Default n_steps raised from 600 to 1200; n_burn from 100 to 200.

Test suite

  • 28 real-template integration tests (tests/test_real_templates.py) now run instead of skipping: FITS files reside in ~/data/legacysurvey/dr10/systematics/0032/, not the directory root.

  • ElasticNet timing budget raised from 120 s to 600 s to accommodate CPU-only machines under load.

  • New tests: tests/test_simulation.py (19 tests) and tests/test_glass_mocks.py (12 tests).

Deprecations / cleanup

  • Five scripts moved to scripts/archive/: build_systematic_maps.py, generate_recommendations.py, generate_results_ls10_summary.py, generate_sample_pages.py, patch_params_from_partials.py. These were development utilities referenced only in the changelog.

Documentation

  • New results page results_simulation_tests.rst covering the full 3 × 3 contamination grid on GLASS and Uchuu mocks at NSIDE 32 and 64.

  • New API pages api/simulation.rst and api/glass_mocks.rst.

  • results_real_template_validation.rst updated to reflect NSIDE = 64 footprint (22 641 pixels) and correct FITS file paths.

  • All five documentation result pages regenerated with fresh figures from the fixed pipeline: mock analysis, systematic tests, validation, progressive contamination, real-template validation.

Package metadata

  • Version bumped to 0.9.5.

  • glass>=2026.1 and corrfunc added as optional extras in pyproject.toml (pip install sys_mapping[glass], pip install sys_mapping[corrfunc]).

0.2.4 (2026-05-11)

ISD-3 speed and accuracy improvements:

  • Speed — normal equations replace per-iteration SVD. iterative_systematics_decontamination now precomputes Xt = X.T once outside the loop. Each iteration forms the weighted normal equations \((\mathbf{X}^\top \mathbf{W} \mathbf{X} + \mathbf{\Lambda})\hat{\boldsymbol\alpha} = \mathbf{X}^\top \mathbf{W} \hat{\boldsymbol\delta}_g\) via a single broadcast-scale + BLAS-3 matmul, then calls lstsq on the resulting \((N_{\rm exp} \times N_{\rm exp})\) system (e.g. 55×55 for \(n_s=5, d=3\)) rather than on the full \((N_{\rm pix} \times N_{\rm exp})\) matrix. Measured speedup: 5–15× per iteration for \(N_{\rm pix} \gtrsim 10^4\).

  • Accuracy — ridge penalty on polynomial-only columns. New lambda_poly: float = 0.0 parameter added to iterative_systematics_decontamination(). When positive it adds a diagonal ridge term \(\lambda_{\rm poly}\) to all polynomial expansion columns (indices \(n_s:\)) while leaving the \(n_s\) linear columns unpenalised. This prevents the polynomial cross-terms from absorbing signal that belongs to the linear coefficients — which are the only ones used in the weight update — resolving the fixed-point inconsistency that caused ISD-3 to converge to implausible solutions. Recommended value: 1e-3 * np.var(delta_g_obs). Default 0.0 (backward-compatible).

  • Interface. run_decontamination() gains isd_lambda_poly: float = 0.0, forwarded to both the pass-1 and pass-2 ISD calls.

  • Documentation. Methods Reference ISD section rewritten to document the normal-equations formulation, the fixed-point inconsistency root cause, and the lambda_poly parameter with its mathematical derivation.

0.2.3 (2026-05-06)

Pipeline refactoring — method-first execution:

  • New orchestration script scripts/run_all_methods_sequential.sh runs the LS10 BGS dataset for each fitting method in speed order: OLS (seconds) → ElasticNet (minutes) → MCMC-add (hours) → MCMC-comb (hours). Sphinx HTML docs are rebuilt after each method phase so intermediate results are immediately browsable.

  • All four analysis scripts gain --only-methods (restrict which methods run in a given invocation) and --force (ignore skip-if-done cache). The shell script exposes FORCE=1 as an environment variable.

  • Skip-if-done cache: run_bin() and run_sample() check for an existing results_partial_{slug}.json (partial run) or results.json (full run) at the start of each bin/sample and return immediately if found. Delete the file or pass --force to recompute.

  • Partial results JSON written after fast-method phases (OLS, ElasticNet, MCMC-add) — lightweight dict with schema_version, methods_run, timestamp_utc, n_sys_*, syst_names, delta_g_std/skew, and per-method sub-dicts. The MCMC-comb phase merges these into the final results.json.

  • Incremental RST updates via sentinel-comment protocol: .. _auto-method-OLS-start: / .. _auto-method-OLS-end: markers are written to results_ls10.rst after each method phase; re-running replaces (not duplicates) the section.

Figures — consistent per-method color palette:

  • New module sys_mapping/plotting.py with shared METHOD_COLORS, METHOD_LINESTYLES, METHOD_MARKERS, METHOD_LABELS, and METHOD_ORDER constants exported from the package. The palette is Tableau-10-inspired and colorblind-safe: OLS (#4E79A7 steel blue), ElasticNet (#59A14F forest green), ISD-1 (#F28E2B amber), ISD-3 (#E15759 coral red), MCMC-add (#B07AA1 purple), MCMC-comb (#76B7B2 teal).

  • run_ls10_analysis.py: 2PCF diagnostic plot updated to use METHOD_COLORS.

0.2.2 (2026-05-06)

  • Bug fix (runtime): scripts/run_mock_analysis_real_templates.py — both ISD blocks computed np.asarray(alpha_hat[-1])[:n_sys] instead of np.asarray(alpha_hat)[:n_sys]. alpha_hat is a 1-D array of shape (n_expanded,); alpha_hat[-1] is a scalar, and slicing it raises IndexError. The surrounding except Exception silently swallowed the error, so ISD results were always missing from this script. Fixed by removing the erroneous [-1] index.

  • Bug fix (consistency): scripts/run_systematic_tests.py — ElasticNet discarded the weight map returned by elasticnet_contamination_fit() and applied additive correction δ_obs − α̂·t instead of weight-based w(p)·(1+δ_obs)−1. Fixed to store the weights in the results dict so that corrected_delta() applies the paper-correct correction consistent with run_validation.py and Weaverdyck & Huterer (2020) Eq. 46-47.

  • Consistency: scripts/run_systematic_tests.py — OLS and MCMC-additive now also use weight-based correction w(p) = 1/(1 + α̂·t(p)), making all six methods consistent with Weaverdyck & Huterer (2020) Eq. 46-47 and with run_validation.py. The _compute_weights() helper is used directly.

  • Documentation: re-ran run_systematic_tests.py, run_validation.py, run_mock_analysis.py, and run_mock_analysis_progressive.py; updated all result tables throughout the docs.

0.2.1 (2026-05-06)

  • Bug fix (performance): iterative_systematics_decontamination was allocating an explicit (n_pix × n_pix) diagonal weight matrix via np.diag(weights) at every iteration. For realistic pixel counts (NSIDE ≥ 64, n_pix ≳ 50 000) this consumed tens of GB and made each ISD iteration orders of magnitude slower than necessary. Fixed by replacing X.T @ np.diag(w) with the equivalent broadcast X.T * w[np.newaxis, :] (regression.py).

  • Bug fix (correctness): The lstsq fallback path inside the ISD loop was scaling both sides by weights instead of sqrt(weights), minimising \(\sum_p w(p)^2 (\delta - \alpha t)^2\) rather than the correct weighted least-squares objective \(\sum_p w(p)(\delta - \alpha t)^2\). Fixed by using sw = sqrt(weights) on both the design matrix and the response (regression.py).

  • Minor: moved X = delta_t_expanded.T outside the ISD iteration loop (it is constant across iterations).

  • Validation: re-ran scripts/run_validation.py at NSIDE=32; all ISD metrics are numerically identical to the pre-fix run (the bugs affected performance, not output for the converging cases in the validation suite). ElasticNet results added to all scenario tables in docs/results_validation.rst.

0.2.0 (2026-05-05)

  • Scripts: converted all notebook-based analyses to standalone Python scripts in scripts/ (build_systematic_maps, run_paper_validation, run_ls10_analysis, run_mock_analysis).

  • Execution: added SLURM job templates in slurm/ and nohup bash wrappers in bash/ for 32-core server execution.

  • Bug fix: corrected template 2PCF estimator and covariance propagation (commit 10c01d0); all prior results should be regenerated.

  • Environment: switched from conda to mamba ($HOME/miniforge3/envs/sys_map).

  • Notebooks: deprecated executed/bloated notebooks moved to notebooks/deprecated/; clean versions kept for interactive exploration.

  • Corrfunc: added measure_two_point_function_corrfunc utility and benchmarking script scripts/benchmark_corrfunc_vs_treecorr.py.

0.1.0 (2026-04-23)

  • Initial release.

  • Joint inference of additive and multiplicative systematics via JAX JIT-compiled likelihoods.

  • emcee MCMC backend (250 walkers, 1500 steps).

  • Noise debiasing (Eq. 21) and two-point function correction (Eq. 15-16).

  • PCA template rotation (Appendix A).

  • Likelihood ratio model selection (Eq. 19).

  • Block bootstrap variance estimation (Sec. 6.2).

  • Five HEALPix systematic template families (Sec. 4.2).

  • TreeCorr Landy-Szalay two-point function measurement.