{"id":"9a49cd3f-0271-435c-9fbd-fbf4a28102b9","arxiv_id":"2506.10526","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":2,"one_line_summary":"A gridded-visibility, FFT-accelerated estimator for the binned 21 cm multi-frequency angular bispectrum and 3D bispectrum is presented and validated against simulations with an input bispectrum model.","lead":"This paper builds a fast computer tool that estimates a three-point statistic (the bispectrum) of the cosmic 21 cm signal from radio telescope visibility data, covering all triangle shapes in one pass. It tests the tool on simulated Murchison Widefield Array observations and finds the recovered bispectrum matches the input model within expected noise, even when some frequency channels are flagged.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The 3D BS recovery rests on frequency ergodicity (Eq. 10), but the validation uses a periodic, statistically homogeneous LoS; a light-cone test is needed before generalizing to real 30.72 MHz data.","rationale":"The paper's central claim is a methods claim: a visibility-based estimator for the binned MABS and 3D 21 cm bispectrum, validated against a known input model. For the MABS component, the estimator does not require frequency ergodicity and the gridding/FFT construction is clearly described. For the 3D BS component, however, the route from MABS to B(k⊥,k∥) is Eq. (10), which is derived under statistical homogeneity and ergodicity along the LoS. The validation deliberately simulates a periodic, stationary 3D box, so it exercises exactly the assumption that makes Eq. (10) exact and provides no information about how the estimator behaves when that assumption fails. This is the weakest load-bearing point because the abstract and summary claim reliable recovery of the 3D BS for observationally relevant bandwidths where the signal is known to evolve across redshift. I considered two other potential concerns: the effective damping factor exp(-π²θ0²ΔU²/3)≈0.89 from gridding, and the computational cost of evaluating Eq. (18) over all frequency and ring triplets. Both are worth checking, but the text does not give enough implementation detail to establish either as a definite internal inconsistency, whereas the ergodicity limitation is explicit in the formalism and absent from the validation. The proposed light-cone test is the single check that would settle whether the 3D BS claim generalizes beyond periodic-box simulations. The reader already identified this as the weakest assumption, and the paper's own limitation statements do not cover it explicitly. I therefore keep the reader's CONDITIONAL verdict unchanged.","tokens_in":22118,"tokens_out":18351,"duration_ms":245918,"concrete_test":"Run the estimator on a light-cone mock: construct a 3D frequency-angle cube from the same f_NG local transform but with an input P(k,z) that changes across the 30.72 MHz band (e.g., splice coeval boxes with a known evolving power spectrum), generate mock visibilities with the actual MWA sampling and flagging, and compare the recovered cylindrical B(k⊥,k∥) with the prediction obtained by numerically applying Eq. (10) to the ensemble-averaged MABS of that light-cone model. Passing criterion: fractional deviations at or below the ~20% level and |Δσ|≤3 in the bins where the periodic-box validation passed; any systematic k∥-dependent excess indicates light-cone bias.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Section 2.2 derives the central 3D BS relation (Eq. 10) from the MABS under the explicit assumption that the signal is ergodic along frequency, i.e., that the MABS depends only on Δν1 and Δν2, not on the absolute frequency ν3. Section 4 validates the estimator in exactly this regime: the input non-Gaussian field is generated in a 2048^3 periodic box with P(k)=k^-2 and then sliced into frequency channels, so the LoS is statistically homogeneous and the 2D Fourier transform in Eq. (19) is exact by construction. For the real MWA band, ν_c=154.25 MHz with B_bw=30.72 MHz corresponds to z≈8.2 and edges z≈7.4 and z≈9.2, where the 21 cm signal evolves substantially. Under such light-cone evolution the MABS depends on absolute frequencies, and Eq. (10) no longer holds exactly; the k∥ modes recovered by the DFT become a weighted mix of true modes. The paper acknowledges that the MABS itself does not require ergodicity, and it lists foregrounds and noise as omissions, but it nowhere tests or quantifies the bias from LoS evolution. Thus the headline claim of a reliable 3D BS estimator is load-bearing on an assumption that the presented validation never stresses. This is a limitation, not an internal inconsistency: the MABS estimator and the periodic-box validation are self-consistent, and the 3D BS claim should be read as conditional on frequency stationarity.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper extends the single-frequency angular-bispectrum estimator of Paper I to multifrequency observations. It defines a binned estimator for the multi-frequency angular bispectrum (MABS) from gridded visibilities using annular rings and FFTs, and then obtains the 3D 21-cm bispectrum by a 2D DFT over frequency separations under the assumption of frequency ergodicity. The estimator is validated with 250 realizations of a non-Gaussian field with a known analytic bispectrum, using an MWA drift-scan baseline distribution, both with and without the observed flagging pattern. The paper reports agreement with analytic predictions at the ≲20% level over k1 in [0.003, 1.258] Mpc^{-1} and concludes that the deviations are mostly consistent with statistical fluctuations.","tokens_in":22405,"tokens_out":27798,"duration_ms":332846,"significance":"The estimator is a useful step toward measuring 21-cm non-Gaussianity: it is presented as the first visibility-based MABS estimator and it covers all triangle shapes in a binned sense. The validation strategy is honest in an important respect: no free parameters are fitted to the target bispectrum, and the analytical predictions are binned over the same discrete k modes as the estimates. The flagging-robustness test is also valuable. However, two issues currently limit the strength of the central claims: as written, the implemented estimator in Eq. (18) omits the non-closure exponential factor that appears in the defining relation Eq. (14), and the 3D-bispectrum part is validated only under a frequency-stationarity assumption that real light-cone data will violate. Both issues are fixable, but they are load-bearing for the unbiasedness and generalizability claims.","major_comments":[{"comment":"Eq. (14) contains the non-closure correction factor exp(π^2 θ0^2 ΔU^2/3), but the implemented binned estimator in Eq. (18) contains only the prefactor 1/A and does not apply this exponential. Section 4 explicitly states that the factor is 0.89 for a typical value of (ΔU)^2 = (ΔU_g)^2/2. If the exponential is not applied, the estimator is biased low by roughly 5–11%, depending on the typical ΔU distribution. This is not negligible relative to the claimed 20% accuracy, and it is a natural explanation for the 4–5σ deviations reported at high k1 in the no-flagging case (Section 5, Figs. 7–8). The authors should include the correction, e.g., as a bin-dependent weight, or quantify its residual effect and show that it is below the statistical error bars.","section":"Section 3.2, Eqs. (14) and (18)"},{"comment":"The 3D BS estimator is built on the assumption that the signal is ergodic along frequency, so that the MABS depends only on (Δν1, Δν2) and Eq. (10) is an exact 2D Fourier relation. The validation in Section 4 uses a periodic 2048^3 box with a statistically homogeneous line of sight, so this assumption is exact by construction and is never stressed. For the real MWA band (ν_c = 154.25 MHz, B_bw = 30.72 MHz, z ≈ 7.4–9.2) the 21 cm signal evolves substantially across the band; under light-cone evolution the MABS depends on absolute frequencies and the recovered k∥ modes become a weighted mix of true modes. Since the headline claim concerns the 3D BS for the EoR, the authors should either add a light-cone simulation test that quantifies the bias, or explicitly restrict the 3D BS claim to a regime where frequency ergodicity holds and discuss the expected bias on the full band.","section":"Section 2.2, Eq. (10); Section 4"},{"comment":"The text states that for the largest k1 bins (0.413–1.258 Mpc^{-1}) and linear triangles, most no-flagging estimates have 4 < Δσ ≤ 5, while the flagged case shows smaller Δσ. With 250 independent realizations, many 4–5σ deviations in a single configuration are not plausibly statistical, and the decrease of Δσ when flagging is added suggests a fixed systematic offset rather than a purely statistical effect. The conclusion that the deviations are 'mostly consistent with the expected statistical fluctuations' needs quantitative support: the authors should report the full Δσ distribution, test whether the analytic binning uses the same effective weights as the estimator (especially near degenerate/linear triangles where B̄ varies rapidly within a bin), and either correct or explicitly exclude the discrepant bins.","section":"Section 5, Figs. 7 and 8"}],"minor_comments":[{"comment":"The cosmological parameters used to evaluate r and r' are not stated, although all k⊥ and k∥ values quoted in Section 5 depend on them.","section":"Section 4"},{"comment":"Calling (μ, t) = (0.55, 0.95) 'equilateral' is imprecise; the exact equilateral configuration is (μ, t) = (0.5, 1). Please refer to this bin as the one nearest to equilateral.","section":"Section 5, Fig. 5 caption"},{"comment":"The statement that U = 250 corresponds to an angular scale of 0.115° appears to be off by a factor of two; for a baseline U, the angular scale is approximately 1/U in radians, i.e., about 0.23° for U = 250.","section":"Section 4"},{"comment":"The claimed computational complexity O(N_c^3 N_t^2 log N_t^2) omits the number of radial rings N_r; the triple-ring combinations in Eq. (18) contribute an additional factor N_r^3, which should be stated.","section":"Section 3.2 and Section 6"},{"comment":"The description of the flagging after collapsing four frequency channels is ambiguous when a collapsed block is only partially flagged; please specify how such blocks are treated.","section":"Section 4"},{"comment":"The definition of Δσ in the captions is written as '|B̄ − [B̄]_d/σ', with the absolute value left open; it should read '|B̄ − [B̄]_d|/σ'.","section":"Figs. 7 and 8 captions"}],"recommendation":"major_revision","confidential_remarks":"The main technical concern is the apparent omission of the non-closure exponential between Eqs. (14) and (18); this is a concrete, fixable issue and may explain the high-k Δσ patterns in Section 5. The light-cone/ergodicity concern is also valid and should be addressed with a dedicated simulation rather than only a caveat. I do not see grounds for rejection; the paper is a solid methods contribution if these two points are handled and the claimed error consistency is made quantitative."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Bottom line: this is a competent methods paper that delivers what it promises—a visibility-based estimator for the binned multi-frequency angular bispectrum and the 3D 21 cm bispectrum, with an FFT acceleration that makes full triangle coverage feasible. The novelty is incremental over Paper I, but the extension to multi-frequency data and the 3D bispectrum is real, and the paper is the first to offer a visibility-based route to all triangle configurations rather than just equilateral/isosceles.\n\nWhat it does well: the estimator is clearly formulated, the binning scheme is sensible, and the validation is self-consistent. The sky is generated from a known non-Gaussian model, the analytic bispectrum is derived separately, and the estimates are compared to identically binned predictions with no fitted parameters. The flagging test is a nice touch; recovering the bispectrum within 20% despite missing channels is a meaningful robustness check. The paper also acknowledges that foregrounds, noise, and systematics are not included.\n\nThe soft spots are real but not fatal. The 3D bispectrum recovery in Eq. (10) assumes the signal is ergodic along frequency, so the MABS depends only on frequency separations. The validation uses a statistically homogeneous, periodic box, so it never tests light-cone evolution. For a 30.72 MHz band centered at z≈8.2, the signal evolves appreciably across the band, and the Fourier relation will mix k∥ modes. The paper notes the assumption but does not quantify the bias. That means the 3D BS claim is conditional on frequency stationarity; the MABS estimator itself does not carry that burden. Second, there is no code or data release; 'available on reasonable request' is weak for a methods paper. Third, the 4–5 sigma deviations at small k1 and squeezed triangles are attributed to beam convolution and binning, but the explanation is brief and those cases deserve a closer look.\n\nWho this is for: anyone working on 21 cm bispectrum estimators from interferometric data. It is a solid methodological step, not a science result. I would send it to peer review; the referee should push for a light-cone test or an explicit estimate of the evolution bias, and for code release. If those are addressed, this will be a useful reference.","headline":"A solid, self-consistent extension of the visibility-based bispectrum estimator to multi-frequency and 3D, but the 3D recovery rests on a frequency-ergodicity assumption that the validation never stresses.","tokens_in":22958,"tokens_out":2562,"would_cite":true,"duration_ms":29537,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":[],"model":"deepseek-v4-flash","headline":"A visibility-based, FFT-accelerated estimator recovers the binned multi-frequency angular bispectrum and the 3D 21 cm bispectrum from radio-interferometric data, validated on simulated MWA observations with deviations below 20 percent.","keywords":["21 cm cosmology","epoch of reionization","bispectrum estimator","radio interferometry","multi-frequency angular bispectrum","Murchison Widefield Array","FFT-based estimator","non-Gaussianity"],"falsifier":"Apply the estimator to a simulated 21 cm light cone whose input bispectrum is known to evolve across the band; if the recovered 3D bispectrum deviates from the true band-averaged values by more than the quoted 20 percent, the ergodicity assumption behind the frequency-to-parallel-wavenumber transform is falsified.","tokens_in":21904,"feed_emoji":"📡","tokens_out":9026,"duration_ms":89110,"temperature":0.7,"pith_summary":"The paper claims to provide the first visibility-based estimator for the binned multi-frequency angular bispectrum (MABS) and the 3D 21 cm bispectrum that covers all triangle configurations. It works directly on gridded interferometric visibilities and uses FFT acceleration so the full computation is feasible on a workstation. The authors validate it with simulated Murchison Widefield Array observations built from a known non-Gaussian input model, both with and without the frequency-flagging pattern of real data. The recovered bispectra agree with analytical predictions within 20 percent, and most deviations are consistent with expected statistical fluctuations. If the claim holds, it gives the epoch-of-reionization community a practical tool for measuring the highly non-Gaussian 21 cm signal rather than only its power spectrum.","feed_headline":"FFT-based estimator recovers the full 21 cm bispectrum","feed_subtitle":"Simulated MWA data match analytic predictions within 20 percent even with flagged channels, opening a window on reionization.","key_machinery":"The central object is the binned MABS estimator, Eq. (18): for three annular rings $(a_1,a_2,a_3)$ in three frequency channels, it forms the product $D(\\ell_1,\\nu_1,\\boldsymbol{\\theta})D(\\ell_2,\\nu_2,\\boldsymbol{\\theta})D(\\ell_3,\\nu_3,\\boldsymbol{\\theta})$, where each $D$ is an inverse FFT of the gridded visibilities restricted to one ring, normalizes by the corresponding product of $I$ functions (inverse FFTs of the weighted sampling), and divides by $A=\\pi\\theta_0^2Q^3/3$. This implements the three-visibility correlation of Eq. (14) that equals the MABS. A 2D DFT over frequency separations, Eq. (19), then maps the MABS to the 3D cylindrical bispectrum $B(k_{1\\perp},k_{2\\perp},k_{3\\perp},k_{1\\parallel},k_{2\\parallel})$, assuming the signal is ergodic along frequency. The FFT structure reduces the cost from $O(N_c^3 N_t^4)$ to $O(N_c^3 N_t^2 \\log N_t^2)$.","core_discovery":"The central claim is that the binned MABS estimator of Eq. (18) — a normalized product of three inverse FFTs of gridded visibilities restricted to annular rings at three frequencies — is an unbiased estimate of the MABS, and that a 2D discrete Fourier transform over frequency separations (Eq. 19) recovers the cylindrical 3D bispectrum under ergodicity along the frequency axis. The validation uses 250 independent realizations of a simulated MWA pointing with a known input bispectrum. The estimated monopole bispectrum agrees with analytical predictions across $0.003\\,\\mathrm{Mpc}^{-1} \\leq k_1 \\leq 1.258\\,\\mathrm{Mpc}^{-1}$ and a wide range of triangle shapes; fractional deviations stay below 20 percent even when the simulated data have exactly the same flagged frequency channels as the actual MWA observations, and most deviations lie within $1\\sigma$ to $3\\sigma$ of the expected statistical fluctuations.","pith_inferences":["Because the validation box is periodic and statistically homogeneous along frequency, the published claims do not yet cover light-cone evolution; applying the estimator to wide-band real data will likely require splitting the band into redshift windows in which ergodicity approximately holds.","Foregrounds are absent from this validation. In real data, spectrally smooth foregrounds concentrated in the wedge could still contribute to the FFT products, so foreground avoidance or subtraction will likely be needed before the estimated bispectrum can be interpreted cosmologically.","Combining this bispectrum estimator with a power-spectrum estimate from the same gridded visibilities could break degeneracies among reionization parameters (for example bubble size versus ionizing efficiency) that the power spectrum alone cannot separate.","The formalism already yields all non-zero multipole moments of the bispectrum; measuring them on real data would directly probe redshift-space distortions and line-of-sight anisotropy, a step the paper leaves for future work."],"forward_implications":["Provides, for the first time, a visibility-based estimator for the binned MABS and the 3D 21 cm bispectrum that covers all triangle configurations, not just equilateral and isosceles shapes.","Brings the computational cost down to $O(N_c^3 N_t^2 \\log N_t^2)$; the full-band MWA simulation analysed here takes about one hour on a 16-core CPU.","Recovers the bispectrum monopole over $k_1 \\in [0.003, 1.258]\\,\\mathrm{Mpc}^{-1}$ for a wide range of triangle shapes, with fractional deviations below 20 percent even under the periodic flagging pattern of actual MWA data.","Provides the practical route to measuring the non-Gaussian 21 cm signal from the epoch of reionization with interferometers such as MWA, complementing power-spectrum measurements."],"supporting_citations":[{"why":"Supplies the single-frequency angular bispectrum estimator whose gridding, ring-binning, and normalization are generalized here to multiple frequencies; it is the methodological foundation being extended.","marker":"Paper I"},{"why":"Establishes the 2D Fourier-transform relation (Eq. 10) between the MABS and the 3D 21 cm bispectrum, the theoretical bridge the estimator relies on.","marker":"Bharadwaj & Pandey 2005"},{"why":"Provides the MWA drift-scan pointing and baseline distribution used to set up the simulated visibility data and the corresponding flagging pattern.","marker":"Patwa et al. (2021)"},{"why":"Defines the multi-frequency angular power spectrum (MAPS) and its frequency-decorrelation behaviour, the two-point analogue whose framework the MABS construction builds upon.","marker":"Datta et al. (2007)"},{"why":"Introduces the $(k_1,\\mu,t)$ parametrization of triangle shape and the spherical-harmonic multipole expansion used to present the monopole bispectrum results.","marker":"Bharadwaj et al. (2020)"},{"why":"Supplies the FFT-based technique that reduces the computational cost of summing closed-triangle bispectrum estimators, which is adopted for the gridded visibility products.","marker":"Sefusatti et al. (2006)"},{"why":"Provides additional FFT bispectrum estimation formalism used to accelerate the computation across frequency channels.","marker":"Scoccimarro (2015)"},{"why":"Gives the discrete version of the multipole integral and the earlier single-frequency analysis that the present paper extends, including the validation procedure.","marker":"Gill & Bharadwaj 2024"}],"fun_headline_variants":["FFT-accelerated estimator maps 21 cm bispectrum","Fast 21 cm bispectrum estimator validated on MWA data","New bispectrum estimator captures reionization signal","Bispectrum estimator robust to flagged channels, aids EoR"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The method assumes the 21 cm signal is statistically uniform along the frequency axis, so that the 3D bispectrum can be obtained by a 2D Fourier transform over frequency separations; if the signal evolves appreciably across the 30.72 MHz band, that transform introduces a bias.","fun_headline_variants_meta":{"raw":{"variants":["FFT-accelerated estimator maps 21 cm bispectrum","Fast 21 cm bispectrum estimator validated on MWA data","New bispectrum estimator captures reionization signal","Bispectrum estimator robust to flagged channels, aids EoR"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000302,"raw_usage":{"total_tokens":1771,"prompt_tokens":1012,"completion_tokens":759,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":628,"completion_tokens_details":{"reasoning_tokens":698}},"tokens_in":628,"tokens_out":759,"duration_ms":7788,"temperature":1.0,"reasoning_tokens":698,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T04:24:29.786760+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Apply the estimator to a simulated 21 cm light cone whose input bispectrum is known to evolve across the band; if the recovered 3D bispectrum deviates from the true band-averaged values by more than the quoted 20 percent, the ergodicity assumption behind the frequency-to-parallel-wavenumber transform is falsified.","supporting_citations":[],"review_version":1}