{"id":"3b5820d9-9b73-4d21-8cc6-d1e4e44ba9dd","arxiv_id":"2412.15858","paper_version":4,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":2,"one_line_summary":"A new Julia package, Vela.jl, provides an independent, parallelized Bayesian pulsar timing and noise analysis implementation with a Python interface, validated against PINT and tempo2.","lead":"Researchers built Vela.jl, a new Julia package for Bayesian analysis of pulsar timing data, with a Python interface. It independently implements the full timing and noise model and was checked against existing tools, agreeing at the ten-nanosecond level.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The validation's load-bearing premise is that posterior agreement with pint.bayesian is not driven by 'cheat' priors; the paper's only support is a prose robustness claim with no quantitative comparison.","rationale":"The reader's weakest_assumption identifies exactly the same load-bearing concern: the validation examples rely on 'cheat' priors and the claimed robustness check is only described in prose. My stress-test agrees with this identification. The central claim of the paper is that Vela.jl provides an independent, efficient implementation of the full non-linear pulsar timing and noise model, and that its Bayesian results agree with pint.bayesian and recover injected parameters. The residual agreement in Section 3.1 supports the deterministic timing model but not the full Bayesian pipeline. The posterior agreement in Section 4.1 is the primary evidence for the Bayesian pipeline, and it is vulnerable to prior-driven circularity. The paper's own Appendix B acknowledges this risk. Because the robustness check is unquantified, the validation is incomplete, not demonstrably wrong. Therefore the verdict should remain CONDITIONAL, pending a quantitative demonstration that the posterior is insensitive to prior width and prior location. I do not see grounds to reject the paper or to declare it unverdictable; the software is released, the code is open, and the residual comparison is a meaningful positive check. The recommended action is to keep the conditional verdict and request the additional validation as a condition of acceptance.","tokens_in":21972,"tokens_out":2326,"duration_ms":22952,"concrete_test":"Re-run the Section 4.1 analysis (PSR J1748-2021E) with two modified prior settings: (a) widen all 'cheat' priors by a factor of 100 (or replace them with broad physically motivated priors), and (b) center the 'cheat' priors on values deliberately offset from the PINT MLE by 5 sigma. Compare the resulting posterior distributions with the published posteriors in Figure 7 using a quantitative metric, such as the shift in posterior mean in units of posterior standard deviation for each parameter. If the posterior shifts significantly or fails to recover injected values in either run, the Figure 7 agreement is prior-driven and the Bayesian validation in Section 4.1 is not probative.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central claim that Vela.jl correctly implements Bayesian non-linear pulsar timing and noise analysis rests on two validation examples (Sections 4.1 and 4.2). Both use 'cheat' priors: uniform distributions centered on PINT maximum-likelihood estimates, with widths 10-40 times the frequentist uncertainties. Appendix B explicitly warns that such priors 'effectively do circular analysis' when the data lacks signal-to-noise. If the likelihood were incorrect but the prior were narrow enough to cover the true parameters, the posterior would still concentrate near the prior center, producing apparent agreement with pint.bayesian and apparent recovery of injected parameters. The paper's response is a single unquantified sentence in each section: 'We have checked that increasing the width of the 'cheat' priors does not appreciably alter the posterior distribution.' No widened-prior posterior plots, summary statistics, or comparisons are provided. Thus the most direct evidence that the likelihood and sampling are correct is not independently established. The residual agreement with PINT and tempo2 (Section 3.1, Figure 1) is a credible check of the deterministic timing model, but it does not validate the noise-model likelihood, the marginalization, or the sampler integration. The J1909-3744 example partially recovers parameters and the paper attributes the discrepancies to DM model inadequacy, which further limits the strength of the Bayesian validation.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper presents Vela.jl, a Julia package for Bayesian pulsar timing and noise analysis, together with a Python binding called pyvela. Vela.jl implements the nonlinear timing and noise model independently, using PINT only for file input, clock corrections, and solar system ephemeris computations. The manuscript describes the package architecture, numerical-precision choices, component types, prior handling, likelihood kernels, and two validation examples: a residual comparison against PINT and tempo2 at the ~10 ns level, and posterior comparisons against pint.bayesian on a small simulated dataset and on a larger simulated J1909-3744 dataset. The paper argues that Vela.jl is an efficient, modular, and reliable alternative to TEMPONEST and a complement to ENTERPRISE.","tokens_in":22195,"tokens_out":5733,"duration_ms":52583,"significance":"If the central claims are borne out, Vela.jl would be a useful independent tool for Bayesian pulsar timing and noise analysis, with the practical advantages of sampler-agnostic interfaces, multi-threading, and a Python binding. The paper has clear strengths: the residual comparison against PINT and tempo2 is an externally grounded check of the deterministic timing model; the design discussion (extended precision, dimensional types, component hierarchy) is informative; and the package is distributed with version control, documentation, and a test suite. The main weakness is that the Bayesian validation relies on 'cheat' priors centered on PINT maximum-likelihood values, and the stated robustness to prior width is not quantitatively documented. Because this concerns the load-bearing evidence for the likelihood and sampler implementation, the claims as presented are defensible but not yet fully supported.","major_comments":[{"comment":"The Bayesian validation examples both use 'cheat' priors, uniform distributions centered on PINT maximum-likelihood estimates with widths 10–40 times the frequentist uncertainties. Appendix B itself warns that such priors 'effectively do circular analysis' when the data do not provide enough signal-to-noise. The only response in each example is the sentence 'We have checked that increasing the width of the 'cheat' priors does not appreciably alter the posterior distribution,' with no widened-prior plots, summary statistics, or quantitative comparisons. A prior that is narrow enough to concentrate near the PINT or true values can produce apparent agreement with pint.bayesian and apparent parameter recovery even if the likelihood or sampler were incorrect. Please add quantitative evidence: for example, posterior medians and credible intervals for prior widths of 10x, 40x, and 100x, overlap or distance metrics between posteriors, or at least one validation run with broad physically motivated priors.","section":"§4.1, §4.2, Appendix B"},{"comment":"The J1909-3744 example does not cleanly validate the likelihood implementation for the correlated noise model. The data are simulated with epoch-wise DM measurements from InPTA DR1, while the fitted model uses a Taylor-series-plus-40-harmonic Gaussian-process DM model; the text admits that some estimated parameters do not agree well with the injected values and attributes this to inadequacy of the DM model. As presented, the example cannot separate model mismatch from an error in the Gaussian-process likelihood, prior transform, or sampler in the 104-dimensional parameter space. The caption's statement that astrometric and binary parameters are consistent within 3σ does not quantify the noise parameters or the Fourier coefficients. Please provide a full parameter-by-parameter coverage or z-score comparison for all free parameters, and ideally an injection-recovery test in which the data are generated from exactly the same Gaussian-process model as the fit.","section":"§4.2, Figure 8"},{"comment":"The ~10 ns residual agreement against PINT and tempo2 (Figure 1) is a credible check of the deterministic timing model, but it does not exercise the noise-model likelihood, the reduced-rank covariance expression (Eq. 5), the ECORR likelihood (Eq. C5), or the Gaussian-process prior transform (Appendix D). Because the abstract claims a full timing and noise analysis package, the Bayesian validation should include at least one comparison that isolates the stochastic likelihood, for example a comparison with ENTERPRISE or TEMPONEST on simulated data with broad priors and a single well-understood red-noise component.","section":"§3.1"}],"minor_comments":[{"comment":"In the sentence defining the split of F0, the text reads 'F_big ≪ F_small'; this inequality appears to be inverted and should presumably be F_big ≫ F_small. Please correct it.","section":"§3.2"},{"comment":"The par and tim file names in the code sample appear with inserted spaces, for example 'N GC 64 40E . par'. If these spaces are present in the actual code, the example will not run; if they are a line-wrapping artifact, please adjust the figure so the code is unambiguous.","section":"Figure 4"},{"comment":"The captions appear as 'T able 1' and 'T able 2' in the compiled text; please check the typesetting so the table labels are rendered correctly.","section":"Tables 1 and 2"},{"comment":"Equation (D7) is difficult to parse in the provided version: the placement of the reference frequency f_yr is ambiguous. Please ensure the equation is typeset unambiguously and define f_yr explicitly in the text.","section":"Equation (D7)"}],"recommendation":"major_revision","confidential_remarks":"This is a useful software-overview paper, and the deterministic validation is solid. The main gap is the unquantified robustness of the Bayesian validation to the 'cheat' priors, which is a fixable evidence problem rather than a demonstrated error. I would not recommend rejection, but the manuscript should not be accepted until quantitative prior-robustness tests and a cleaner validation of the stochastic likelihood are provided."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"You should know this paper is a software overview for Vela.jl, a Julia package for Bayesian pulsar timing and noise analysis, with a Python binding. The genuinely new thing is the independent implementation of the full non-linear timing model, not the underlying math, which follows Lentati and van Haasteren. What the paper does well is validate the deterministic part of the code: residuals agree with PINT and tempo2 at the ~10 ns level, which is real evidence that the timing model is correctly implemented. The code is open source, versioned on Zenodo, and the documentation is thorough. The author also deserves credit for disclosing that the second example does not recover all injected parameters and attributing that to DM model inadequacy rather than hiding it.\n\nThe soft spot is the Bayesian validation. Both examples use 'cheat' priors: uniform distributions centered on PINT maximum-likelihood values with widths 10–40 times the frequentist uncertainties. Appendix B rightly warns this can become circular analysis when the data lacks SNR. The paper says in one sentence that widening the priors does not appreciably alter the posterior, but no widened-prior plots or summary statistics are shown. That is a legitimate concern, though not fatal: for parameters like F0 and coordinates, 10–40x uncertainties is likely broad enough that the likelihood dominates, and the first example is a small dataset where pint.bayesian agrees with Vela.jl. But the reader cannot verify that from the paper. The efficiency claim is also not benchmarked against tempo2 or ENTERPRISE, which is minor but worth noting.\n\nIn my reading, the central correctness claim holds up, but the validation would be stronger with a quantitative prior-width test. This is the kind of paper that should go to peer review, and the referee should ask for that robustness check. I would not cite it myself in the next year unless I needed a Python-callable non-linear Bayesian timing tool, but colleagues doing PTA noise analysis will likely find it useful. Bring it to reading group if you are interested in how modern pulsar timing software is engineered; otherwise it is a solid, honest contribution that deserves a fair referee.","headline":"A solid software paper for a new Bayesian pulsar timing package, with one honest caveat: the cheat-prior robustness check is asserted in prose but not shown, so the Bayesian validation is weaker than it could be.","tokens_in":22764,"tokens_out":2057,"would_cite":false,"duration_ms":19438,"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":"Independent pulsar-timing code reproduces existing results to 10 ns","keywords":["pulsar timing","Bayesian inference","noise analysis","Julia","pulsar timing arrays","nested sampling","Markov chain Monte Carlo","narrowband timing"],"falsifier":"Run the same two simulated datasets through Vela.jl with broad priors anchored to physical ranges rather than to the reference fitter's point estimates (e.g., flat priors tens to thousands of times wider than the quoted uncertainties), and check whether the posterior medians and the reference-fitter agreement survive; if they shift beyond the quoted uncertainties, the circular-analysis concern is confirmed.","tokens_in":21721,"feed_emoji":"🕰️","tokens_out":12910,"duration_ms":107257,"temperature":0.7,"pith_summary":"This paper presents Vela.jl, a Julia package for Bayesian pulsar timing and noise analysis, along with pyvela, a Python binding. It aims to establish that a from-scratch implementation of the full non-linear timing and noise model can be efficient, parallel, and accurate enough for routine single-pulsar use and for the large datasets of pulsar timing arrays. The evidence comes in two forms: timing residuals that agree with two established pulsar-timing packages at the ~10 ns level, and Bayesian posterior runs on simulated datasets that agree with an established fitter and recover injected parameters. A practical benefit is that Bayesian noise characterization no longer has to be separated from the full timing-model fit, and can be run with whatever sampler the user prefers.","feed_headline":"Independent pulsar-timing code reproduces existing results to 10 ns","feed_subtitle":"Runs the full non-linear timing and noise model in parallel, so Bayesian fits can scale to pulsar timing arrays.","key_machinery":"The load-bearing machinery is a modular, from-scratch implementation of the timing and noise model: TOA delays and phase corrections are represented as Component objects, each with a correction routine, and the log-likelihood is assembled by a Kernel object that evaluates the Gaussian likelihood in equation (4). Computational cost is controlled by a reduced-rank covariance model together with the Woodbury and matrix-determinant lemmas, and by an ECORR block decomposition that evaluates the correlated-noise likelihood in linear time. Numerically, TOA values and rotational phases use the Double64 extended-precision representation, while other quantities use ordinary double precision; the rotational frequency is stored as a sum of two doubles to keep the parameter type uniform. Red-noise Fourier coefficients are reparameterised by their prior standard deviations so they are a priori unit-normal, which avoids hard-to-sample funnel geometries.","core_discovery":"On its own terms, this paper claims that Vela.jl delivers an independent, efficient, parallelized implementation of the full non-linear pulsar timing and noise model, with a Python binding called pyvela. The evidence is threefold: timing residuals computed by Vela.jl agree with those from PINT and tempo2 to within about 10 ns for an identical model; a small simulated dataset yields posterior distributions that agree with pint.bayesian; and a larger simulated binary-pulsar dataset with injected dispersion-measure variations recovers astrometric and binary parameters within 3 sigma. The larger example also shows a limitation the paper concedes: some dispersion-measure parameters do not recover their injected values, attributed to the DM model not capturing short-timescale variations.","pith_inferences":["If the prior-widening check were reported quantitatively, the cheat-prior validation could be retired for these datasets; until then, the posterior agreement shown should be read as conditional on those priors.","The linear-time likelihood evaluation opens a practical benchmark: running per-pulsar noise characterization for an entire pulsar timing array on a single workstation, which the paper does not itself demonstrate.","The sampler-agnostic interface makes model comparison by Bayesian evidence a direct next step, since nested sampling can be run without porting Vela.jl into a specific engine.","If the planned wideband and photon-domain timing are added, the same component-and-kernel architecture could unify narrowband, wideband, and high-energy pulsar timing analysis in one package."],"forward_implications":["Bayesian inference over the full non-linear timing and noise model becomes available with any MCMC or nested sampler, including Python samplers through pyvela.","The reported ~10 ns residual agreement provides an independent numerical cross-check of existing timing-model implementations, at roughly the level at which those implementations already agree with each other.","Multi-threaded, linear-time likelihood evaluation makes full Bayesian noise characterization practical on the hundreds-of-TOAs datasets typical of current pulsar timing array pulsars.","Single-pulsar PTA analyses no longer need to choose between linearized analytic marginalisation and a full non-linear sampler locked to one sampling engine.","The same modular component-and-kernel design can carry future wideband and photon-domain timing without restructuring the likelihood machinery."],"supporting_citations":[{"why":"Describes PINT, which supplies file parsing, clock corrections, and ephemerides, and is the source of the J1748–2021E example dataset and of the comparison posterior in the first validation.","marker":"(Luo et al. 2021)"},{"why":"Documents PINT's current interfaces, including the simulation tool that generated the two test datasets and the pint.bayesian fitter used as the comparison in Section 4.1.","marker":"(Susobhanan et al. 2024)"},{"why":"Defines the full non-linear timing and noise likelihood with Fourier-series Gaussian processes that Vela.jl independently implements, and the TEMPONEST approach it aims to offer an alternative to.","marker":"(Lentati et al. 2014)"},{"why":"Supplies the tempo2 timing-model equations and residual definition underlying equations (1)–(3), and is one of the two packages Vela.jl is compared against in Figure 1.","marker":"(Hobbs et al. 2006)"},{"why":"Provides the reduced-rank approximation of the covariance matrix and the Woodbury and matrix-determinant formulas that make the likelihood evaluation cheap.","marker":"(van Haasteren & Vallisneri 2014)"},{"why":"Gives the ECORR noise representation and block-diagonal formulas used in Appendix C for the linear-time correlated-noise likelihood.","marker":"(Johnson et al. 2024)"},{"why":"Provides the emcee ensemble sampler used in both Bayesian examples, demonstrating Vela.jl's compatibility with standard MCMC samplers.","marker":"(Foreman-Mackey et al. 2013)"},{"why":"The circular-analysis warning the paper cites when explaining why 'cheat' priors can be invalid if the data lack signal-to-noise, grounding the weakest assumption.","marker":"(Kriegeskorte et al. 2009)"},{"why":"Supplies the Double64 extended-precision type used to represent TOA values and rotational phases, enabling the ~10 ns residual agreement.","marker":"(Sarnoff et al. 2022)"},{"why":"Provides the InPTA DR1 narrowband data subset and epoch-wise DM measurements used to construct the J1909–3744 simulation.","marker":"(Tarafdar et al. 2022)"}],"fun_headline_variants":["Vela.jl: Bayesian pulsar timing with parallel speed","New Julia package nails pulsar timing to 10 ns","Vela.jl scales Bayesian pulsar analysis to arrays","Independent timing code passes pulsar tests","Bayesian pulsar timing gets a fast new tool"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the validation posteriors are not steered by 'cheat' priors centered on the very estimates being reproduced -- a risk the paper acknowledges and says it checked by widening the priors, but only asserts in prose, without reporting the numbers.","fun_headline_variants_meta":{"raw":{"variants":["Vela.jl: Bayesian pulsar timing with parallel speed","New Julia package nails pulsar timing to 10 ns","Vela.jl scales Bayesian pulsar analysis to arrays","Independent timing code passes pulsar tests","Bayesian pulsar timing gets a fast new tool"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000151,"raw_usage":{"total_tokens":1130,"prompt_tokens":808,"completion_tokens":322,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":424,"completion_tokens_details":{"reasoning_tokens":245}},"tokens_in":424,"tokens_out":322,"duration_ms":3119,"temperature":1.0,"reasoning_tokens":245,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-11T11:01:20.838742+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Run the same two simulated datasets through Vela.jl with broad priors anchored to physical ranges rather than to the reference fitter's point estimates (e.g., flat priors tens to thousands of times wider than the quoted uncertainties), and check whether the posterior medians and the reference-fitter agreement survive; if they shift beyond the quoted uncertainties, the circular-analysis concern is confirmed.","supporting_citations":[{"cited_title":"2014, Monthly Notices of the Royal Astronomical Society, 446, 1170, 10.1093/mnras/stu2157","cited_arxiv_id":null,"evidence_quote":"Provides the reduced-rank approximation of the covariance matrix and the Woodbury and matrix-determinant formulas that make the likelihood evaluation cheap."},{"cited_title":"2022, DoubleFloats , 1.2.2","cited_arxiv_id":null,"evidence_quote":"Supplies the Double64 extended-precision type used to represent TOA values and rotational phases, enabling the ~10 ns residual agreement."}],"review_version":1}