pith. sign in

arxiv: 2510.05215 · v2 · submitted 2025-10-06 · 🌌 astro-ph.CO

QML-FAST -- A Fast Code for low-ell Tomographic Maximum Likelihood Power Spectrum Estimation

Yurii Kvasiuk , Anderson Lai , Moritz M\"unchmeyer , Kendrick M. Smith This is my paper

Pith reviewed 2026-05-18 09:04 UTC · model grok-4.3

classification 🌌 astro-ph.CO
keywords power spectrum estimationquadratic maximum likelihoodtomographic analysislarge-scale structurecosmic microwave backgroundkSZ reconstruction
0
0 comments X

The pith

A fast QML estimator performs fully optimal power spectrum estimation for multiple correlated fields at low multipoles.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The paper presents QML-FAST, a fast implementation of the quadratic maximum likelihood estimator for power spectra of multiple correlated scalar fields on the sphere. The code supports arbitrary binning in redshift and multipoles along with cross-correlations and uses a pixel-wise covariance model to achieve full optimality. Optimizations make the estimator tractable for low-ℓ tomographic analyses with many fields, delivering significant precision gains at large scales over the pseudo-C_ℓ method while running in about an hour on standard hardware.

Core claim

The authors present a novel implementation of the quadratic maximum likelihood estimator for power spectra of multiple correlated scalar fields. The estimator uses a pixel-wise covariance model for a fully optimal analysis, supports arbitrary binning in redshift and multipoles, and includes cross-correlations. Through a number of optimizations, the estimator and covariance matrix become computationally tractable for low-ℓ analyses, as demonstrated by constructing estimators on 40 correlated fields up to N_side=32 in about an hour on a single 24-core CPU using less than 256 Gb RAM. Validation on simulations shows significant gains in precision at large scales compared to the pseudo-C_ℓ method

What carries the argument

Quadratic maximum likelihood estimator with pixel-wise covariance model and numerical optimizations for low-ℓ tractability.

If this is right

  • Yields unbiased optimal estimates of auto- and cross-power spectra at large angular scales for tomographic data.
  • Enables efficient computation for analyses with up to 40 correlated fields at N_side=32.
  • Provides higher precision than pseudo-C_ℓ methods at low multipoles for applications like kSZ reconstruction.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

  • The optimizations could extend to modestly higher resolution with further work.
  • Similar techniques may improve other quadratic estimators in cosmology.
  • Demonstrated performance suggests readiness for real survey data.

Load-bearing premise

The numerical optimizations needed to make the computation tractable at N_side=32 with 40 fields preserve the exact optimality and unbiasedness of the quadratic maximum likelihood estimator.

What would settle it

Comparison of power spectra recovered from simulated maps with known inputs to check for bias or loss of optimality at low ℓ.

Figures

Figures reproduced from arXiv: 2510.05215 by Anderson Lai, Kendrick M. Smith, Moritz M\"unchmeyer, Yurii Kvasiuk.

Figure 1
Figure 1. Figure 1: An Nside = 16 binary mask downgraded from Nside = 32 and applied to the kSZ analysis in [9], the unmasked pixels have a value of 1. Mode removal. We investigate the robustness of the mode deprojection described in Sec. 2.5 by removing the QML estimations below ℓ = 5 through three separate methods. 1. Calculate the complete Fisher matrix Fℓℓ′ and the estimator Qℓ, then remove the rows and columns of the Fis… view at source ↗
Figure 2
Figure 2. Figure 2: Comparison of mode-deprojection methods. Left: The full-sky QML estimates and the associated [PITH_FULL_IMAGE:figures/full_fig_p008_2.png] view at source ↗
Figure 3
Figure 3. Figure 3: Comparison of (un-)biasedness of mode-deprojection methods. The QML estimates of the [PITH_FULL_IMAGE:figures/full_fig_p009_3.png] view at source ↗
Figure 4
Figure 4. Figure 4: Influence of the fiducial power on the estimator optimality. The average of 5000 QML estimates, [PITH_FULL_IMAGE:figures/full_fig_p009_4.png] view at source ↗
Figure 5
Figure 5. Figure 5: Test of the pseudo-Cℓ pipeline. Pseudo-Cℓ estimates on the autopower C a,a ℓ from averaging 105 pseudo-Cℓ outputs. The pseudo-Cℓ power spectrum is calculated using PyMaster without apodization to the mask. The red curve shows the fiducial power spectrum convolved with the bandpower window function. There is residual bias beyond ℓ = 46. The spikiness of the fiducial power is due to the window function. The … view at source ↗
Figure 6
Figure 6. Figure 6: Comparison of pseudo-Cℓ and QML error bars. Left: A comparison between pseudo-Cℓ estimate and QML estimate of the autopower C a,a ℓ for the same input masked maps. The pseudo-Cℓ power spectrum is evaluated with the PyMaster pipeline starting at ℓ = 1. ℓ = 0 has been deprojected with a pseudo-inverse covariance for QML output. The grey shaded region highlights the modes that fail to converge for both method… view at source ↗
Figure 7
Figure 7. Figure 7: Influence of fiducial cross-powers on estimator error bars. Left: The ratio of QML error bar to the [PITH_FULL_IMAGE:figures/full_fig_p012_7.png] view at source ↗
Figure 8
Figure 8. Figure 8: Small sky fraction mask. An Nside = 16 binary mask covering a connected patch in the sky with sky fraction is 0.055, the restricted minimum angular size of the unmasked region leads to a singular Fisher matrix. between an initial guess on the fiducial power and the one that generates the observed maps can be addressed through an iterative process, where the fiducial power for the next QML estimator is upda… view at source ↗
Figure 9
Figure 9. Figure 9: Bandpowered analysis. Left: The Fisher element corresponds to the [PITH_FULL_IMAGE:figures/full_fig_p013_9.png] view at source ↗
Figure 10
Figure 10. Figure 10: Iterative QML application. Left: The best-fit modeling parameter [PITH_FULL_IMAGE:figures/full_fig_p014_10.png] view at source ↗
read the original abstract

We present a novel implementation for the quadratic maximum likelihood (QML) power spectrum estimator for multiple correlated scalar fields on the sphere. Our estimator supports arbitrary binning in redshift and multipoles $\ell$ and includes cross-correlations of redshift bins. It implements a fully optimal analysis with a pixel-wise covariance model. We implement a number of optimizations which make the estimator and associated covariance matrix computationally tractable for a low-$\ell$ analysis, suitable for example for kSZ velocity reconstruction or primordial non-Gaussianity from scale-dependent bias analyses. We validate our estimator extensively on simulations and compare its features and precision with the common pseudo-$C_\ell$ method, showing significant gains at large scales. We make our code publicly available. In a companion paper, we apply the estimator to kSZ velocity reconstruction using data from ACT and DESI Legacy Survey and construct full set of QML estimators on 40 correlated fields up to $N_{\text{side}}= 32$ in timescale of an hour on a single 24-core CPU requiring $<256\ \mathrm{Gb}$ RAM, demonstrating the performance of the code.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

2 major / 2 minor

Summary. The manuscript presents QML-FAST, a publicly released code implementing the quadratic maximum likelihood (QML) estimator for tomographic power spectra of multiple correlated scalar fields on the sphere. It supports arbitrary redshift and multipole binning (including cross-bin correlations), adopts a pixel-wise covariance model for a claimed fully optimal analysis, and introduces optimizations to render the high-dimensional covariance (~5×10^5 elements) tractable at low ℓ. Validation on simulations is reported to show significant gains relative to the pseudo-C_ℓ estimator at large scales, with computational performance demonstrated for 40 fields at N_side=32 (processed in ~1 hour on a 24-core CPU with <256 GB RAM) in a companion kSZ application.

Significance. If the optimizations preserve the exact unbiasedness and minimum-variance properties of the textbook QML estimator, the work would provide a practically useful tool for large-scale cosmological analyses such as kSZ velocity reconstruction and scale-dependent bias measurements for primordial non-Gaussianity. The public code release and reported wall-clock performance constitute concrete strengths that lower the barrier for optimal analyses of correlated fields.

major comments (2)
  1. [§2 and §3 (optimizations and covariance construction)] The central claim of a 'fully optimal analysis' (abstract and §2) rests on the pixel-wise covariance inversion remaining exact. The description of 'a number of optimizations' for tractability at N_side=32 with 40 fields must explicitly state whether these consist of exact symmetry reductions/sparse factorizations or tolerance-controlled approximations (iterative solvers, low-rank updates, early stopping). Any surrogate inverse can shift the recovered C_ℓ at the lowest multipoles in a manner not removed by standard debiasing; quantitative tests of this bias (with error bars) are required in the validation section.
  2. [Validation section / simulation results] Table or figure presenting the simulation validation (presumably §4 or §5) reports only qualitative 'significant gains' over pseudo-C_ℓ. Load-bearing for the optimality claim are the actual bias and variance ratios at ℓ ≲ 10, including error bars on the lowest bin; without these numbers the cross-method comparison cannot be assessed.
minor comments (2)
  1. [Abstract] The abstract would benefit from one or two quantitative metrics (e.g., fractional improvement in variance or bias level at ℓ=2–5) to substantiate the 'significant gains' statement.
  2. [§2] Notation for the multi-field covariance matrix and the binning operator should be introduced with an explicit equation early in §2 to aid readability.

Simulated Author's Rebuttal

2 responses · 0 unresolved

We thank the referee for their careful and constructive review of our manuscript. Their comments highlight important points regarding clarity on the exactness of our optimizations and the need for quantitative validation metrics. We address each major comment below and will incorporate revisions to strengthen the paper.

read point-by-point responses

  1. Referee: [§2 and §3 (optimizations and covariance construction)] The central claim of a 'fully optimal analysis' (abstract and §2) rests on the pixel-wise covariance inversion remaining exact. The description of 'a number of optimizations' for tractability at N_side=32 with 40 fields must explicitly state whether these consist of exact symmetry reductions/sparse factorizations or tolerance-controlled approximations (iterative solvers, low-rank updates, early stopping). Any surrogate inverse can shift the recovered C_ℓ at the lowest multipoles in a manner not removed by standard debiasing; quantitative tests of this bias (with error bars) are required in the validation section.

    Authors: We agree that explicit clarification is needed to support the 'fully optimal' claim. Our optimizations rely on exact symmetry reductions of the spherical harmonic basis and sparse matrix factorizations that preserve the exact pixel-wise covariance inversion without introducing approximations, iterative solvers, or low-rank updates. We will revise §§2 and 3 to state this explicitly and will add quantitative bias tests (with error bars) in the validation section demonstrating that recovered C_ℓ values at ℓ ≲ 10 show no systematic shift relative to the input. revision: yes

  2. Referee: [Validation section / simulation results] Table or figure presenting the simulation validation (presumably §4 or §5) reports only qualitative 'significant gains' over pseudo-C_ℓ. Load-bearing for the optimality claim are the actual bias and variance ratios at ℓ ≲ 10, including error bars on the lowest bin; without these numbers the cross-method comparison cannot be assessed.

    Authors: We concur that qualitative statements alone are insufficient for a rigorous assessment. We will expand the validation section (and associated table/figure) to report explicit bias and variance ratios at ℓ ≲ 10, including error bars on the lowest multipole bin, enabling direct quantitative comparison between QML-FAST and the pseudo-C_ℓ estimator. revision: yes

Circularity Check

0 steps flagged

No circularity: implementation and validation of existing QML estimator

full rationale

The paper presents a code implementation of the standard quadratic maximum likelihood estimator for tomographic power spectra, with numerical optimizations for tractability at low ℓ and N_side=32. It validates performance on simulations and compares to pseudo-C_ℓ, but contains no derivation chain, no fitted parameters renamed as predictions, and no load-bearing self-citations that reduce the central claims to inputs by construction. The work is self-contained as a computational tool release with external simulation benchmarks.

Axiom & Free-Parameter Ledger

0 free parameters · 1 axioms · 0 invented entities

The central contribution is a software implementation of a known estimator class; no new free parameters, axioms, or invented entities are introduced beyond standard spherical harmonic and covariance machinery.

axioms (1)
  • domain assumption The pixel-space covariance model for the observed fields is known and correctly specified.
    Required for the quadratic maximum likelihood estimator to be optimal; stated in the abstract as 'fully optimal analysis with a pixel-wise covariance model'.

pith-pipeline@v0.9.0 · 5745 in / 1388 out tokens · 27934 ms · 2026-05-18T09:04:27.863991+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

15 extracted references · 15 canonical work pages · 4 internal anchors

  1. [1]

    The Pseudo-$C_l$ method: Cosmic microwave background anisotropy power spectrum statistics for high precision cosmology

    Benjamin D. Wandelt, Eric Hivon, and Krzysztof M. Gorski. “The pseudo-c l method: cosmic microwave background anisotropy power spectrum statistics for high precision cosmology”. In:Phys. Rev. D64 (2001), p. 083003. arXiv:astro-ph/0008111. 14

    work page internal anchor Pith review Pith/arXiv arXiv 2001

  2. [2]

    How to measure CMB power spectra without losing information

    Max Tegmark. “How to measure CMB power spectra without losing information”. In:Physical Review D 55.10 (May 1997), pp. 5895–5907.issn: 1089-4918

    work page 1997

  3. [3]

    How to measure CMB polarization power spectra without losing information

    Max Tegmark and Angelica de Oliveira-Costa. “How to measure CMB polarization power spectra without losing information”. In:Physical Review D64.6 (Aug. 2001).issn: 1089-4918

    work page 2001

  4. [4]

    OPTIMIZED LARGE-SCALE CMB LIKELIHOOD AND QUADRATIC MAXIMUM LIKELIHOOD POWER SPECTRUM ESTIMATION

    E. Gjerløw et al. “OPTIMIZED LARGE-SCALE CMB LIKELIHOOD AND QUADRATIC MAXIMUM LIKELIHOOD POWER SPECTRUM ESTIMATION”. In:The Astrophysical Journal Supplement Series 221.1 (Oct. 2015), p. 5.issn: 1538-4365

    work page 2015

  5. [5]

    Quadratic estimator for CMB cross-correlation

    S. Vanneste et al. “Quadratic estimator for CMB cross-correlation”. In:Phys. Rev. D98.10 (2018), p. 103526. arXiv:1807.02484 [astro-ph.CO]

    work page internal anchor Pith review Pith/arXiv arXiv 2018

  6. [6]

    ECLIPSE: a fast Quadratic Maximum Likelihood estimator for CMB intensity and polarization power spectra

    J.D. Bilbao-Ahedo et al. “ECLIPSE: a fast Quadratic Maximum Likelihood estimator for CMB intensity and polarization power spectra”. In:Journal of Cosmology and Astroparticle Physics2021.07 (July 2021), p. 034.issn: 1475-7516

    work page 2021

  7. [7]

    Testing quadratic maximum likelihood estimators for forthcoming Stage-IV weak lensing surveys

    Alessandro Maraio, Alex Hall, and Andy Taylor. “Testing quadratic maximum likelihood estimators for forthcoming Stage-IV weak lensing surveys”. In:Monthly Notices of the Royal Astronomical Society520.4 (Feb. 2023), pp. 4836–4852.issn: 1365-2966

    work page 2023

  8. [8]

    Quadratic estimators for unwindowed power spectrum of galaxy-galaxy weak lensing and its application toP gm(k) estimation

    Taisei Terawaki, Masahiro Takada, and Takanori Taniguchi. “Quadratic estimators for unwindowed power spectrum of galaxy-galaxy weak lensing and its application toP gm(k) estimation”. In: (July 2025). arXiv: 2507.18789 [astro-ph.CO]

    work page arXiv 2025

  9. [9]

    KSZ Velocity Reconstruction with ACT and DESI-LS using a Tomographic QML Power Spectrum Estimator

    Anderson C. M. Lai, Yurii Kvasiuk, and Moritz M¨ unchmeyer. “KSZ Velocity Reconstruction with ACT and DESI-LS using a Tomographic QML Power Spectrum Estimator”. In: (June 2025). arXiv:2506.21684 [astro-ph.CO]

    work page arXiv 2025

  10. [10]

    A method for extracting maximum resolution power spectra from microwave sky maps

    M. Tegmark. “A method for extracting maximum resolution power spectra from microwave sky maps”. In:Monthly Notices of the Royal Astronomical Society280.1 (May 1996), pp. 299–308.issn: 1365-2966

    work page 1996

  11. [11]

    opt einsum - A Python package for optimizing contraction order for einsum-like expressions

    Daniel G. a. Smith and Johnnie Gray. “opt einsum - A Python package for optimizing contraction order for einsum-like expressions”. In:Journal of Open Source Software3.26 (2018), p. 753

    work page 2018

  12. [12]

    Numba: A llvm-based python jit compiler

    Siu Kwan Lam, Antoine Pitrou, and Stanley Seibert. “Numba: A llvm-based python jit compiler”. In: Proceedings of the Second Workshop on the LL VM Compiler Infrastructure in HPC. 2015, pp. 1–6

    work page 2015

  13. [13]

    healpy: equal area pixelization and spherical harmonics transforms for data on the sphere in Python

    Andrea Zonca et al. “healpy: equal area pixelization and spherical harmonics transforms for data on the sphere in Python”. In:Journal of Open Source Software4.35 (Mar. 2019), p. 1298

    work page 2019

  14. [14]

    M., Hivon, E., Banday, A

    K. M. G´ orski et al. “HEALPix: A Framework for High-Resolution Discretization and Fast Analysis of Data Distributed on the Sphere”. In:ApJ622 (Apr. 2005), pp. 759–771. eprint:arXiv:astro-ph/0409513

    work page internal anchor Pith review Pith/arXiv arXiv 2005

  15. [15]

    Erminia Calabrese et al.The Atacama Cosmology Telescope: DR6 Constraints on Extended Cosmological Models. 2025. arXiv:2503.14454 [astro-ph.CO]. A Variance due to mismodeling of the signal covariance Here, we show that the QML estimator, built withC f id ̸=C true is suboptimal. Let us •assume we mismodel the signal covariance, but not the noise, so that th...

    work page internal anchor Pith review Pith/arXiv arXiv 2025