{"id":"e4c7eee4-186c-42c0-85ff-cb0780d497ee","arxiv_id":"2505.06459","paper_version":2,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":3,"one_line_summary":"A two-step Bayesian PINN that sets its likelihood variance to a residual error bound improves calibration on cosmological ODEs inside the training region, but is unstable outside it.","lead":"The paper adds uncertainty estimates to physics-informed neural networks with a two-step scheme: train a deterministic solution bundle, then train a Bayesian network using residual error bounds as the noise level. It tests the approach on cosmological models and uses the resulting uncertainty bands for parameter estimation from Hubble data.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The central step σ_Like(x,λ)=B(x,λ) treats a deterministic worst-case error bound as a Gaussian noise scale; unless the BNN posterior mean equals the bounded network, and unless errors are Gaussian at that scale, the reported calibration gains are enclosure artifacts rather than calibrated UQ.","rationale":"Reader's weakest assumption is correct; I agree and sharpen it. The key issue is not that B is large but that B is a deterministic certificate, so the only way to claim 'calibration' is to show the Gaussian scale and the bound coincide. The paper provides no such derivation, and Appendix I tests a different distribution. The proposed check uses the analytic solution for ΛCDM, so it is inexpensive and decisive. I do not recommend rejecting the paper: the two-step training is a plausible conservative-UQ recipe and the public repository makes reproduction feasible. The verdict should stay conditional on the coverage and scale validation.","tokens_in":28971,"tokens_out":7594,"duration_ms":81609,"concrete_test":"On the ΛCDM bundle with analytic solution, freeze the deterministic network uθdet as the point predictor and form predictive intervals with σ=B only (no BNN retraining). Compute empirical coverage at nominal 50%, 68%, 95%, and 99% levels over 200 held-out (z,Ωm) points. If the 68% interval covers near 100% of points, the EB contribution is a worst-case enclosure, not a calibrated Gaussian scale; if coverage matches nominal levels, the bound behaves like a statistical dispersion. Then repeat the same coverage computation for the actual BNN predictions and report max_x |E[uθ](x)-uθdet(x)|/B(x); if this ratio is not small, the bound is being applied to a predictor it does not certify.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Section 4.2 sets σ_Like(x,λ)=B(x,λ), where B in Eq. (11) is a guaranteed (or intended) upper bound on |u(x,λ)−uθdet(x,λ)|. A worst-case bound is not a standard deviation: a Gaussian with sigma=B places only 68% of its mass within one sigma, while the bound implies the true solution lies within one sigma 100% of the time. Thus EB predictive intervals are deterministic enclosures, not calibrated probability bands; miscalibration area can improve simply because intervals are made wide enough to contain the deterministic error, independent of whether the Gaussian error model is correct. The paper's own diagnostic (Appendix I) checks the distribution of solutions, not the distribution of errors δ=u−uθdet, so it does not validate the Gaussian assumption. There is also an object mismatch: B bounds u−uθdet, but the likelihood in Eq. (17) is centered at uθ and fitted to uθdet. The bound transfers to the BNN predictive distribution only if E[uθ|D]≈uθdet everywhere; this is never checked, and Table 2 already shows BBB+2S+EB with median RE 0.047 for ΛCDM and substantial drift in CPL/HS in the OOD tables. Whenever the posterior mean drifts, B no longer applies to the predictive interval. Finally, Algorithm 1 computes εi as a maximum over K discrete grid points, not a true supremum, so the certified-bound interpretation is weaker than stated. These gaps make the central calibration claim unverified.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper proposes a two-step Bayesian uncertainty quantification method for Physics-Informed Neural Networks (PINNs). In the first step, a deterministic PINN is trained as a Solution Bundle over equation parameters; in the second step, a Bayesian Neural Network (BNN) is trained on the deterministic network's outputs, with a Gaussian likelihood whose standard deviation is either homoscedastic or set to a residual-based deterministic error bound B(x,lambda) from Liu et al. The method is evaluated on forward problems for four cosmological models (ΛCDM, CPL, Quintessence, and HS) using miscalibration area, median relative error, and related metrics, and on inverse parameter estimation with Cosmic Chronometer data. The central claim is that using error bounds to define the heteroscedastic likelihood variance significantly improves uncertainty quantification, particularly calibration, while also providing a posterior predictive distribution for inverse problems.","tokens_in":29338,"tokens_out":3768,"duration_ms":40675,"significance":"If the central claim were established, the paper would be useful: combining rigorous PINN error bounds with BNN predictive distributions is a natural and potentially valuable idea, and the systematic comparison of NLM, BBB, and HMC within a two-step framework, together with open-source code and cosmological applications, would be a service to the community. However, the current manuscript does not validate the key modeling assumption that a deterministic worst-case error bound can be used as the standard deviation of a Gaussian likelihood. The reported calibration improvements may therefore be an artifact of making predictive intervals wide enough to contain the deterministic error, rather than evidence of well-calibrated probabilistic predictions. The paper's strengths are its clear problem formulation, the breadth of experiments, and the inclusion of inverse-problem applications, but the load-bearing inference step needs substantial additional support before the main claim can be accepted.","major_comments":[{"comment":"Setting σ_Like(x,lambda)=B(x,lambda) equates a deterministic worst-case error bound with the standard deviation of a Gaussian likelihood. A Gaussian with scale B assigns only about 68% of its mass within one standard deviation, whereas the bound in Eq. (11) guarantees |u-uθ_det|≤B with (intended) certainty. Thus the EB predictive intervals are deterministic enclosures, not calibrated probability bands, and the reported reduction in miscalibration area can be forced by construction: wider intervals trivially cover more of the true solution. The manuscript does not derive this identification, does not test the distribution of the errors δ=u-uθ_det that enter Eq. (17), and does not compare against a homoscedastic likelihood whose total width is matched to B. Without such checks, the calibration gains reported in Section 7 are not evidence for the Gaussian error model.","section":"§4.2, Eq. (17)-(18) and Eq. (11)"},{"comment":"There is an object mismatch between the object bounded by B and the object used in the likelihood. Eq. (11) bounds u(x,lambda)-uθ_det(x,lambda), but the likelihood in Eq. (17) is centered at the BNN output uθ(x,lambda) and fitted to uθ_det(x,lambda). The bound transfers to the BNN predictive distribution only if the posterior mean E[uθ|D] closely approximates uθ_det everywhere; this is never checked. Table 2 shows substantial drift in several EB configurations (e.g., BBB+2S+EB has median RE 0.047 for ΛCDM while NLM+2S+EB has 0.002), and Tables 4 and 13 show large OOD drift and extreme outliers (e.g., CPL HMC+2S+EB reaches Q100=3.88e+26 in Table 13). Whenever the BNN posterior mean drifts from uθ_det, the bound B no longer applies to the predictive interval, so the calibration claim is unsupported for those configurations.","section":"§4.1 and §4.2, Eq. (17), Table 2"},{"comment":"Algorithm 1 computes ε_i as the maximum of |r(τ)| over K discrete grid points in each partition, not a true supremum over the continuous interval. Consequently the quantity produced is an empirical approximate bound, not the certified bound stated in Eqs. (30)-(37). This weakens the interpretation of the EB intervals as guaranteed enclosures and affects the miscalibration-area analysis, since the claimed coverage property is the only formal justification for using B as a likelihood scale. The authors should either compute certified bounds (e.g., using interval arithmetic or Lipschitz-based supremum estimates) or explicitly state that B is an approximate bound and quantify the discretization error.","section":"Appendix B.3, Algorithm 1"},{"comment":"The validation of the Gaussian likelihood assumption is not targeted at the quantity that matters. Appendix I examines the distribution of solutions uθ_det across training runs and overlays Gaussians with the sample mean and variance, but the likelihood in Eq. (17) concerns the conditional distribution of the error δ=u-uθ_det (or, in the two-step setup, the residual of the BNN fit). A histogram of solutions does not test whether the errors are Gaussian with scale B, nor whether the scale varies with x and lambda as claimed. The authors should provide direct diagnostics such as standardized residuals or quantile-quantile plots for δ/B, or at least explicitly acknowledge that the Gaussian assumption is an unverified modeling choice.","section":"Appendix I"},{"comment":"The caption of Table 2 states that metrics are computed in the 'Testing Region', but this region is never defined. Tables 3 and 4 separately consider training and OOD regions, and the choice of region can materially change the reported results; the catastrophic outliers in Table 13 for CPL and HS (e.g., 9.01e+08 and 3.88e+26) are not visible in Table 2, suggesting that the testing region excludes the most extreme points. Since the central claim that 'error bounds significantly enhance calibration' is based on Table 2, the authors must define the testing region precisely, report results on a fixed held-out set that includes all evaluation points, and explain how the OOD extremes are handled in the reported miscalibration areas.","section":"Table 2 and Section 7"}],"minor_comments":[{"comment":"The phrase 'Testing Region' is undefined; please specify the exact domain and parameter ranges used for the evaluation.","section":"Table 2 caption"},{"comment":"The expression for e^{P(z)} is printed twice with an apparent typo; the first line after 'we can then simplify' repeats the previous equation rather than showing the simplified form.","section":"Appendix G.2, Eq. (51)"},{"comment":"The heading 'ADITIONAL RESULTS' should be 'ADDITIONAL RESULTS'.","section":"Appendix J heading"},{"comment":"The reference to Hoffman and Gelman lacks publication details (year and venue), and several URLs in the reference list are incomplete or contain line breaks that may be artifacts of formatting.","section":"References"},{"comment":"The baseline residual-likelihood method is not implemented for NLM, and the paper does not discuss how the absence of NLM in the baseline affects the comparison; a brief comment would improve clarity.","section":"Section 4.3 and Table 2"}],"recommendation":"major_revision","confidential_remarks":"The central idea is attractive and the experimental scope is appropriate, but the key step σ_Like=B is currently an unvalidated identification of a worst-case bound with a statistical dispersion scale. I would be willing to see a revision that adds (i) direct checks of the error distribution and empirical coverage of the predictive intervals, (ii) a comparison against non-EB intervals of matched width, (iii) verification that the posterior mean of the BNN is close to uθ_det where B is applied, and (iv) a precise definition of the testing region. If those checks do not support the calibration claim, the paper should be reframed as presenting a heuristic uncertainty heuristic rather than a calibrated UQ method."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Pablo and colleagues propose training a deterministic PINN solution bundle, then training a BNN on its outputs with the likelihood standard deviation set to a residual-based error bound B(x,λ). That specific construction is new as far as I know, and the cosmological testbed (ΛCDM, CPL, Quintessence, HS) is a reasonable place to see whether it works. Credit where due: they ship code, they report a range of uncertainty metrics, and they are frank about OOD failures and numerical blow-ups in CPL/HS. The two-step recipe itself is a useful engineering pattern.\n\nThe load-bearing step is where I get off the train. Setting σ_Like = B treats a worst-case deterministic error bound as a Gaussian noise scale. A bound is not a standard deviation; a Gaussian with σ=B gives 68% coverage at one sigma, while a bound implies 100%. The paper never checks the empirical distribution of η = u − u_θ_det. Appendix I looks at the distribution of solutions during training, not at the error process. There is also an object mismatch: B bounds u − u_θ_det, but the likelihood is centered on the BNN output u_θ and fitted to u_θ_det. That transfer works only if the BNN posterior mean is close to u_θ_det everywhere. Table 2 already shows a median RE of 0.047 for BBB+2S+EB on ΛCDM, and the OOD tables show the posterior mean drifting well outside the bound's regime. So the reported calibration gains are not as cleanly attributable to 'using the error bound' as the paper claims.\n\nOther soft spots: the 'testing region' in Table 2 is never defined; the training/OOD split in Tables 3–4 is clearer but paints a mixed picture, with some catastrophic CPL/HS results. All metrics appear to come from a single run, with no error bars or sensitivity analysis over σ_Like, σ_prior, or the partition parameters N and K.\n\nI don't think the idea is dead. The heuristic may well be useful: error bounds are cheap to compute and encode where the deterministic net is bad, so using them to inflate predictive variance in those regions is plausible. But the paper needs to either derive the mapping under explicit assumptions or, failing that, run a proper empirical calibration study of the error distribution and report coverage curves. As it stands, the central claim is unverified.\n\nWho is this for? Researchers working on UQ in scientific ML, especially those using PINNs for ODEs and parameter estimation. I'd bring it to a reading group because the methodological question is a good one. It deserves peer review, but a serious referee should demand major revision on the σ_Like=B justification.","headline":"A promising but unproven heuristic: setting the likelihood variance equal to a deterministic PINN error bound can improve calibration, yet the paper never validates that the bound behaves like a Gaussian noise scale.","tokens_in":29846,"tokens_out":4961,"would_cite":false,"duration_ms":50136,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":[],"model":"deepseek-v4-flash","headline":"The paper's central assertion is that using residual-based error bounds to define a heteroscedastic likelihood variance in a two-step Bayesian training scheme substantially improves calibration of physics-informed neural network…","keywords":["physics-informed neural networks","uncertainty quantification","Bayesian neural networks","error bounds","solution bundles","heteroscedasticity","cosmological models","inverse problems"],"falsifier":"Compute the empirical distribution of the deterministic error $\\eta(x,\\lambda)=u_{\\mathrm{true}}-\\tilde{u}_{\\theta_{\\mathrm{det}}}$ across a dense grid for the ΛCDM and CPL models, and compare it with $\\mathcal{N}(0,B(x,\\lambda)^2)$. If the true coverage of the reported 68% and 95% predictive intervals deviates substantially from nominal levels, or if replacing $B$ by a fraction $cB$ with $c<1$ does not change the miscalibration area, the central identification of the error bound as a likelihood variance is not doing the claimed work.","tokens_in":28786,"feed_emoji":"🎯","tokens_out":5570,"duration_ms":49014,"temperature":0.7,"pith_summary":"This paper claims that uncertainty estimates for physics-informed neural networks become substantially better calibrated if the noise scale of the Bayesian likelihood is tied to a worst-case residual error bound rather than a constant. The proposed two-step procedure first trains a deterministic Solution Bundle, a network that takes equation parameters as extra inputs, and then trains a Bayesian network on that deterministic solution, setting the likelihood standard deviation $\\sigma_{\\mathrm{Like}}(x,\\lambda)$ equal to the error bound $B(x,\\lambda)$. On four cosmological differential-equation models, the method reduces miscalibration area relative to residual-likelihood baselines and produces parameter posteriors for inverse problems that largely agree with literature values. The sympathetic reading is that this gives practitioners calibrated predictive distributions for PINN solutions without external data.","feed_headline":"Error bounds turn PINN uncertainties into calibrated intervals","feed_subtitle":"A two-step scheme sets Bayesian noise to a worst-case residual bound, fixing overconfident predictions for cosmological ODEs.","key_machinery":"The load-bearing object is the error bound $B(x,\\lambda)$: a deterministic scalar function computed from the differential equation's coefficients and the PINN's residual, satisfying $\\|\\eta(x,\\lambda)\\nleq B(x,\\lambda)$ for the approximation error $\\eta$. The paper's pivotal move is the substitution $\\sigma_{\\mathrm{Like}}(x,\\lambda)=B(x,\\lambda)$, turning a worst-case deterministic bound into the standard deviation of a Gaussian likelihood in the Bayesian second step. The second carrier is the Solution Bundle, a network $u_\\theta(x,\\lambda)$ trained over a range of equation parameters $\\lambda$ in one pass, which supplies the target outputs for the Bayesian network and makes inverse problems tractable because no retraining per parameter is needed.","core_discovery":"The paper's central discovery is the identification of the error bound as a heteroscedastic variance for Bayesian inference: $\\sigma_{\\mathrm{Like}}(x,\\lambda)=B(x,\\lambda)$, where $B$ is a deterministic, residual-based upper bound on the solution error. In the two-step scheme, the first step trains a deterministic Solution Bundle by minimizing residuals; the second step uses its outputs as noise-free targets for a Bayesian neural network whose likelihood variance is either constant or equal to $B$. The paper argues that the error-bound choice consistently lowers miscalibration area compared with residual-likelihood baselines and with a homoscedastic variance, while also enabling inverse-problem parameter estimation through marginalization over the predictive distribution. The authors present the result as a trade-off: accuracy is sacrificed relative to the deterministic network, but the calibrated uncertainty is the intended product.","pith_inferences":["I would expect the calibration improvement to be partly a consequence of the bound's conservativeness: $B$ is a worst-case upper bound, so intervals must cover the deterministic error by construction; a natural test is to compare against a heteroscedastic variance estimated from empirical error quantiles, which should give similar coverage with narrower intervals if the error distribution is well ","The method implicitly assumes the deterministic error is approximately Gaussian at scale $B$; the paper's own histograms show some models are noticeably non-Gaussian, so applying the same recipe to strongly non-Gaussian error distributions would require a different likelihood family.","Because $B$ is available for certain linear and weakly nonlinear ODEs, the technique should transfer to other areas where residual-based PINN bounds exist, such as linear PDEs, provided the bound is evaluated on the parameter-augmented domain."],"forward_implications":["If the central claim is right, PINN practitioners can obtain calibrated predictive intervals for forward solutions of ODE systems using only residuals and equation structure, with no external solution data.","The same predictive distributions can be fed into a Bayesian parameter-estimation loop, giving posterior means and error bars for cosmological parameters such as $\\Omega_{m,0}$ and $H_0$.","The error-bound variance is architecture-independent, so the calibration benefit should transfer to any Bayesian network used in the second step, including neural linear models, Bayes-by-backpropagation, and Hamiltonian Monte Carlo.","Calibration gains concentrate in the training region; out-of-distribution behavior is mixed, so the method is a statement about in-distribution uncertainty quantification rather than extrapolation.","A trade-off appears: intervals become wider and numerically noisier, particularly for the more complex cosmological models, so calibration at the cost of sharpness is a deliberate choice."],"supporting_citations":[{"why":"Supplies the residual-based error bound for linear ODEs that is used as the heteroscedastic likelihood standard deviation.","marker":"[Liu et al., 2022]"},{"why":"Introduces Solution Bundles, the parameter-augmented networks that produce the deterministic targets for the Bayesian second step.","marker":"[Flamant et al., 2020]"},{"why":"Establishes the two-step Bayesian training recipe of treating a deterministic network's outputs as targets for a Bayesian neural network.","marker":"[Graf et al., 2021]"},{"why":"Provides the cosmological PINN training method and the reparameterization used for the CPL model.","marker":"[Chantada et al., 2022]"},{"why":"Documents the decoupling of residuals from solution errors that motivates using error bounds instead of residuals for uncertainty.","marker":"[Wang et al., 2022]"},{"why":"Supplies the miscalibration area metric used to compare uncertainty quality.","marker":"[Chung et al., 2021]"},{"why":"Provides the MCMC sampler used to approximate posterior distributions of equation parameters in inverse problems.","marker":"[Foreman-Mackey et al., 2013]"}],"fun_headline_variants":["Error bounds calibrate PINN uncertainty","Two-step PINNs turn error bounds into confidence","Using residual bounds to fix PINN overconfidence","Calibrated PINN uncertainties via error-bound variance"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the worst-case deterministic error bound $B(x,\\lambda)$ can be treated as the standard deviation of a Gaussian likelihood; if the actual errors are not Gaussian with that scale, the reported calibration improvement is an artifact of intervals that are forced wide enough to cover the deterministic error by construction.","fun_headline_variants_meta":{"raw":{"variants":["Error bounds calibrate PINN uncertainty","Two-step PINNs turn error bounds into confidence","Using residual bounds to fix PINN overconfidence","Calibrated PINN uncertainties via error-bound variance"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000142,"raw_usage":{"total_tokens":1108,"prompt_tokens":822,"completion_tokens":286,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":438,"completion_tokens_details":{"reasoning_tokens":228}},"tokens_in":438,"tokens_out":286,"duration_ms":3301,"temperature":1.0,"reasoning_tokens":228,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-15T22:42:40.573809+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Compute the empirical distribution of the deterministic error $\\eta(x,\\lambda)=u_{\\mathrm{true}}-\\tilde{u}_{\\theta_{\\mathrm{det}}}$ across a dense grid for the ΛCDM and CPL models, and compare it with $\\mathcal{N}(0,B(x,\\lambda)^2)$. If the true coverage of the reported 68% and 95% predictive intervals deviates substantially from nominal levels, or if replacing $B$ by a fraction $cB$ with $c<1$ does not change the miscalibration area, the central identification of the error bound as a likelihood variance is not doing the claimed work.","supporting_citations":[],"review_version":1}